Skip to content

Commit cffe10e

Browse files
rouaultrouault
authored
Merge pull request OSGeo#3741 from rouault/mercator_spherical
Implement EPSG:1026 Mercator (Spherical) method
2 parents 4c3b426 + 2ba6637 commit cffe10e

File tree

11 files changed

+214
-11
lines changed

11 files changed

+214
-11
lines changed

‎docs/source/usage/ellipsoids.rst

Lines changed: 11 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -65,23 +65,29 @@ the ellipsoid into a sphere with features defined by the ellipsoid.
6565
Ellipsoid spherification parameters
6666
+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
6767

68-
.. option:: +R_A=<value>
68+
.. option:: +R_A
6969

7070
A sphere with the same surface area as the ellipsoid.
7171

72-
.. option:: +R_V=<value>
72+
.. option:: +R_V
7373

7474
A sphere with the same volume as the ellipsoid.
7575

76-
.. option:: +R_a=<value>
76+
.. option:: +R_C
77+
78+
.. versionadded:: 9.3.0
79+
80+
A sphere whose radius is the radius of the conformal sphere at :math:`\phi_0`.
81+
82+
.. option:: +R_a
7783

7884
A sphere with :math:`R = (a + b)/2` (arithmetic mean).
7985

80-
.. option:: +R_g=<value>
86+
.. option:: +R_g
8187

8288
A sphere with :math:`R = \sqrt{ab}` (geometric mean).
8389

84-
.. option:: +R_h=<value>
90+
.. option:: +R_h
8591

8692
A sphere with :math:`R = 2ab/(a+b)` (harmonic mean).
8793

‎include/proj/coordinateoperation.hpp

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1292,6 +1292,11 @@ class PROJ_GCC_DLL Conversion : public SingleOperation {
12921292
const common::Angle &centerLong, const common::Length &falseEasting,
12931293
const common::Length &falseNorthing);
12941294

