Ignore:
Timestamp:
2008-12-23T15:07:05+01:00 (11 years ago)
Author:
stoecker
Message:

removed usage of tab stops

Location:
trunk/src/org/openstreetmap/josm/data/projection
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/org/openstreetmap/josm/data/projection/Epsg4326.java

    r1032 r1169  
    1212public class Epsg4326 implements Projection {
    1313
    14         public EastNorth latlon2eastNorth(LatLon p) {
    15                 return new EastNorth(p.lon(), p.lat());
    16         }
     14    public EastNorth latlon2eastNorth(LatLon p) {
     15        return new EastNorth(p.lon(), p.lat());
     16    }
    1717
    18         public LatLon eastNorth2latlon(EastNorth p) {
    19                 return new LatLon(p.north(), p.east());
    20         }
     18    public LatLon eastNorth2latlon(EastNorth p) {
     19        return new LatLon(p.north(), p.east());
     20    }
    2121
    22         @Override public String toString() {
    23                 return tr("EPSG:4326");
    24         }
     22    @Override public String toString() {
     23        return tr("EPSG:4326");
     24    }
    2525
    2626    public String getCacheDirectoryName() {
     
    2828    }
    2929
    30         public double scaleFactor() {
    31             return 1.0/360;
     30    public double scaleFactor() {
     31        return 1.0/360;
    3232    }
    3333
    34         @Override public boolean equals(Object o) {
    35                 return o instanceof Epsg4326;
    36         }
     34    @Override public boolean equals(Object o) {
     35        return o instanceof Epsg4326;
     36    }
    3737
    38         @Override public int hashCode() {
    39                 return Epsg4326.class.hashCode();
    40         }
     38    @Override public int hashCode() {
     39        return Epsg4326.class.hashCode();
     40    }
    4141}
  • trunk/src/org/openstreetmap/josm/data/projection/Lambert.java

    r865 r1169  
    1414
    1515public class Lambert implements Projection {
    16         /**
    17         * Lambert I, II, III, and IV projection exponents
    18         */
    19         public static final double n[] = { 0.7604059656, 0.7289686274, 0.6959127966, 0.6712679322 };
    20 
    21         /**
    22         * Lambert I, II, III, and IV projection constants
    23         */
    24         public static final double c[] = { 11603796.98, 11745793.39, 11947992.52, 12136281.99 };
    25 
    26         /**
    27         * Lambert I, II, III, and IV false east
    28         */
    29         public static final double Xs[] = { 600000.0, 600000.0, 600000.0, 234.358 };
    30 
    31         /**
    32         * Lambert I, II, III, and IV false north
    33         */
    34         public static final double Ys[] = { 5657616.674, 6199695.768, 6791905.085, 7239161.542 };
    35 
    36         /**
    37         * Lambert I, II, III, and IV longitudinal offset to Greenwich meridian
    38         */
    39         public static final double lg0 = 0.04079234433198; // 2deg20'14.025"
    40 
    41         /**
    42         * precision in iterative schema
    43         */
    44 
    45         public static final double epsilon = 1e-11;
    46 
    47         /**
    48         * France is divided in 4 Lambert projection zones (1,2,3 + 4th for Corsica)
    49         */
    50         public static final double cMaxLatZone1 = Math.toRadians(57 * 0.9);
    51 
    52         public static final double zoneLimits[] = { Math.toRadians(53.5 * 0.9), // between Zone 1 and Zone 2 (in grad *0.9)
    53                         Math.toRadians(50.5 * 0.9), // between Zone 2 and Zone 3
    54                         Math.toRadians(47.51963 * 0.9), // between Zone 3 and Zone 4
    55                         Math.toRadians(46.17821 * 0.9) };// lowest latitude of Zone 4
    56 
    57         public static final double cMinLonZones = Math.toRadians(-4.9074074074074059 * 0.9);
    58 
    59         public static final double cMaxLonZones = Math.toRadians(10.2 * 0.9);
    60 
    61         /**
    62         *  Because josm cannot work correctly if two zones are displayed, we allow some overlapping
    63         */
    64         public static final double cMaxOverlappingZones = Math.toRadians(1.5 * 0.9);
    65 
    66         public static int layoutZone = -1;
    67 
    68         private static int currentZone = 0;
    69 
    70         private static boolean dontDisplayErrors = false;
    71 
    72         /**
    73         * @param p  WGS84 lat/lon (ellipsoid GRS80) (in degree)
    74         * @return eastnorth projection in Lambert Zone (ellipsoid Clark)
    75         */
    76         public EastNorth latlon2eastNorth(LatLon p) {
    77                 // translate ellipsoid GRS80 (WGS83) => Clark
    78                 LatLon geo = GRS802Clark(p);
    79                 double lt = geo.lat(); // in radian
    80                 double lg = geo.lon();
    81 
    82                 // check if longitude and latitude are inside the french Lambert zones
    83                 currentZone = 0;
    84                 boolean outOfLambertZones = false;
    85                 if (lt >= zoneLimits[3] && lt <= cMaxLatZone1 && lg >= cMinLonZones && lg <= cMaxLonZones) {
    86                         // zone I
    87                         if (lt > zoneLimits[0])
    88                                 currentZone = 0;
    89                         // zone II
    90                         else if (lt > zoneLimits[1])
    91                                 currentZone = 1;
    92                         // zone III
    93                         else if (lt > zoneLimits[2])
    94                                 currentZone = 2;
    95                         // zone III or IV
    96                         else if (lt > zoneLimits[3])
    97                                 // Note: zone IV is dedicated to Corsica island and extends from 47.8 to
    98                                 // 45.9 degrees of latitude. There is an overlap with zone III that can be
    99                                 // solved only with longitude (covers Corsica if lon > 7.2 degree)
    100                                 if (lg < Math.toRadians(8 * 0.9))
    101                                         currentZone = 2;
    102                                 else
    103                                         currentZone = 3;
    104                 } else {
    105                         outOfLambertZones = true; // possible when MAX_LAT is used
    106                         if (p.lat() != 0 && Math.abs(p.lat()) != Projection.MAX_LAT
    107                                         && p.lon() != 0 && Math.abs(p.lon()) != Projection.MAX_LON
    108                                         && dontDisplayErrors == false) {
    109                                 JOptionPane.showMessageDialog(Main.parent,
    110                                                 tr("The projection \"{0}\" is designed for\n"
    111                                                 + "latitudes between 46.1° and 57° only.\n"
    112                                                 + "Use another projection system if you are not using\n"
    113                                                 + "a french WMS server.\n"
    114                                                 + "Do not upload any data after this message.", this.toString()));
    115                                 dontDisplayErrors = true;
    116                         }
    117                 }
    118                 if (!outOfLambertZones) {
    119                         if (layoutZone == -1) {
    120                                 layoutZone = currentZone;
    121                                 dontDisplayErrors = false;
    122                         } else if (layoutZone != currentZone) {
    123                                 if ((currentZone < layoutZone && Math.abs(zoneLimits[currentZone] - lt) > cMaxOverlappingZones)
    124                                                 || (currentZone > layoutZone && Math.abs(zoneLimits[layoutZone] - lt) > cMaxOverlappingZones)) {
    125                                         JOptionPane.showMessageDialog(Main.parent,
    126                                                                         tr("IMPORTANT : data positionned far away from\n"
    127                                                                                         + "the current Lambert zone limits.\n"
    128                                                                                         + "Do not upload any data after this message.\n"
    129                                                                                         + "Undo your last action, Save your work \n"
    130                                                                                         + "and Start a new layer on the new zone."));
    131                                         layoutZone = -1;
    132                                         dontDisplayErrors = true;
    133                                 } else {
    134                                         System.out.println("temporarily extend Lambert zone " + layoutZone + " projection at lat,lon:"
    135                                                         + lt + "," + lg);
    136                                 }
    137                         }
    138                 }
    139                 if (layoutZone == -1) {
    140                         return ConicProjection(lt, lg, Xs[currentZone], Ys[currentZone], c[currentZone], n[currentZone]);
    141                 } // else
    142                 return ConicProjection(lt, lg, Xs[layoutZone], Ys[layoutZone], c[layoutZone], n[layoutZone]);
    143         }
    144 
    145         public LatLon eastNorth2latlon(EastNorth p) {
    146                 LatLon geo;
    147                 if (layoutZone == -1)
    148                         // possible until the Lambert zone is determined by latlon2eastNorth() with a valid LatLon
    149                         geo = Geographic(p, Xs[currentZone], Ys[currentZone], c[currentZone], n[currentZone]);
    150                 else
    151                         geo = Geographic(p, Xs[layoutZone], Ys[layoutZone], c[layoutZone], n[layoutZone]);
    152                 // translate ellipsoid Clark => GRS80 (WGS83)
    153                 LatLon wgs = Clark2GRS80(geo);
    154                 return new LatLon(Math.toDegrees(wgs.lat()), Math.toDegrees(wgs.lon()));
    155         }
    156 
    157         @Override public String toString() {
    158                 return tr("Lambert Zone (France)");
    159         }
    160 
    161         public String getCacheDirectoryName() {
    162                 return "lambert";
    163         }
    164 
    165         public double scaleFactor() {
    166                 return 1.0 / 360;
    167         }
    168 
    169         @Override
    170         public boolean equals(Object o) {
    171                 return o instanceof Lambert;
    172         }
    173 
    174         @Override
    175         public int hashCode() {
    176                 return Lambert.class.hashCode();
    177         }
    178 
    179         /**
    180         * Initializes from geographic coordinates. Note that reference ellipsoid
    181         * used by Lambert is the Clark ellipsoid.
    182         *
    183         * @param lat latitude in grad
    184         * @param lon longitude in grad
    185         * @param Xs  false east (coordinate system origin) in meters
    186         * @param Ys  false north (coordinate system origin) in meters
    187         * @param c   projection constant
    188         * @param n   projection exponent
    189         * @return EastNorth projected coordinates in meter
    190         */
    191         private EastNorth ConicProjection(double lat, double lon, double Xs, double Ys, double c, double n) {
    192                 double eslt = Ellipsoid.clarke.e * Math.sin(lat);
    193                 double l = Math.log(Math.tan(Math.PI / 4.0 + (lat / 2.0))
    194                                 * Math.pow((1.0 - eslt) / (1.0 + eslt), Ellipsoid.clarke.e / 2.0));
    195                 double east = Xs + c * Math.exp(-n * l) * Math.sin(n * (lon - lg0));
    196                 double north = Ys - c * Math.exp(-n * l) * Math.cos(n * (lon - lg0));
    197                 return new EastNorth(east, north);
    198         }
    199 
    200         /**
    201         * Initializes from projected coordinates (conic projection). Note that
    202         * reference ellipsoid used by Lambert is Clark
    203         *
    204         * @param coord projected coordinates pair in meters
    205         * @param Xs    false east (coordinate system origin) in meters
    206         * @param Ys    false north (coordinate system origin) in meters
    207         * @param c     projection constant
    208         * @param n     projection exponent
    209         * @return LatLon in radian
    210         */
    211         private LatLon Geographic(EastNorth eastNorth, double Xs, double Ys, double c, double n) {
    212                 double dx = eastNorth.east() - Xs;
    213                 double dy = Ys - eastNorth.north();
    214                 double R = Math.sqrt(dx * dx + dy * dy);
    215                 double gamma = Math.atan(dx / dy);
    216                 double l = -1.0 / n * Math.log(Math.abs(R / c));
    217                 l = Math.exp(l);
    218                 double lon = lg0 + gamma / n;
    219                 double lat = 2.0 * Math.atan(l) - Math.PI / 2.0;
    220                 double delta = 1.0;
    221                 while (delta > epsilon) {
    222                         double eslt = Ellipsoid.clarke.e * Math.sin(lat);
    223                         double nlt = 2.0 * Math.atan(Math.pow((1.0 + eslt) / (1.0 - eslt), Ellipsoid.clarke.e / 2.0) * l) - Math.PI
    224                                         / 2.0;
    225                         delta = Math.abs(nlt - lat);
    226                         lat = nlt;
    227                 }
    228                 return new LatLon(lat, lon); // in radian
    229         }
    230 
    231         /**
    232         * Translate latitude/longitude in WGS84, (ellipsoid GRS80) to Lambert
    233         * geographic, (ellipsoid Clark)
    234         *
    235         * @param wgs
    236         * @return
    237         */
    238         private LatLon GRS802Clark(LatLon wgs) {
    239                 double lat = Math.toRadians(wgs.lat()); // degree to radian
    240                 double lon = Math.toRadians(wgs.lon());
    241                 // WGS84 geographic => WGS84 cartesian
    242                 double N = Ellipsoid.GRS80.a / (Math.sqrt(1.0 - Ellipsoid.GRS80.e2 * Math.sin(lat) * Math.sin(lat)));
    243                 double X = (N/* +height */) * Math.cos(lat) * Math.cos(lon);
    244                 double Y = (N/* +height */) * Math.cos(lat) * Math.sin(lon);
    245                 double Z = (N * (1.0 - Ellipsoid.GRS80.e2)/* + height */) * Math.sin(lat);
    246                 // WGS84 => Lambert ellipsoide similarity
    247                 X += 168.0;
    248                 Y += 60.0;
    249                 Z += -320.0;
    250                 // Lambert cartesian => Lambert geographic
    251                 return Geographic(X, Y, Z, Ellipsoid.clarke);
    252         }
    253 
    254         /**
    255         * @param lambert
    256         * @return
    257         */
    258         private LatLon Clark2GRS80(LatLon lambert) {
    259                 double lat = lambert.lat(); // in radian
    260                 double lon = lambert.lon();
    261                 // Lambert geographic => Lambert cartesian
    262                 double N = Ellipsoid.clarke.a / (Math.sqrt(1.0 - Ellipsoid.clarke.e2 * Math.sin(lat) * Math.sin(lat)));
    263                 double X = (N/* +height */) * Math.cos(lat) * Math.cos(lon);
    264                 double Y = (N/* +height */) * Math.cos(lat) * Math.sin(lon);
    265                 double Z = (N * (1.0 - Ellipsoid.clarke.e2)/* + height */) * Math.sin(lat);
    266                 // Lambert => WGS84 ellipsoide similarity
    267                 X += -168.0;
    268                 Y += -60.0;
    269                 Z += 320.0;
    270                 // WGS84 cartesian => WGS84 geographic
    271                 return Geographic(X, Y, Z, Ellipsoid.GRS80);
    272         }
    273 
    274         /**
    275         * initializes from cartesian coordinates
    276         *
    277         * @param X
    278         *            1st coordinate in meters
    279         * @param Y
    280         *            2nd coordinate in meters
    281         * @param Z
    282         *            3rd coordinate in meters
    283         * @param ell
    284         *            reference ellipsoid
    285         */
    286         private LatLon Geographic(double X, double Y, double Z, Ellipsoid ell) {
    287                 double norm = Math.sqrt(X * X + Y * Y);
    288                 double lg = 2.0 * Math.atan(Y / (X + norm));
    289                 double lt = Math.atan(Z / (norm * (1.0 - (ell.a * ell.e2 / Math.sqrt(X * X + Y * Y + Z * Z)))));
    290                 double delta = 1.0;
    291                 while (delta > epsilon) {
    292                         double s2 = Math.sin(lt);
    293                         s2 *= s2;
    294                         double l = Math.atan((Z / norm)
    295                                 / (1.0 - (ell.a * ell.e2 * Math.cos(lt) / (norm * Math.sqrt(1.0 - ell.e2 * s2)))));
    296                         delta = Math.abs(l - lt);
    297                         lt = l;
    298                 }
    299                 double s2 = Math.sin(lt);
    300                 s2 *= s2;
    301                 // h = norm / Math.cos(lt) - ell.a / Math.sqrt(1.0 - ell.e2 * s2);
    302                 return new LatLon(lt, lg);
    303         }
     16    /**
     17    * Lambert I, II, III, and IV projection exponents
     18    */
     19    public static final double n[] = { 0.7604059656, 0.7289686274, 0.6959127966, 0.6712679322 };
     20
     21    /**
     22    * Lambert I, II, III, and IV projection constants
     23    */
     24    public static final double c[] = { 11603796.98, 11745793.39, 11947992.52, 12136281.99 };
     25
     26    /**
     27    * Lambert I, II, III, and IV false east
     28    */
     29    public static final double Xs[] = { 600000.0, 600000.0, 600000.0, 234.358 };
     30
     31    /**
     32    * Lambert I, II, III, and IV false north
     33    */
     34    public static final double Ys[] = { 5657616.674, 6199695.768, 6791905.085, 7239161.542 };
     35
     36    /**
     37    * Lambert I, II, III, and IV longitudinal offset to Greenwich meridian
     38    */
     39    public static final double lg0 = 0.04079234433198; // 2deg20'14.025"
     40
     41    /**
     42    * precision in iterative schema
     43    */
     44
     45    public static final double epsilon = 1e-11;
     46
     47    /**
     48    * France is divided in 4 Lambert projection zones (1,2,3 + 4th for Corsica)
     49    */
     50    public static final double cMaxLatZone1 = Math.toRadians(57 * 0.9);
     51
     52    public static final double zoneLimits[] = { Math.toRadians(53.5 * 0.9), // between Zone 1 and Zone 2 (in grad *0.9)
     53            Math.toRadians(50.5 * 0.9), // between Zone 2 and Zone 3
     54            Math.toRadians(47.51963 * 0.9), // between Zone 3 and Zone 4
     55            Math.toRadians(46.17821 * 0.9) };// lowest latitude of Zone 4
     56
     57    public static final double cMinLonZones = Math.toRadians(-4.9074074074074059 * 0.9);
     58
     59    public static final double cMaxLonZones = Math.toRadians(10.2 * 0.9);
     60
     61    /**
     62    *  Because josm cannot work correctly if two zones are displayed, we allow some overlapping
     63    */
     64    public static final double cMaxOverlappingZones = Math.toRadians(1.5 * 0.9);
     65
     66    public static int layoutZone = -1;
     67
     68    private static int currentZone = 0;
     69
     70    private static boolean dontDisplayErrors = false;
     71
     72    /**
     73    * @param p  WGS84 lat/lon (ellipsoid GRS80) (in degree)
     74    * @return eastnorth projection in Lambert Zone (ellipsoid Clark)
     75    */
     76    public EastNorth latlon2eastNorth(LatLon p) {
     77        // translate ellipsoid GRS80 (WGS83) => Clark
     78        LatLon geo = GRS802Clark(p);
     79        double lt = geo.lat(); // in radian
     80        double lg = geo.lon();
     81
     82        // check if longitude and latitude are inside the french Lambert zones
     83        currentZone = 0;
     84        boolean outOfLambertZones = false;
     85        if (lt >= zoneLimits[3] && lt <= cMaxLatZone1 && lg >= cMinLonZones && lg <= cMaxLonZones) {
     86            // zone I
     87            if (lt > zoneLimits[0])
     88                currentZone = 0;
     89            // zone II
     90            else if (lt > zoneLimits[1])
     91                currentZone = 1;
     92            // zone III
     93            else if (lt > zoneLimits[2])
     94                currentZone = 2;
     95            // zone III or IV
     96            else if (lt > zoneLimits[3])
     97                // Note: zone IV is dedicated to Corsica island and extends from 47.8 to
     98                // 45.9 degrees of latitude. There is an overlap with zone III that can be
     99                // solved only with longitude (covers Corsica if lon > 7.2 degree)
     100                if (lg < Math.toRadians(8 * 0.9))
     101                    currentZone = 2;
     102                else
     103                    currentZone = 3;
     104        } else {
     105            outOfLambertZones = true; // possible when MAX_LAT is used
     106            if (p.lat() != 0 && Math.abs(p.lat()) != Projection.MAX_LAT
     107                    && p.lon() != 0 && Math.abs(p.lon()) != Projection.MAX_LON
     108                    && dontDisplayErrors == false) {
     109                JOptionPane.showMessageDialog(Main.parent,
     110                        tr("The projection \"{0}\" is designed for\n"
     111                        + "latitudes between 46.1° and 57° only.\n"
     112                        + "Use another projection system if you are not using\n"
     113                        + "a french WMS server.\n"
     114                        + "Do not upload any data after this message.", this.toString()));
     115                dontDisplayErrors = true;
     116            }
     117        }
     118        if (!outOfLambertZones) {
     119            if (layoutZone == -1) {
     120                layoutZone = currentZone;
     121                dontDisplayErrors = false;
     122            } else if (layoutZone != currentZone) {
     123                if ((currentZone < layoutZone && Math.abs(zoneLimits[currentZone] - lt) > cMaxOverlappingZones)
     124                        || (currentZone > layoutZone && Math.abs(zoneLimits[layoutZone] - lt) > cMaxOverlappingZones)) {
     125                    JOptionPane.showMessageDialog(Main.parent,
     126                                    tr("IMPORTANT : data positionned far away from\n"
     127                                            + "the current Lambert zone limits.\n"
     128                                            + "Do not upload any data after this message.\n"
     129                                            + "Undo your last action, Save your work \n"
     130                                            + "and Start a new layer on the new zone."));
     131                    layoutZone = -1;
     132                    dontDisplayErrors = true;
     133                } else {
     134                    System.out.println("temporarily extend Lambert zone " + layoutZone + " projection at lat,lon:"
     135                            + lt + "," + lg);
     136                }
     137            }
     138        }
     139        if (layoutZone == -1) {
     140            return ConicProjection(lt, lg, Xs[currentZone], Ys[currentZone], c[currentZone], n[currentZone]);
     141        } // else
     142        return ConicProjection(lt, lg, Xs[layoutZone], Ys[layoutZone], c[layoutZone], n[layoutZone]);
     143    }
     144
     145    public LatLon eastNorth2latlon(EastNorth p) {
     146        LatLon geo;
     147        if (layoutZone == -1)
     148            // possible until the Lambert zone is determined by latlon2eastNorth() with a valid LatLon
     149            geo = Geographic(p, Xs[currentZone], Ys[currentZone], c[currentZone], n[currentZone]);
     150        else
     151            geo = Geographic(p, Xs[layoutZone], Ys[layoutZone], c[layoutZone], n[layoutZone]);
     152        // translate ellipsoid Clark => GRS80 (WGS83)
     153        LatLon wgs = Clark2GRS80(geo);
     154        return new LatLon(Math.toDegrees(wgs.lat()), Math.toDegrees(wgs.lon()));
     155    }
     156
     157    @Override public String toString() {
     158        return tr("Lambert Zone (France)");
     159    }
     160
     161    public String getCacheDirectoryName() {
     162        return "lambert";
     163    }
     164
     165    public double scaleFactor() {
     166        return 1.0 / 360;
     167    }
     168
     169    @Override
     170    public boolean equals(Object o) {
     171        return o instanceof Lambert;
     172    }
     173
     174    @Override
     175    public int hashCode() {
     176        return Lambert.class.hashCode();
     177    }
     178
     179    /**
     180    * Initializes from geographic coordinates. Note that reference ellipsoid
     181    * used by Lambert is the Clark ellipsoid.
     182    *
     183    * @param lat latitude in grad
     184    * @param lon longitude in grad
     185    * @param Xs  false east (coordinate system origin) in meters
     186    * @param Ys  false north (coordinate system origin) in meters
     187    * @param c   projection constant
     188    * @param n   projection exponent
     189    * @return EastNorth projected coordinates in meter
     190    */
     191    private EastNorth ConicProjection(double lat, double lon, double Xs, double Ys, double c, double n) {
     192        double eslt = Ellipsoid.clarke.e * Math.sin(lat);
     193        double l = Math.log(Math.tan(Math.PI / 4.0 + (lat / 2.0))
     194                * Math.pow((1.0 - eslt) / (1.0 + eslt), Ellipsoid.clarke.e / 2.0));
     195        double east = Xs + c * Math.exp(-n * l) * Math.sin(n * (lon - lg0));
     196        double north = Ys - c * Math.exp(-n * l) * Math.cos(n * (lon - lg0));
     197        return new EastNorth(east, north);
     198    }
     199
     200    /**
     201    * Initializes from projected coordinates (conic projection). Note that
     202    * reference ellipsoid used by Lambert is Clark
     203    *
     204    * @param coord projected coordinates pair in meters
     205    * @param Xs    false east (coordinate system origin) in meters
     206    * @param Ys    false north (coordinate system origin) in meters
     207    * @param c     projection constant
     208    * @param n     projection exponent
     209    * @return LatLon in radian
     210    */
     211    private LatLon Geographic(EastNorth eastNorth, double Xs, double Ys, double c, double n) {
     212        double dx = eastNorth.east() - Xs;
     213        double dy = Ys - eastNorth.north();
     214        double R = Math.sqrt(dx * dx + dy * dy);
     215        double gamma = Math.atan(dx / dy);
     216        double l = -1.0 / n * Math.log(Math.abs(R / c));
     217        l = Math.exp(l);
     218        double lon = lg0 + gamma / n;
     219        double lat = 2.0 * Math.atan(l) - Math.PI / 2.0;
     220        double delta = 1.0;
     221        while (delta > epsilon) {
     222            double eslt = Ellipsoid.clarke.e * Math.sin(lat);
     223            double nlt = 2.0 * Math.atan(Math.pow((1.0 + eslt) / (1.0 - eslt), Ellipsoid.clarke.e / 2.0) * l) - Math.PI
     224                    / 2.0;
     225            delta = Math.abs(nlt - lat);
     226            lat = nlt;
     227        }
     228        return new LatLon(lat, lon); // in radian
     229    }
     230
     231    /**
     232    * Translate latitude/longitude in WGS84, (ellipsoid GRS80) to Lambert
     233    * geographic, (ellipsoid Clark)
     234    *
     235    * @param wgs
     236    * @return
     237    */
     238    private LatLon GRS802Clark(LatLon wgs) {
     239        double lat = Math.toRadians(wgs.lat()); // degree to radian
     240        double lon = Math.toRadians(wgs.lon());
     241        // WGS84 geographic => WGS84 cartesian
     242        double N = Ellipsoid.GRS80.a / (Math.sqrt(1.0 - Ellipsoid.GRS80.e2 * Math.sin(lat) * Math.sin(lat)));
     243        double X = (N/* +height */) * Math.cos(lat) * Math.cos(lon);
     244        double Y = (N/* +height */) * Math.cos(lat) * Math.sin(lon);
     245        double Z = (N * (1.0 - Ellipsoid.GRS80.e2)/* + height */) * Math.sin(lat);
     246        // WGS84 => Lambert ellipsoide similarity
     247        X += 168.0;
     248        Y += 60.0;
     249        Z += -320.0;
     250        // Lambert cartesian => Lambert geographic
     251        return Geographic(X, Y, Z, Ellipsoid.clarke);
     252    }
     253
     254    /**
     255    * @param lambert
     256    * @return
     257    */
     258    private LatLon Clark2GRS80(LatLon lambert) {
     259        double lat = lambert.lat(); // in radian
     260        double lon = lambert.lon();
     261        // Lambert geographic => Lambert cartesian
     262        double N = Ellipsoid.clarke.a / (Math.sqrt(1.0 - Ellipsoid.clarke.e2 * Math.sin(lat) * Math.sin(lat)));
     263        double X = (N/* +height */) * Math.cos(lat) * Math.cos(lon);
     264        double Y = (N/* +height */) * Math.cos(lat) * Math.sin(lon);
     265        double Z = (N * (1.0 - Ellipsoid.clarke.e2)/* + height */) * Math.sin(lat);
     266        // Lambert => WGS84 ellipsoide similarity
     267        X += -168.0;
     268        Y += -60.0;
     269        Z += 320.0;
     270        // WGS84 cartesian => WGS84 geographic
     271        return Geographic(X, Y, Z, Ellipsoid.GRS80);
     272    }
     273
     274    /**
     275    * initializes from cartesian coordinates
     276    *
     277    * @param X
     278    *            1st coordinate in meters
     279    * @param Y
     280    *            2nd coordinate in meters
     281    * @param Z
     282    *            3rd coordinate in meters
     283    * @param ell
     284    *            reference ellipsoid
     285    */
     286    private LatLon Geographic(double X, double Y, double Z, Ellipsoid ell) {
     287        double norm = Math.sqrt(X * X + Y * Y);
     288        double lg = 2.0 * Math.atan(Y / (X + norm));
     289        double lt = Math.atan(Z / (norm * (1.0 - (ell.a * ell.e2 / Math.sqrt(X * X + Y * Y + Z * Z)))));
     290        double delta = 1.0;
     291        while (delta > epsilon) {
     292            double s2 = Math.sin(lt);
     293            s2 *= s2;
     294            double l = Math.atan((Z / norm)
     295                / (1.0 - (ell.a * ell.e2 * Math.cos(lt) / (norm * Math.sqrt(1.0 - ell.e2 * s2)))));
     296            delta = Math.abs(l - lt);
     297            lt = l;
     298        }
     299        double s2 = Math.sin(lt);
     300        s2 *= s2;
     301        // h = norm / Math.cos(lt) - ell.a / Math.sqrt(1.0 - ell.e2 * s2);
     302        return new LatLon(lt, lg);
     303    }
    304304
    305305}
  • trunk/src/org/openstreetmap/josm/data/projection/Mercator.java

    r1032 r1169  
    1010 * Implement Mercator Projection code, coded after documentation
    1111 * from wikipedia.
    12  * 
    13  * The center of the mercator projection is always the 0 grad 
     12 *
     13 * The center of the mercator projection is always the 0 grad
    1414 * coordinate.
    15  * 
     15 *
    1616 * @author imi
    1717 */
    1818public class Mercator implements Projection {
    1919
    20         public EastNorth latlon2eastNorth(LatLon p) {
    21                 return new EastNorth(
    22                         p.lon()*Math.PI/180,
    23                         Math.log(Math.tan(Math.PI/4+p.lat()*Math.PI/360)));
    24         }
     20    public EastNorth latlon2eastNorth(LatLon p) {
     21        return new EastNorth(
     22            p.lon()*Math.PI/180,
     23            Math.log(Math.tan(Math.PI/4+p.lat()*Math.PI/360)));
     24    }
    2525
    26         public LatLon eastNorth2latlon(EastNorth p) {
    27                 return new LatLon(
    28                         Math.atan(Math.sinh(p.north()))*180/Math.PI,
    29                         p.east()*180/Math.PI);
    30         }
     26    public LatLon eastNorth2latlon(EastNorth p) {
     27        return new LatLon(
     28            Math.atan(Math.sinh(p.north()))*180/Math.PI,
     29            p.east()*180/Math.PI);
     30    }
    3131
    32         @Override public String toString() {
    33                 return tr("Mercator");
    34         }
     32    @Override public String toString() {
     33        return tr("Mercator");
     34    }
    3535
    3636    public String getCacheDirectoryName() {
     
    3838    }
    3939
    40         public double scaleFactor() {
    41             return 1/Math.PI/2;
     40    public double scaleFactor() {
     41        return 1/Math.PI/2;
    4242    }
    4343
    44         @Override public boolean equals(Object o) {
    45                 return o instanceof Mercator;
    46         }
     44    @Override public boolean equals(Object o) {
     45        return o instanceof Mercator;
     46    }
    4747
    48         @Override public int hashCode() {
    49                 return Mercator.class.hashCode();
    50         }
     48    @Override public int hashCode() {
     49        return Mercator.class.hashCode();
     50    }
    5151}
  • trunk/src/org/openstreetmap/josm/data/projection/Projection.java

    r656 r1169  
    66
    77/**
    8  * Classes implementing this are able to convert lat/lon values to 
     8 * Classes implementing this are able to convert lat/lon values to
    99 * planar screen coordinates.
    10  * 
     10 *
    1111 * @author imi
    1212 */
    1313public interface Projection {
    1414
    15         /**
    16          * Maximum latitude representable.
    17          */
    18         public static final double MAX_LAT = 85.05112877980659; // Mercator squares the world
    19        
    20         /**
    21          * Maximum longditude representable.
    22          */
    23         public static final double MAX_LON = 180;
     15    /**
     16     * Maximum latitude representable.
     17     */
     18    public static final double MAX_LAT = 85.05112877980659; // Mercator squares the world
    2419
    25         /**
    26          * Minimum difference in location to not be represented as the same position.
    27         */
    28         public static final double MAX_SERVER_PRECISION = 1e12;
     20    /**
     21     * Maximum longditude representable.
     22    */
     23    public static final double MAX_LON = 180;
    2924
    30         /**
    31          * List of all available projections.
    32          */
    33         public static Projection[] allProjections = new Projection[]{
    34                 new Epsg4326(),
    35                 new Mercator(),
    36                 new Lambert()
    37         };
    38        
    39         /**
    40          * Convert from lat/lon to northing/easting.
    41          *
    42          * @param p             The geo point to convert. x/y members of the point are filled.
    43          */
    44         EastNorth latlon2eastNorth(LatLon p);
    45        
    46         /**
    47          * Convert from norting/easting to lat/lon.
    48          *
    49          * @param p             The geo point to convert. lat/lon members of the point are filled.
    50          */
    51         LatLon eastNorth2latlon(EastNorth p);
     25    /**
     26     * Minimum difference in location to not be represented as the same position.
     27     */
     28    public static final double MAX_SERVER_PRECISION = 1e12;
    5229
    53         /**
    54          * Describe the projection converter in one or two words.
    55          */
    56         String toString();
    57    
     30    /**
     31     * List of all available projections.
     32     */
     33    public static Projection[] allProjections = new Projection[]{
     34        new Epsg4326(),
     35        new Mercator(),
     36        new Lambert()
     37    };
     38
     39    /**
     40     * Convert from lat/lon to northing/easting.
     41     *
     42     * @param p     The geo point to convert. x/y members of the point are filled.
     43     */
     44    EastNorth latlon2eastNorth(LatLon p);
     45
     46    /**
     47     * Convert from norting/easting to lat/lon.
     48     *
     49     * @param p     The geo point to convert. lat/lon members of the point are filled.
     50     */
     51    LatLon eastNorth2latlon(EastNorth p);
     52
     53    /**
     54     * Describe the projection converter in one or two words.
     55     */
     56    String toString();
     57
    5858    /**
    5959     * Get a filename compatible string (for the cache directory)
    6060     */
    6161    String getCacheDirectoryName();
    62    
     62
    6363    /**
    64      * The factor to multiply with an easting coordinate to get from "easting 
     64     * The factor to multiply with an easting coordinate to get from "easting
    6565     * units per pixel" to "meters per pixel"
    6666     */
Note: See TracChangeset for help on using the changeset viewer.