Understand shapefiles in different projections (notably OSGB)
[potlatch2.git] / com / gradoservice / proj4as / proj / ProjOrtho.as
1 /*******************************************************************************\r
2 NAME                             ORTHOGRAPHIC \r
3 \r
4 PURPOSE:        Transforms input longitude and latitude to Easting and\r
5                 Northing for the Orthographic projection.  The\r
6                 longitude and latitude must be in radians.  The Easting\r
7                 and Northing values will be returned in meters.\r
8 \r
9 PROGRAMMER              DATE\r
10 ----------              ----\r
11 T. Mittan               Mar, 1993\r
12 \r
13 ALGORITHM REFERENCES\r
14 \r
15 1.  Snyder, John P., "Map Projections--A Working Manual", U.S. Geological\r
16     Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United\r
17     State Government Printing Office, Washington D.C., 1987.\r
18 \r
19 2.  Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",\r
20     U.S. Geological Survey Professional Paper 1453 , United State Government\r
21     Printing Office, Washington D.C., 1989.\r
22 *******************************************************************************/\r
23 \r
24 package com.gradoservice.proj4as.proj\r
25 {\r
26         import com.gradoservice.proj4as.ProjPoint;\r
27         import com.gradoservice.proj4as.ProjConstants;\r
28         import com.gradoservice.proj4as.Datum;\r
29                         \r
30         public class ProjOrtho extends AbstractProjProjection\r
31         {\r
32                 \r
33                 private var cos_p14:Number;\r
34                 private var sin_p14:Number;\r
35                 \r
36                 public function ProjOrtho(data:ProjParams)\r
37                 {\r
38                         super(data);\r
39                 }\r
40                 \r
41                   override public function init():void\r
42                   {\r
43                     //double temp;                      /* temporary variable           */\r
44                 \r
45                     /* Place parameters in static storage for common use\r
46                       -------------------------------------------------*/;\r
47                     this.sin_p14=Math.sin(this.lat0);\r
48                     this.cos_p14=Math.cos(this.lat0);   \r
49                   }\r
50                 \r
51                 \r
52                   /* Orthographic forward equations--mapping lat,long to x,y\r
53                     ---------------------------------------------------*/\r
54                   override public function forward(p:ProjPoint):ProjPoint\r
55                   {\r
56                     var sinphi:Number, cosphi:Number;   /* sin and cos value                            */\r
57                     var dlon:Number;            /* delta longitude value                        */\r
58                     var coslon:Number;          /* cos of longitude                             */\r
59                     var ksp:Number;             /* scale factor                                 */\r
60                     var g:Number;               \r
61                     var lon:Number=p.x;\r
62                     var lat:Number=p.y; \r
63                     var x:Number,y:Number;\r
64                     /* Forward equations\r
65                       -----------------*/\r
66                     dlon = ProjConstants.adjust_lon(lon - this.long0);\r
67                 \r
68                     sinphi=Math.sin(lat);\r
69                     cosphi=Math.cos(lat);       \r
70                 \r
71                     coslon = Math.cos(dlon);\r
72                     g = this.sin_p14 * sinphi + this.cos_p14 * cosphi * coslon;\r
73                     ksp = 1.0;\r
74                     if ((g > 0) || (Math.abs(g) <= ProjConstants.EPSLN)) {\r
75                       x = this.a * ksp * cosphi * Math.sin(dlon);\r
76                       y = this.y0 + this.a * ksp * (this.cos_p14 * sinphi - this.sin_p14 * cosphi * coslon);\r
77                     } else {\r
78                       trace("orthoFwdPointError");\r
79                     }\r
80                     p.x=x;\r
81                     p.y=y;\r
82                     return p;\r
83                   }\r
84                 \r
85                 \r
86                  override public function inverse(p:ProjPoint):ProjPoint\r
87                  {\r
88                     var rh:Number;              /* height above ellipsoid                       */\r
89                     var x:Number,y:Number,z:Number;             /* angle                                        */\r
90                     var sinz:Number,cosz:Number,cosi:Number;    /* sin of z and cos of z                        */\r
91                     var temp:Number;\r
92                     var con:Number;\r
93                     var lon:Number , lat:Number;\r
94                     /* Inverse equations\r
95                       -----------------*/\r
96                     p.x -= this.x0;\r
97                     p.y -= this.y0;\r
98                     rh = Math.sqrt(p.x * p.x + p.y * p.y);\r
99                     if (rh > this.a + .0000001) {\r
100                       trace("orthoInvDataError");\r
101                     }\r
102                     z = ProjConstants.asinz(rh / this.a);\r
103                 \r
104                     sinz=Math.sin(z);\r
105                     cosi=Math.cos(z);\r
106                 \r
107                     lon = this.long0;\r
108                     if (Math.abs(rh) <= ProjConstants.EPSLN) {\r
109                       lat = this.lat0; \r
110                     }\r
111                     lat = ProjConstants.asinz(cosz * this.sin_p14 + (y * sinz * this.cos_p14)/rh);\r
112                     con = Math.abs(lat0) - ProjConstants.HALF_PI;\r
113                     if (Math.abs(con) <= ProjConstants.EPSLN) {\r
114                        if (this.lat0 >= 0) {\r
115                           lon = ProjConstants.adjust_lon(this.long0 + Math.atan2(p.x, -p.y));\r
116                        } else {\r
117                           lon = ProjConstants.adjust_lon(this.long0 -Math.atan2(-p.x, p.y));\r
118                        }\r
119                     }\r
120                     con = cosz - this.sin_p14 * Math.sin(lat);\r
121                     if ((Math.abs(con) >= ProjConstants.EPSLN) || (Math.abs(x) >= ProjConstants.EPSLN)) {\r
122                        lon = ProjConstants.adjust_lon(this.long0 + Math.atan2((p.x * sinz * this.cos_p14), (con * rh)));\r
123                     }\r
124                     p.x=lon;\r
125                     p.y=lat;\r
126                     return p;\r
127                   }\r
128                 \r
129                 \r
130         }\r
131 }