1295+
PROJ_DLL static ConversionNNPtr createMercatorSpherical(
1296+
const util::PropertyMap &properties, const common::Angle &centerLat,
1297+
const common::Angle &centerLong, const common::Length &falseEasting,
1298+
const common::Length &falseNorthing);
1299+
12951300
PROJ_DLL static ConversionNNPtr
12961301
createMollweide(const util::PropertyMap &properties,
12971302
const common::Angle &centerLong,

‎scripts/reference_exported_symbols.txt

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -589,6 +589,7 @@ osgeo::proj::operation::Conversion::createLambertConicConformal_2SP_Michigan(osg
589589
osgeo::proj::operation::Conversion::createLambertConicConformal_2SP(osgeo::proj::util::PropertyMap const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Length const&, osgeo::proj::common::Length const&)
590590
osgeo::proj::operation::Conversion::createLambertCylindricalEqualArea(osgeo::proj::util::PropertyMap const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Length const&, osgeo::proj::common::Length const&)
591591
osgeo::proj::operation::Conversion::createLambertCylindricalEqualAreaSpherical(osgeo::proj::util::PropertyMap const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Length const&, osgeo::proj::common::Length const&)
592+
osgeo::proj::operation::Conversion::createMercatorSpherical(osgeo::proj::util::PropertyMap const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Length const&, osgeo::proj::common::Length const&)
592593
osgeo::proj::operation::Conversion::createMercatorVariantA(osgeo::proj::util::PropertyMap const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Scale const&, osgeo::proj::common::Length const&, osgeo::proj::common::Length const&)
593594
osgeo::proj::operation::Conversion::createMercatorVariantB(osgeo::proj::util::PropertyMap const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Length const&, osgeo::proj::common::Length const&)
594595
osgeo::proj::operation::Conversion::createMillerCylindrical(osgeo::proj::util::PropertyMap const&, osgeo::proj::common::Angle const&, osgeo::proj::common::Length const&, osgeo::proj::common::Length const&)

‎src/ell_set.cpp

Lines changed: 27 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -49,16 +49,21 @@ int pj_ellipsoid(PJ *P) {
4949
5050
Spherification parameters supported are:
5151
R_A, which gives a sphere with the same surface area as the
52-
ellipsoid R_V, which gives a sphere with the same volume as the ellipsoid
52+
ellipsoid R_V, which gives a sphere with the same volume as the
53+
ellipsoid
5354
5455
R_a, which gives a sphere with R = (a + b)/2 (arithmetic mean)
5556
R_g, which gives a sphere with R = sqrt(a*b) (geometric mean)
5657
R_h, which gives a sphere with R = 2*a*b/(a+b) (harmonic mean)
5758
5859
R_lat_a=phi, which gives a sphere with R being the arithmetic mean
59-
of of the corresponding ellipsoid at latitude phi. R_lat_g=phi, which gives
60-
a sphere with R being the geometric mean of of the corresponding ellipsoid
61-
at latitude phi.
60+
of of the corresponding ellipsoid at latitude phi.
61+
R_lat_g=phi, which gives
62+
a sphere with R being the geometric mean of of the corresponding
63+
ellipsoid at latitude phi.
64+
65+
R_C, which gives a sphere with the radius of the conformal sphere
66+
at phi0.
6267
6368
If R is given as size parameter, any shape and spherification parameters
6469
given are ignored.
@@ -349,8 +354,8 @@ static const double RV6 = 55 / 1296.;
349354
/***************************************************************************************/
350355
static int ellps_spherification(PJ *P) {
351356
/***************************************************************************************/
352-
const char *keys[] = {"R_A", "R_V", "R_a", "R_g",
353-
"R_h", "R_lat_a", "R_lat_g"};
357+
const char *keys[] = {"R_A", "R_V", "R_a", "R_g",
358+
"R_h", "R_lat_a", "R_lat_g", "R_C"};
354359
size_t len, i;
355360
paralist *par = nullptr;
356361

@@ -428,6 +433,22 @@ static int ellps_spherification(PJ *P) {
428433
else /* geometric */
429434
P->a *= sqrt(1 - P->es) / t;
430435
break;
436+
437+
/* R_C - a sphere with R = radius of conformal sphere, taken at a
438+
* latitude that is phi0 (note: at least for mercator. for other
439+
* projection methods, this could be phi1)
440+
* Formula from IOGP Publication 373-7-2 – Geomatics Guidance Note number 7,
441+
* part 2 1.1 Ellipsoid parameters
442+
*/
443+
case 7:
444+
t = sin(P->phi0);
445+
t = 1 - P->es * t * t;
446+
if (t == 0.) {
447+
proj_log_error(P, _("Invalid eccentricity"));
448+
return proj_errno_set(P, PROJ_ERR_INVALID_OP_ILLEGAL_ARG_VALUE);
449+
}
450+
P->a *= sqrt(1 - P->es) / t;
451+
break;
431452
}
432453

433454
if (P->a <= 0.) {

‎src/iso19111/io.cpp

Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -11180,7 +11180,26 @@ PROJStringParser::Private::buildProjectedCRS(int iStep,
1118011180
}
1118111181
}
1118211182
} else if (hasParamValue(step, "lat_ts")) {
11183+
if (hasParamValue(step, "R_C") &&
11184+
!geodCRS->ellipsoid()->isSphere() &&
11185+
getAngularValue(getParamValue(step, "lat_ts")) != 0) {
11186+
throw ParsingException("lat_ts != 0 not supported for "
11187+
"spherical Mercator on an ellipsoid");
11188+
}
1118311189
mapping = getMapping(EPSG_CODE_METHOD_MERCATOR_VARIANT_B);
11190+
} else if (hasParamValue(step, "R_C")) {
11191+
const auto &k = getParamValueK(step);
11192+
if (!k.empty() && getNumericValue(k) != 1.0) {
11193+
if (geodCRS->ellipsoid()->isSphere()) {
11194+
mapping = getMapping(EPSG_CODE_METHOD_MERCATOR_VARIANT_A);
11195+
} else {
11196+
throw ParsingException(
11197+
"k_0 != 1 not supported for spherical Mercator on an "
11198+
"ellipsoid");
11199+
}
11200+
} else {
11201+
mapping = getMapping(EPSG_CODE_METHOD_MERCATOR_SPHERICAL);
11202+
}
1118411203
} else {
1118511204
mapping = getMapping(EPSG_CODE_METHOD_MERCATOR_VARIANT_A);
1118611205
}

‎src/iso19111/operation/conversion.cpp

Lines changed: 32 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1751,6 +1751,38 @@ ConversionNNPtr Conversion::createPopularVisualisationPseudoMercator(
17511751

17521752
// ---------------------------------------------------------------------------
17531753

1754+
/** \brief Instantiate a conversion based on the
1755+
* <a href="../../../operations/projections/merc.html">
1756+
* Mercator</a> projection method, using its spherical formulation
1757+
*
1758+
* When used with an ellipsoid, the radius used is the radius of the conformal
1759+
* sphere at centerLat.
1760+
*
1761+
* This method is defined as
1762+
* <a
1763+
* href="https://epsg.org/coord-operation-method_1026/Mercator-Spherical.html">
1764+
* EPSG:1026</a>.
1765+
*
1766+
* @param properties See \ref general_properties of the conversion. If the name
1767+
* is not provided, it is automatically set.
1768+
* @param centerLat See \ref center_latitude . Usually 0
1769+
* @param centerLong See \ref center_longitude . Usually 0
1770+
* @param falseEasting See \ref false_easting . Usually 0
1771+
* @param falseNorthing See \ref false_northing . Usually 0
1772+
* @return a new Conversion.
1773+
* @since 9.3
1774+
*/
1775+
ConversionNNPtr Conversion::createMercatorSpherical(
1776+
const util::PropertyMap &properties, const common::Angle &centerLat,
1777+
const common::Angle &centerLong, const common::Length &falseEasting,
1778+
const common::Length &falseNorthing) {
1779+
return create(
1780+
properties, EPSG_CODE_METHOD_MERCATOR_SPHERICAL,
1781+
createParams(centerLat, centerLong, falseEasting, falseNorthing));
1782+
}
1783+
1784+
// ---------------------------------------------------------------------------
1785+
17541786
/** \brief Instantiate a conversion based on the
17551787
* <a href="../../../operations/projections/moll.html">
17561788
* Mollweide</a> projection method.

‎src/iso19111/operation/parammappings.cpp

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -770,6 +770,9 @@ static const MethodMapping projectionMethodMappings[] = {
770770
// handled manually
771771
"webmerc", nullptr, paramsNatOrigin},
772772

773+
{EPSG_NAME_METHOD_MERCATOR_SPHERICAL, EPSG_CODE_METHOD_MERCATOR_SPHERICAL,
774+
nullptr, "merc", "R_C", paramsNatOrigin},
775+
773776
{PROJ_WKT2_NAME_METHOD_MOLLWEIDE, 0, "Mollweide", "moll", nullptr,
774777
paramsLongNatOrigin},
775778

@@ -938,6 +941,7 @@ const struct MethodNameCode methodNameCodes[] = {
938941
METHOD_NAME_CODE(KROVAK),
939942
METHOD_NAME_CODE(LAMBERT_AZIMUTHAL_EQUAL_AREA),
940943
METHOD_NAME_CODE(POPULAR_VISUALISATION_PSEUDO_MERCATOR),
944+
METHOD_NAME_CODE(MERCATOR_SPHERICAL),
941945
METHOD_NAME_CODE(MERCATOR_VARIANT_A),
942946
METHOD_NAME_CODE(MERCATOR_VARIANT_B),
943947
METHOD_NAME_CODE(OBLIQUE_STEREOGRAPHIC),

‎src/proj_constants.h

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -187,6 +187,9 @@
187187
"Popular Visualisation Pseudo Mercator"
188188
#define EPSG_CODE_METHOD_POPULAR_VISUALISATION_PSEUDO_MERCATOR 1024
189189

190+
#define EPSG_NAME_METHOD_MERCATOR_SPHERICAL "Mercator (Spherical)"
191+
#define EPSG_CODE_METHOD_MERCATOR_SPHERICAL 1026
192+
190193
#define PROJ_WKT2_NAME_METHOD_MOLLWEIDE "Mollweide"
191194

192195
#define PROJ_WKT2_NAME_METHOD_NATURAL_EARTH "Natural Earth"

‎test/gie/builtins.gie

Lines changed: 16 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -3749,6 +3749,22 @@ direction inverse
37493749
accept 0 1e-15
37503750
expect 0 57.295779513e-15
37513751

3752+
-------------------------------------------------------------------------------
3753+
# Mercator Spherical
3754+
operation +proj=merc +R_C +ellps=WGS84 +lat_0=45
3755+
-------------------------------------------------------------------------------
3756+
accept 2 49
3757+
expect 221892.515234695253 6253822.971804938279
3758+
3759+
-------------------------------------------------------------------------------
3760+
# Mercator Spherical
3761+
operation +proj=merc +R_C +a=6371007 +b=6371007
3762+
-------------------------------------------------------------------------------
3763+
# Test point from IOGP guidance note 7.2
3764+
tolerance 1cm
3765+
accept -100.33333333333333 24.381786944444446
3766+
expect -11156569.90 2796869.94
3767+
37523768

37533769
===============================================================================
37543770
# Miller Oblated Stereographic

‎test/unit/test_io.cpp

Lines changed: 17 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10618,6 +10618,23 @@ TEST(io, projparse_cea_ellipsoidal_with_k_0) {
1061810618

1061910619
// ---------------------------------------------------------------------------
1062010620

10621+
TEST(io, projparse_merc_spherical_on_ellipsoid) {
10622+
std::string input("+proj=merc +R_C +lat_0=1 +lon_0=2 +x_0=3 +y_0=4 "
10623+
"+ellps=WGS84 +units=m +no_defs +type=crs");
10624+
auto obj = PROJStringParser().createFromPROJString(input);
10625+
auto crs = nn_dynamic_pointer_cast<ProjectedCRS>(obj);
10626+
ASSERT_TRUE(crs != nullptr);
10627+
EXPECT_EQ(crs->derivingConversion()->method()->getEPSGCode(),
10628+
EPSG_CODE_METHOD_MERCATOR_SPHERICAL);
10629+
EXPECT_EQ(
10630+
crs->exportToPROJString(
10631+
PROJStringFormatter::create(PROJStringFormatter::Convention::PROJ_4)
10632+
.get()),
10633+
input);
10634+
}
10635+
10636+
// ---------------------------------------------------------------------------
10637+
1062110638
TEST(io, projparse_geos_sweep_x) {
1062210639
auto obj = PROJStringParser().createFromPROJString(
1062310640
"+proj=geos +sweep=x +h=1 +type=crs");

0 commit comments

Comments
 (0)