001// License: GPL. For details, see LICENSE file.
002package org.openstreetmap.josm.data.projection;
003
004import java.util.Collections;
005import java.util.HashMap;
006import java.util.Map;
007import java.util.function.DoubleUnaryOperator;
008
009import org.openstreetmap.josm.data.Bounds;
010import org.openstreetmap.josm.data.ProjectionBounds;
011import org.openstreetmap.josm.data.coor.EastNorth;
012import org.openstreetmap.josm.data.coor.LatLon;
013import org.openstreetmap.josm.data.projection.datum.Datum;
014import org.openstreetmap.josm.data.projection.proj.Proj;
015import org.openstreetmap.josm.tools.Utils;
016
017/**
018 * Implementation of the Projection interface that represents a coordinate reference system and delegates
019 * the real projection and datum conversion to other classes.
020 *
021 * It handles false easting and northing, central meridian and general scale factor before calling the
022 * delegate projection.
023 *
024 * Forwards lat/lon values to the real projection in units of radians.
025 *
026 * The fields are named after Proj.4 parameters.
027 *
028 * Subclasses of AbstractProjection must set ellps and proj to a non-null value.
029 * In addition, either datum or nadgrid has to be initialized to some value.
030 */
031public abstract class AbstractProjection implements Projection {
032
033    protected Ellipsoid ellps;
034    protected Datum datum;
035    protected Proj proj;
036    protected double x0;            /* false easting (in meters) */
037    protected double y0;            /* false northing (in meters) */
038    protected double lon0;          /* central meridian */
039    protected double pm;            /* prime meridian */
040    protected double k0 = 1.0;      /* general scale factor */
041    protected double toMeter = 1.0; /* switch from meters to east/north coordinate units */
042
043    private volatile ProjectionBounds projectionBoundsBox;
044
045    /**
046     * Get the base ellipsoid that this projection uses.
047     * @return The {@link Ellipsoid}
048     */
049    public final Ellipsoid getEllipsoid() {
050        return ellps;
051    }
052
053    /**
054     * Gets the datum this projection is based on.
055     * @return The datum
056     */
057    public final Datum getDatum() {
058        return datum;
059    }
060
061    /**
062     * Replies the projection (in the narrow sense)
063     * @return The projection object
064     */
065    public final Proj getProj() {
066        return proj;
067    }
068
069    /**
070     * Gets an east offset that gets applied when converting the coordinate
071     * @return The offset to apply in meter
072     */
073    public final double getFalseEasting() {
074        return x0;
075    }
076
077    /**
078     * Gets an north offset that gets applied when converting the coordinate
079     * @return The offset to apply in meter
080     */
081    public final double getFalseNorthing() {
082        return y0;
083    }
084
085    /**
086     * Gets the meridian that this projection is centered on.
087     * @return The longitude of the meridian.
088     */
089    public final double getCentralMeridian() {
090        return lon0;
091    }
092
093    public final double getScaleFactor() {
094        return k0;
095    }
096
097    /**
098     * Get the factor that converts meters to intended units of east/north coordinates.
099     *
100     * For projected coordinate systems, the semi-major axis of the ellipsoid is
101     * always given in meters, which means the preliminary projection result will
102     * be in meters as well. This factor is used to convert to the intended units
103     * of east/north coordinates (e.g. feet in the US).
104     *
105     * For geographic coordinate systems, the preliminary "projection" result will
106     * be in degrees, so there is no reason to convert anything and this factor
107     * will by 1 by default.
108     *
109     * @return factor that converts meters to intended units of east/north coordinates
110     */
111    public final double getToMeter() {
112        return toMeter;
113    }
114
115    @Override
116    public EastNorth latlon2eastNorth(LatLon ll) {
117        ll = datum.fromWGS84(ll);
118        double[] en = proj.project(Math.toRadians(ll.lat()), Math.toRadians(LatLon.normalizeLon(ll.lon() - lon0 - pm)));
119        return new EastNorth((ellps.a * k0 * en[0] + x0) / toMeter, (ellps.a * k0 * en[1] + y0) / toMeter);
120    }
121
122    @Override
123    public LatLon eastNorth2latlon(EastNorth en) {
124        return eastNorth2latlon(en, LatLon::normalizeLon);
125    }
126
127    @Override
128    public LatLon eastNorth2latlonClamped(EastNorth en) {
129        LatLon ll = eastNorth2latlon(en, lon -> Utils.clamp(lon, -180, 180));
130        Bounds bounds = getWorldBoundsLatLon();
131        return new LatLon(Utils.clamp(ll.lat(), bounds.getMinLat(), bounds.getMaxLat()),
132                Utils.clamp(ll.lon(), bounds.getMinLon(), bounds.getMaxLon()));
133    }
134
135    private LatLon eastNorth2latlon(EastNorth en, DoubleUnaryOperator normalizeLon) {
136        double[] latlonRad = proj.invproject((en.east() * toMeter - x0) / ellps.a / k0, (en.north() * toMeter - y0) / ellps.a / k0);
137        double lon = Math.toDegrees(latlonRad[1]) + lon0 + pm;
138        LatLon ll = new LatLon(Math.toDegrees(latlonRad[0]), normalizeLon.applyAsDouble(lon));
139        return datum.toWGS84(ll);
140    }
141
142    @Override
143    public Map<ProjectionBounds, Projecting> getProjectingsForArea(ProjectionBounds area) {
144        if (proj.lonIsLinearToEast()) {
145            //FIXME: Respect datum?
146            // wrap the wrold around
147            Bounds bounds = getWorldBoundsLatLon();
148            double minEast = latlon2eastNorth(bounds.getMin()).east();
149            double maxEast = latlon2eastNorth(bounds.getMax()).east();
150            double dEast = maxEast - minEast;
151            if ((area.minEast < minEast || area.maxEast > maxEast) && dEast > 0) {
152                // We could handle the dEast < 0 case but we don't need it atm.
153                int minChunk = (int) Math.floor((area.minEast - minEast) / dEast);
154                int maxChunk = (int) Math.floor((area.maxEast - minEast) / dEast);
155                HashMap<ProjectionBounds, Projecting> ret = new HashMap<>();
156                for (int chunk = minChunk; chunk <= maxChunk; chunk++) {
157                    ret.put(new ProjectionBounds(Math.max(area.minEast, minEast + chunk * dEast), area.minNorth,
158                            Math.min(area.maxEast, maxEast + chunk * dEast), area.maxNorth),
159                            new ShiftedProjecting(this, new EastNorth(-chunk * dEast, 0)));
160                }
161                return ret;
162            }
163        }
164
165        return Collections.singletonMap(area, this);
166    }
167
168    @Override
169    public double getDefaultZoomInPPD() {
170        // this will set the map scaler to about 1000 m
171        return 10;
172    }
173
174    /**
175     * @return The EPSG Code of this CRS, null if it doesn't have one.
176     */
177    public abstract Integer getEpsgCode();
178
179    /**
180     * Default implementation of toCode().
181     * Should be overridden, if there is no EPSG code for this CRS.
182     */
183    @Override
184    public String toCode() {
185        return "EPSG:" + getEpsgCode();
186    }
187
188    protected static final double convertMinuteSecond(double minute, double second) {
189        return (minute/60.0) + (second/3600.0);
190    }
191
192    protected static final double convertDegreeMinuteSecond(double degree, double minute, double second) {
193        return degree + (minute/60.0) + (second/3600.0);
194    }
195
196    @Override
197    public final ProjectionBounds getWorldBoundsBoxEastNorth() {
198        ProjectionBounds result = projectionBoundsBox;
199        if (result == null) {
200            synchronized (this) {
201                result = projectionBoundsBox;
202                if (result == null) {
203                    Bounds b = getWorldBoundsLatLon();
204                    // add 4 corners
205                    result = new ProjectionBounds(latlon2eastNorth(b.getMin()));
206                    result.extend(latlon2eastNorth(b.getMax()));
207                    result.extend(latlon2eastNorth(new LatLon(b.getMinLat(), b.getMaxLon())));
208                    result.extend(latlon2eastNorth(new LatLon(b.getMaxLat(), b.getMinLon())));
209                    // and trace along the outline
210                    double dLon = (b.getMaxLon() - b.getMinLon()) / 1000;
211                    double dLat = (b.getMaxLat() - b.getMinLat()) / 1000;
212                    for (double lon = b.getMinLon(); lon < b.getMaxLon(); lon += dLon) {
213                        result.extend(latlon2eastNorth(new LatLon(b.getMinLat(), lon)));
214                        result.extend(latlon2eastNorth(new LatLon(b.getMaxLat(), lon)));
215                    }
216                    for (double lat = b.getMinLat(); lat < b.getMaxLat(); lat += dLat) {
217                        result.extend(latlon2eastNorth(new LatLon(lat, b.getMinLon())));
218                        result.extend(latlon2eastNorth(new LatLon(lat, b.getMaxLon())));
219                    }
220                    projectionBoundsBox = result;
221                }
222            }
223        }
224        return projectionBoundsBox;
225    }
226
227    @Override
228    public Projection getBaseProjection() {
229        return this;
230    }
231}