Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: OpenGIS Simple Features Reference Implementation
4 : * Purpose: OGRSpatialReference translation to/from ESRI .prj definitions.
5 : * Author: Frank Warmerdam, warmerdam@pobox.com
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2000, Frank Warmerdam
9 : * Copyright (c) 2007-2013, Even Rouault <even dot rouault at spatialys.com>
10 : * Copyright (c) 2013, Kyle Shannon <kyle at pobox dot com>
11 : *
12 : * Permission is hereby granted, free of charge, to any person obtaining a
13 : * copy of this software and associated documentation files (the "Software"),
14 : * to deal in the Software without restriction, including without limitation
15 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
16 : * and/or sell copies of the Software, and to permit persons to whom the
17 : * Software is furnished to do so, subject to the following conditions:
18 : *
19 : * The above copyright notice and this permission notice shall be included
20 : * in all copies or substantial portions of the Software.
21 : *
22 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
23 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
24 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
25 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
26 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
27 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
28 : * DEALINGS IN THE SOFTWARE.
29 : ****************************************************************************/
30 :
31 : #include "cpl_port.h"
32 : #include "ogr_spatialref.h"
33 : #include "ogr_srs_esri_names.h"
34 :
35 : #include <cmath>
36 : #include <climits>
37 : #include <cstddef>
38 : #include <cstdio>
39 : #include <cstdlib>
40 : #include <cstring>
41 : #include <algorithm>
42 : #include <limits>
43 : #include <string>
44 :
45 : #include "cpl_conv.h"
46 : #include "cpl_csv.h"
47 : #include "cpl_error.h"
48 : #include "cpl_multiproc.h"
49 : #include "cpl_string.h"
50 : #include "cpl_vsi.h"
51 : #include "ogr_core.h"
52 : #include "ogr_p.h"
53 : #include "ogr_srs_api.h"
54 :
55 : /* -------------------------------------------------------------------- */
56 : /* Table relating USGS and ESRI state plane zones. */
57 : /* -------------------------------------------------------------------- */
58 : constexpr int anUsgsEsriZones[] = {
59 : 101, 3101, 102, 3126, 201, 3151, 202, 3176, 203, 3201, 301, 3226,
60 : 302, 3251, 401, 3276, 402, 3301, 403, 3326, 404, 3351, 405, 3376,
61 : 406, 3401, 407, 3426, 501, 3451, 502, 3476, 503, 3501, 600, 3526,
62 : 700, 3551, 901, 3601, 902, 3626, 903, 3576, 1001, 3651, 1002, 3676,
63 : 1101, 3701, 1102, 3726, 1103, 3751, 1201, 3776, 1202, 3801, 1301, 3826,
64 : 1302, 3851, 1401, 3876, 1402, 3901, 1501, 3926, 1502, 3951, 1601, 3976,
65 : 1602, 4001, 1701, 4026, 1702, 4051, 1703, 6426, 1801, 4076, 1802, 4101,
66 : 1900, 4126, 2001, 4151, 2002, 4176, 2101, 4201, 2102, 4226, 2103, 4251,
67 : 2111, 6351, 2112, 6376, 2113, 6401, 2201, 4276, 2202, 4301, 2203, 4326,
68 : 2301, 4351, 2302, 4376, 2401, 4401, 2402, 4426, 2403, 4451, 2500, 0,
69 : 2501, 4476, 2502, 4501, 2503, 4526, 2600, 0, 2601, 4551, 2602, 4576,
70 : 2701, 4601, 2702, 4626, 2703, 4651, 2800, 4676, 2900, 4701, 3001, 4726,
71 : 3002, 4751, 3003, 4776, 3101, 4801, 3102, 4826, 3103, 4851, 3104, 4876,
72 : 3200, 4901, 3301, 4926, 3302, 4951, 3401, 4976, 3402, 5001, 3501, 5026,
73 : 3502, 5051, 3601, 5076, 3602, 5101, 3701, 5126, 3702, 5151, 3800, 5176,
74 : 3900, 0, 3901, 5201, 3902, 5226, 4001, 5251, 4002, 5276, 4100, 5301,
75 : 4201, 5326, 4202, 5351, 4203, 5376, 4204, 5401, 4205, 5426, 4301, 5451,
76 : 4302, 5476, 4303, 5501, 4400, 5526, 4501, 5551, 4502, 5576, 4601, 5601,
77 : 4602, 5626, 4701, 5651, 4702, 5676, 4801, 5701, 4802, 5726, 4803, 5751,
78 : 4901, 5776, 4902, 5801, 4903, 5826, 4904, 5851, 5001, 6101, 5002, 6126,
79 : 5003, 6151, 5004, 6176, 5005, 6201, 5006, 6226, 5007, 6251, 5008, 6276,
80 : 5009, 6301, 5010, 6326, 5101, 5876, 5102, 5901, 5103, 5926, 5104, 5951,
81 : 5105, 5976, 5201, 6001, 5200, 6026, 5200, 6076, 5201, 6051, 5202, 6051,
82 : 5300, 0, 5400, 0};
83 :
84 : /************************************************************************/
85 : /* ESRIToUSGSZone() */
86 : /* */
87 : /* Convert ESRI style state plane zones to USGS style state */
88 : /* plane zones. */
89 : /************************************************************************/
90 :
91 0 : static int ESRIToUSGSZone(int nESRIZone)
92 :
93 : {
94 : // anUsgsEsriZones is a series of ints where 2 consecutive integers
95 : // are used to map from USGS to ESRI state plane zones.
96 : // TODO(schwehr): Would be better as a std::map.
97 0 : const int nPairs = sizeof(anUsgsEsriZones) / (2 * sizeof(int));
98 :
99 0 : for (int i = 0; i < nPairs; i++)
100 : {
101 0 : if (anUsgsEsriZones[i * 2 + 1] == nESRIZone)
102 0 : return anUsgsEsriZones[i * 2];
103 : }
104 :
105 0 : return 0;
106 : }
107 :
108 : /************************************************************************/
109 : /* OSRImportFromESRI() */
110 : /************************************************************************/
111 :
112 : /**
113 : * \brief Import coordinate system from ESRI .prj format(s).
114 : *
115 : * This function is the same as the C++ method
116 : * OGRSpatialReference::importFromESRI().
117 : */
118 8 : OGRErr OSRImportFromESRI(OGRSpatialReferenceH hSRS, char **papszPrj)
119 :
120 : {
121 8 : VALIDATE_POINTER1(hSRS, "OSRImportFromESRI", OGRERR_FAILURE);
122 :
123 8 : return reinterpret_cast<OGRSpatialReference *>(hSRS)->importFromESRI(
124 8 : papszPrj);
125 : }
126 :
127 : /************************************************************************/
128 : /* OSR_GDV() */
129 : /* */
130 : /* Fetch a particular parameter out of the parameter list, or */
131 : /* the indicated default if it isn't available. This is a */
132 : /* helper function for importFromESRI(). */
133 : /************************************************************************/
134 :
135 120 : static double OSR_GDV(char **papszNV, const char *pszField,
136 : double dfDefaultValue)
137 :
138 : {
139 120 : if (papszNV == nullptr || papszNV[0] == nullptr)
140 0 : return dfDefaultValue;
141 :
142 120 : if (STARTS_WITH_CI(pszField, "PARAM_"))
143 : {
144 88 : int iLine = 0; // Used after for loop.
145 698 : for (; papszNV[iLine] != nullptr &&
146 698 : !STARTS_WITH_CI(papszNV[iLine], "Paramet");
147 : iLine++)
148 : {
149 : }
150 :
151 387 : for (int nOffset = atoi(pszField + 6);
152 387 : papszNV[iLine] != nullptr && nOffset > 0; iLine++)
153 : {
154 299 : if (strlen(papszNV[iLine]) > 0)
155 299 : nOffset--;
156 : }
157 :
158 88 : while (papszNV[iLine] != nullptr && strlen(papszNV[iLine]) == 0)
159 0 : iLine++;
160 :
161 88 : if (papszNV[iLine] != nullptr)
162 : {
163 88 : char *const pszLine = papszNV[iLine];
164 :
165 : // Trim comments.
166 6453 : for (int i = 0; pszLine[i] != '\0'; i++)
167 : {
168 6365 : if (pszLine[i] == '/' && pszLine[i + 1] == '*')
169 87 : pszLine[i] = '\0';
170 : }
171 :
172 88 : double dfValue = 0.0;
173 88 : char **papszTokens = CSLTokenizeString(papszNV[iLine]);
174 88 : if (CSLCount(papszTokens) == 3)
175 : {
176 : // http://agdcftp1.wr.usgs.gov/pub/projects/lcc/akcan_lcc/akcan.tar.gz
177 : // contains weird values for the second. Ignore it and
178 : // the result looks correct.
179 56 : double dfSecond = CPLAtof(papszTokens[2]);
180 56 : if (dfSecond < 0.0 || dfSecond >= 60.0)
181 0 : dfSecond = 0.0;
182 :
183 56 : dfValue = std::abs(CPLAtof(papszTokens[0])) +
184 56 : CPLAtof(papszTokens[1]) / 60.0 + dfSecond / 3600.0;
185 :
186 56 : if (CPLAtof(papszTokens[0]) < 0.0)
187 15 : dfValue *= -1;
188 : }
189 32 : else if (CSLCount(papszTokens) > 0)
190 : {
191 32 : dfValue = CPLAtof(papszTokens[0]);
192 : }
193 : else
194 : {
195 0 : dfValue = dfDefaultValue;
196 : }
197 :
198 88 : CSLDestroy(papszTokens);
199 :
200 88 : return dfValue;
201 : }
202 :
203 0 : return dfDefaultValue;
204 : }
205 :
206 32 : int iLine = 0; // Used after for loop.
207 156 : for (; papszNV[iLine] != nullptr &&
208 152 : !EQUALN(papszNV[iLine], pszField, strlen(pszField));
209 : iLine++)
210 : {
211 : }
212 :
213 32 : if (papszNV[iLine] == nullptr)
214 4 : return dfDefaultValue;
215 :
216 28 : return CPLAtof(papszNV[iLine] + strlen(pszField));
217 : }
218 :
219 : /************************************************************************/
220 : /* OSR_GDS() */
221 : /************************************************************************/
222 :
223 106 : static CPLString OSR_GDS(char **papszNV, const char *pszField,
224 : const char *pszDefaultValue)
225 :
226 : {
227 106 : if (papszNV == nullptr || papszNV[0] == nullptr)
228 0 : return pszDefaultValue;
229 :
230 106 : int iLine = 0; // Used after for loop.
231 309 : for (; papszNV[iLine] != nullptr &&
232 302 : !EQUALN(papszNV[iLine], pszField, strlen(pszField));
233 : iLine++)
234 : {
235 : }
236 :
237 106 : if (papszNV[iLine] == nullptr)
238 7 : return pszDefaultValue;
239 :
240 99 : char **papszTokens = CSLTokenizeString(papszNV[iLine]);
241 :
242 : CPLString osResult =
243 198 : CSLCount(papszTokens) > 1 ? papszTokens[1] : pszDefaultValue;
244 :
245 99 : CSLDestroy(papszTokens);
246 99 : return osResult;
247 : }
248 :
249 : /************************************************************************/
250 : /* importFromESRI() */
251 : /************************************************************************/
252 :
253 : /**
254 : * \brief Import coordinate system from ESRI .prj format(s).
255 : *
256 : * This function will read the text loaded from an ESRI .prj file, and
257 : * translate it into an OGRSpatialReference definition. This should support
258 : * many (but by no means all) old style (Arc/Info 7.x) .prj files, as well
259 : * as the newer pseudo-OGC WKT .prj files. Note that new style .prj files
260 : * are in OGC WKT format, but require some manipulation to correct datum
261 : * names, and units on some projection parameters. This is addressed within
262 : * importFromESRI() by an automatic call to morphFromESRI().
263 : *
264 : * Currently only GEOGRAPHIC, UTM, STATEPLANE, GREATBRITIAN_GRID, ALBERS,
265 : * EQUIDISTANT_CONIC, TRANSVERSE (mercator), POLAR, MERCATOR and POLYCONIC
266 : * projections are supported from old style files.
267 : *
268 : * At this time there is no equivalent exportToESRI() method. Writing old
269 : * style .prj files is not supported by OGRSpatialReference. However the
270 : * morphToESRI() and exportToWkt() methods can be used to generate output
271 : * suitable to write to new style (Arc 8) .prj files.
272 : *
273 : * This function is the equivalent of the C function OSRImportFromESRI().
274 : *
275 : * @param papszPrj NULL terminated list of strings containing the definition.
276 : *
277 : * @return OGRERR_NONE on success or an error code in case of failure.
278 : */
279 :
280 910 : OGRErr OGRSpatialReference::importFromESRI(char **papszPrj)
281 :
282 : {
283 910 : if (papszPrj == nullptr || papszPrj[0] == nullptr)
284 29 : return OGRERR_CORRUPT_DATA;
285 :
286 : /* -------------------------------------------------------------------- */
287 : /* ArcGIS and related products now use a variant of Well Known */
288 : /* Text. Try to recognize this and ingest it. WKT is usually */
289 : /* all on one line, but we will accept multi-line formats and */
290 : /* concatenate. */
291 : /* -------------------------------------------------------------------- */
292 881 : if (STARTS_WITH_CI(papszPrj[0], "GEOGCS") ||
293 536 : STARTS_WITH_CI(papszPrj[0], "PROJCS") ||
294 40 : STARTS_WITH_CI(papszPrj[0], "LOCAL_CS")
295 : // Also accept COMPD_CS, even if it is unclear that it is valid
296 : // traditional ESRI WKT. But people might use such PRJ file
297 : // See https://github.com/OSGeo/gdal/issues/1881
298 40 : || STARTS_WITH_CI(papszPrj[0], "COMPD_CS"))
299 : {
300 1682 : std::string osWKT(papszPrj[0]);
301 843 : for (int i = 1; papszPrj[i] != nullptr; i++)
302 : {
303 2 : osWKT += papszPrj[i];
304 : }
305 841 : return importFromWkt(osWKT.c_str());
306 : }
307 :
308 : /* -------------------------------------------------------------------- */
309 : /* Operate on the basis of the projection name. */
310 : /* -------------------------------------------------------------------- */
311 80 : CPLString osProj = OSR_GDS(papszPrj, "Projection", "");
312 40 : bool bDatumApplied = false;
313 :
314 40 : if (EQUAL(osProj, ""))
315 : {
316 5 : CPLDebug("OGR_ESRI", "Can't find Projection");
317 5 : return OGRERR_CORRUPT_DATA;
318 : }
319 35 : else if (EQUAL(osProj, "GEOGRAPHIC"))
320 : {
321 : // Nothing to do.
322 : }
323 32 : else if (EQUAL(osProj, "utm"))
324 : {
325 12 : const double dfOsrGdv = OSR_GDV(papszPrj, "zone", 0.0);
326 12 : if (dfOsrGdv > 0 && dfOsrGdv < 61)
327 : {
328 12 : const double dfYShift = OSR_GDV(papszPrj, "Yshift", 0.0);
329 :
330 12 : SetUTM(static_cast<int>(dfOsrGdv), dfYShift == 0.0);
331 : }
332 : else
333 : {
334 0 : const double dfCentralMeridian = OSR_GDV(papszPrj, "PARAM_1", 0.0);
335 0 : const double dfRefLat = OSR_GDV(papszPrj, "PARAM_2", 0.0);
336 0 : if (dfCentralMeridian >= -180.0 && dfCentralMeridian <= 180.0)
337 : {
338 0 : const int nZone = static_cast<int>(
339 0 : (dfCentralMeridian + 183.0) / 6.0 + 0.0000001);
340 0 : SetUTM(nZone, dfRefLat >= 0.0);
341 : }
342 : }
343 : }
344 20 : else if (EQUAL(osProj, "STATEPLANE"))
345 : {
346 4 : const double dfZone = OSR_GDV(papszPrj, "zone", 0.0);
347 :
348 4 : if (dfZone < std::numeric_limits<int>::min() ||
349 4 : dfZone > std::numeric_limits<int>::max() || CPLIsNan(dfZone))
350 : {
351 0 : CPLError(CE_Failure, CPLE_AppDefined, "zone out of range: %f",
352 : dfZone);
353 0 : return OGRERR_CORRUPT_DATA;
354 : }
355 :
356 4 : int nZone = static_cast<int>(dfZone);
357 :
358 4 : if (nZone != 0)
359 0 : nZone = ESRIToUSGSZone(nZone);
360 : else
361 : {
362 4 : const double dfFipszone = OSR_GDV(papszPrj, "fipszone", 0.0);
363 :
364 4 : if (dfFipszone < std::numeric_limits<int>::min() ||
365 8 : dfFipszone > std::numeric_limits<int>::max() ||
366 4 : CPLIsNan(dfFipszone))
367 : {
368 0 : CPLError(CE_Failure, CPLE_AppDefined,
369 : "fipszone out of range: %f", dfFipszone);
370 0 : return OGRERR_CORRUPT_DATA;
371 : }
372 :
373 4 : nZone = static_cast<int>(dfFipszone);
374 : }
375 :
376 4 : if (nZone != 0)
377 : {
378 4 : bDatumApplied = true;
379 4 : if (EQUAL(OSR_GDS(papszPrj, "Datum", "NAD83"), "NAD27"))
380 0 : SetStatePlane(nZone, FALSE);
381 : else
382 4 : SetStatePlane(nZone, TRUE);
383 : }
384 : }
385 32 : else if (EQUAL(osProj, "GREATBRITIAN_GRID") ||
386 16 : EQUAL(osProj, "GREATBRITAIN_GRID"))
387 : {
388 0 : const char *pszWkt =
389 : "PROJCS[\"OSGB 1936 / British National Grid\","
390 : "GEOGCS[\"OSGB 1936\",DATUM[\"OSGB_1936\","
391 : "SPHEROID[\"Airy 1830\",6377563.396,299.3249646]],"
392 : "PRIMEM[\"Greenwich\",0],UNIT[\"degree\",0.0174532925199433]],"
393 : "PROJECTION[\"Transverse_Mercator\"],"
394 : "PARAMETER[\"latitude_of_origin\",49],"
395 : "PARAMETER[\"central_meridian\",-2],"
396 : "PARAMETER[\"scale_factor\",0.999601272],"
397 : "PARAMETER[\"false_easting\",400000],"
398 : "PARAMETER[\"false_northing\",-100000],UNIT[\"metre\",1]]";
399 :
400 0 : bDatumApplied = true;
401 0 : importFromWkt(pszWkt);
402 : }
403 16 : else if (EQUAL(osProj, "ALBERS"))
404 : {
405 13 : SetACEA(OSR_GDV(papszPrj, "PARAM_1", 0.0),
406 : OSR_GDV(papszPrj, "PARAM_2", 0.0),
407 : OSR_GDV(papszPrj, "PARAM_4", 0.0),
408 : OSR_GDV(papszPrj, "PARAM_3", 0.0),
409 : OSR_GDV(papszPrj, "PARAM_5", 0.0),
410 : OSR_GDV(papszPrj, "PARAM_6", 0.0));
411 : }
412 3 : else if (EQUAL(osProj, "LAMBERT"))
413 : {
414 0 : SetLCC(OSR_GDV(papszPrj, "PARAM_1", 0.0),
415 : OSR_GDV(papszPrj, "PARAM_2", 0.0),
416 : OSR_GDV(papszPrj, "PARAM_4", 0.0),
417 : OSR_GDV(papszPrj, "PARAM_3", 0.0),
418 : OSR_GDV(papszPrj, "PARAM_5", 0.0),
419 : OSR_GDV(papszPrj, "PARAM_6", 0.0));
420 : }
421 3 : else if (EQUAL(osProj, "LAMBERT_AZIMUTHAL"))
422 : {
423 0 : SetLAEA(OSR_GDV(papszPrj, "PARAM_2", 0.0),
424 : OSR_GDV(papszPrj, "PARAM_1", 0.0),
425 : OSR_GDV(papszPrj, "PARAM_3", 0.0),
426 : OSR_GDV(papszPrj, "PARAM_4", 0.0));
427 : }
428 3 : else if (EQUAL(osProj, "EQUIDISTANT_CONIC"))
429 : {
430 1 : const double dfStdPCount = OSR_GDV(papszPrj, "PARAM_1", 0.0);
431 : // TODO(schwehr): What is a reasonable range for StdPCount?
432 1 : if (dfStdPCount < 0 || dfStdPCount > std::numeric_limits<int>::max() ||
433 0 : CPLIsNan(dfStdPCount))
434 : {
435 1 : CPLError(CE_Failure, CPLE_AppDefined, "StdPCount out of range: %lf",
436 : dfStdPCount);
437 1 : return OGRERR_CORRUPT_DATA;
438 : }
439 0 : const int nStdPCount = static_cast<int>(dfStdPCount);
440 :
441 0 : if (nStdPCount == 1)
442 : {
443 0 : SetEC(OSR_GDV(papszPrj, "PARAM_2", 0.0),
444 : OSR_GDV(papszPrj, "PARAM_2", 0.0),
445 : OSR_GDV(papszPrj, "PARAM_4", 0.0),
446 : OSR_GDV(papszPrj, "PARAM_3", 0.0),
447 : OSR_GDV(papszPrj, "PARAM_5", 0.0),
448 : OSR_GDV(papszPrj, "PARAM_6", 0.0));
449 : }
450 : else
451 : {
452 0 : SetEC(OSR_GDV(papszPrj, "PARAM_2", 0.0),
453 : OSR_GDV(papszPrj, "PARAM_3", 0.0),
454 : OSR_GDV(papszPrj, "PARAM_5", 0.0),
455 : OSR_GDV(papszPrj, "PARAM_4", 0.0),
456 : OSR_GDV(papszPrj, "PARAM_5", 0.0),
457 : OSR_GDV(papszPrj, "PARAM_7", 0.0));
458 : }
459 : }
460 2 : else if (EQUAL(osProj, "TRANSVERSE"))
461 : {
462 1 : SetTM(OSR_GDV(papszPrj, "PARAM_3", 0.0),
463 : OSR_GDV(papszPrj, "PARAM_2", 0.0),
464 : OSR_GDV(papszPrj, "PARAM_1", 0.0),
465 : OSR_GDV(papszPrj, "PARAM_4", 0.0),
466 : OSR_GDV(papszPrj, "PARAM_5", 0.0));
467 : }
468 1 : else if (EQUAL(osProj, "POLAR"))
469 : {
470 0 : SetPS(OSR_GDV(papszPrj, "PARAM_2", 0.0),
471 : OSR_GDV(papszPrj, "PARAM_1", 0.0), 1.0,
472 : OSR_GDV(papszPrj, "PARAM_3", 0.0),
473 : OSR_GDV(papszPrj, "PARAM_4", 0.0));
474 : }
475 1 : else if (EQUAL(osProj, "MERCATOR"))
476 : {
477 1 : SetMercator2SP(OSR_GDV(papszPrj, "PARAM_2", 0.0), 0.0,
478 : OSR_GDV(papszPrj, "PARAM_1", 0.0),
479 : OSR_GDV(papszPrj, "PARAM_3", 0.0),
480 : OSR_GDV(papszPrj, "PARAM_4", 0.0));
481 : }
482 0 : else if (EQUAL(osProj, SRS_PT_MERCATOR_AUXILIARY_SPHERE))
483 : {
484 : // This is EPSG:3875 Pseudo Mercator. We might as well import it from
485 : // the EPSG spec.
486 0 : importFromEPSG(3857);
487 0 : bDatumApplied = true;
488 : }
489 0 : else if (EQUAL(osProj, "POLYCONIC"))
490 : {
491 0 : SetPolyconic(OSR_GDV(papszPrj, "PARAM_2", 0.0),
492 : OSR_GDV(papszPrj, "PARAM_1", 0.0),
493 : OSR_GDV(papszPrj, "PARAM_3", 0.0),
494 : OSR_GDV(papszPrj, "PARAM_4", 0.0));
495 : }
496 : else
497 : {
498 0 : CPLDebug("OGR_ESRI", "Unsupported projection: %s", osProj.c_str());
499 0 : SetLocalCS(osProj);
500 : }
501 :
502 : /* -------------------------------------------------------------------- */
503 : /* Try to translate the datum/spheroid. */
504 : /* -------------------------------------------------------------------- */
505 34 : if (!IsLocal() && !bDatumApplied)
506 : {
507 60 : const CPLString osDatum = OSR_GDS(papszPrj, "Datum", "");
508 :
509 55 : if (EQUAL(osDatum, "NAD27") || EQUAL(osDatum, "NAD83") ||
510 55 : EQUAL(osDatum, "WGS84") || EQUAL(osDatum, "WGS72"))
511 : {
512 20 : SetWellKnownGeogCS(osDatum);
513 : }
514 10 : else if (EQUAL(osDatum, "EUR") || EQUAL(osDatum, "ED50"))
515 : {
516 0 : SetWellKnownGeogCS("EPSG:4230");
517 : }
518 10 : else if (EQUAL(osDatum, "GDA94"))
519 : {
520 9 : SetWellKnownGeogCS("EPSG:4283");
521 : }
522 : else
523 : {
524 2 : CPLString osSpheroid = OSR_GDS(papszPrj, "Spheroid", "");
525 :
526 2 : if (EQUAL(osSpheroid, "INT1909") ||
527 1 : EQUAL(osSpheroid, "INTERNATIONAL1909"))
528 : {
529 0 : OGRSpatialReference oGCS;
530 0 : oGCS.importFromEPSG(4022);
531 0 : CopyGeogCSFrom(&oGCS);
532 : }
533 1 : else if (EQUAL(osSpheroid, "AIRY"))
534 : {
535 0 : OGRSpatialReference oGCS;
536 0 : oGCS.importFromEPSG(4001);
537 0 : CopyGeogCSFrom(&oGCS);
538 : }
539 1 : else if (EQUAL(osSpheroid, "CLARKE1866"))
540 : {
541 0 : OGRSpatialReference oGCS;
542 0 : oGCS.importFromEPSG(4008);
543 0 : CopyGeogCSFrom(&oGCS);
544 : }
545 1 : else if (EQUAL(osSpheroid, "GRS80"))
546 : {
547 0 : OGRSpatialReference oGCS;
548 0 : oGCS.importFromEPSG(4019);
549 0 : CopyGeogCSFrom(&oGCS);
550 : }
551 1 : else if (EQUAL(osSpheroid, "KRASOVSKY") ||
552 2 : EQUAL(osSpheroid, "KRASSOVSKY") ||
553 1 : EQUAL(osSpheroid, "KRASSOWSKY"))
554 : {
555 0 : OGRSpatialReference oGCS;
556 0 : oGCS.importFromEPSG(4024);
557 0 : CopyGeogCSFrom(&oGCS);
558 : }
559 1 : else if (EQUAL(osSpheroid, "Bessel"))
560 : {
561 0 : OGRSpatialReference oGCS;
562 0 : oGCS.importFromEPSG(4004);
563 0 : CopyGeogCSFrom(&oGCS);
564 : }
565 : else
566 : {
567 1 : bool bFoundParameters = false;
568 2 : for (int iLine = 0; papszPrj[iLine] != nullptr; iLine++)
569 : {
570 2 : if (STARTS_WITH_CI(papszPrj[iLine], "Parameters"))
571 : {
572 2 : char **papszTokens = CSLTokenizeString(
573 1 : papszPrj[iLine] + strlen("Parameters"));
574 1 : if (CSLCount(papszTokens) == 2)
575 : {
576 1 : OGRSpatialReference oGCS;
577 1 : const double dfSemiMajor = CPLAtof(papszTokens[0]);
578 1 : const double dfSemiMinor = CPLAtof(papszTokens[1]);
579 : const double dfInvFlattening =
580 1 : OSRCalcInvFlattening(dfSemiMajor, dfSemiMinor);
581 1 : oGCS.SetGeogCS("unknown", "unknown", "unknown",
582 : dfSemiMajor, dfInvFlattening);
583 1 : CopyGeogCSFrom(&oGCS);
584 1 : bFoundParameters = true;
585 : }
586 1 : CSLDestroy(papszTokens);
587 1 : break;
588 : }
589 : }
590 1 : if (!bFoundParameters)
591 : {
592 : // If unknown, default to WGS84 so there is something there.
593 0 : SetWellKnownGeogCS("WGS84");
594 : }
595 : }
596 : }
597 : }
598 :
599 : /* -------------------------------------------------------------------- */
600 : /* Linear units translation */
601 : /* -------------------------------------------------------------------- */
602 34 : if (IsLocal() || IsProjected())
603 : {
604 31 : const double dfOldUnits = GetLinearUnits();
605 62 : const CPLString osValue = OSR_GDS(papszPrj, "Units", "");
606 62 : CPLString osOldAuth;
607 : {
608 31 : const char *pszOldAuth = GetAuthorityCode(nullptr);
609 31 : if (pszOldAuth)
610 4 : osOldAuth = pszOldAuth;
611 : }
612 :
613 31 : if (EQUAL(osValue, ""))
614 0 : SetLinearUnitsAndUpdateParameters(SRS_UL_METER, 1.0);
615 31 : else if (EQUAL(osValue, "FEET"))
616 2 : SetLinearUnitsAndUpdateParameters(SRS_UL_US_FOOT,
617 : CPLAtof(SRS_UL_US_FOOT_CONV));
618 29 : else if (CPLAtof(osValue) != 0.0)
619 1 : SetLinearUnitsAndUpdateParameters("user-defined",
620 1 : 1.0 / CPLAtof(osValue));
621 : else
622 28 : SetLinearUnitsAndUpdateParameters(osValue, 1.0);
623 :
624 : // Reinstall authority if linear units value has not changed (bug #1697)
625 31 : const double dfNewUnits = GetLinearUnits();
626 35 : if (IsProjected() && !osOldAuth.empty() && dfOldUnits != 0.0 &&
627 4 : std::abs(dfNewUnits / dfOldUnits - 1) < 1e-8)
628 : {
629 1 : SetAuthority("PROJCS", "EPSG", atoi(osOldAuth));
630 : }
631 : }
632 :
633 34 : return OGRERR_NONE;
634 : }
635 :
636 : /************************************************************************/
637 : /* FindCodeFromDict() */
638 : /* */
639 : /* Find the code from a dict file. */
640 : /************************************************************************/
641 0 : static int FindCodeFromDict(const char *pszDictFile, const char *CSName,
642 : char *code)
643 : {
644 : /* -------------------------------------------------------------------- */
645 : /* Find and open file. */
646 : /* -------------------------------------------------------------------- */
647 0 : const char *pszFilename = CPLFindFile("gdal", pszDictFile);
648 0 : if (pszFilename == nullptr)
649 0 : return OGRERR_UNSUPPORTED_SRS;
650 :
651 0 : VSILFILE *fp = VSIFOpenL(pszFilename, "rb");
652 0 : if (fp == nullptr)
653 0 : return OGRERR_UNSUPPORTED_SRS;
654 :
655 : /* -------------------------------------------------------------------- */
656 : /* Process lines. */
657 : /* -------------------------------------------------------------------- */
658 0 : OGRErr eErr = OGRERR_UNSUPPORTED_SRS;
659 0 : const char *pszLine = nullptr;
660 :
661 0 : while ((pszLine = CPLReadLineL(fp)) != nullptr)
662 : {
663 0 : if (pszLine[0] == '#')
664 0 : continue;
665 :
666 0 : if (strstr(pszLine, CSName))
667 : {
668 0 : const char *pComma = strchr(pszLine, ',');
669 0 : if (pComma)
670 : {
671 0 : strncpy(code, pszLine, pComma - pszLine);
672 0 : code[pComma - pszLine] = '\0';
673 0 : eErr = OGRERR_NONE;
674 : }
675 0 : break;
676 : }
677 : }
678 :
679 : /* -------------------------------------------------------------------- */
680 : /* Cleanup */
681 : /* -------------------------------------------------------------------- */
682 0 : VSIFCloseL(fp);
683 :
684 0 : return eErr;
685 : }
686 :
687 : /************************************************************************/
688 : /* ImportFromESRIStatePlaneWKT() */
689 : /* */
690 : /* Search a ESRI State Plane WKT and import it. */
691 : /************************************************************************/
692 :
693 3 : OGRErr OGRSpatialReference::ImportFromESRIStatePlaneWKT(int code,
694 : const char *datumName,
695 : const char *unitsName,
696 : int pcsCode,
697 : const char *csName)
698 : {
699 : // If the CS name is known.
700 3 : if (code == 0 && !datumName && !unitsName && pcsCode == 32767 && csName)
701 : {
702 0 : char codeS[10] = {};
703 0 : if (FindCodeFromDict("esri_StatePlane_extra.wkt", csName, codeS) !=
704 : OGRERR_NONE)
705 0 : return OGRERR_FAILURE;
706 0 : return importFromDict("esri_StatePlane_extra.wkt", codeS);
707 : }
708 :
709 3 : int searchCode = -1;
710 3 : if (unitsName == nullptr)
711 0 : unitsName = "";
712 :
713 : // Find state plane prj str by pcs code only.
714 3 : if (code == 0 && !datumName && pcsCode != 32767)
715 : {
716 0 : int unitCode = 1;
717 0 : if (EQUAL(unitsName, "international_feet"))
718 0 : unitCode = 3;
719 0 : else if (strstr(unitsName, "feet") || strstr(unitsName, "foot"))
720 0 : unitCode = 2;
721 :
722 0 : for (int i = 0; statePlanePcsCodeToZoneCode[i] != 0; i += 2)
723 : {
724 0 : if (pcsCode == statePlanePcsCodeToZoneCode[i])
725 : {
726 0 : searchCode = statePlanePcsCodeToZoneCode[i + 1];
727 0 : const int unitIndex = searchCode % 10;
728 0 : if ((unitCode == 1 && !(unitIndex == 0 || unitIndex == 1)) ||
729 0 : (unitCode == 2 &&
730 0 : !(unitIndex == 2 || unitIndex == 3 || unitIndex == 4)) ||
731 0 : (unitCode == 3 && !(unitIndex == 5 || unitIndex == 6)))
732 : {
733 0 : searchCode -= unitIndex;
734 0 : switch (unitIndex)
735 : {
736 0 : case 0:
737 : case 3:
738 : case 5:
739 0 : if (unitCode == 2)
740 0 : searchCode += 3;
741 0 : else if (unitCode == 3)
742 0 : searchCode += 5;
743 0 : break;
744 0 : case 1:
745 : case 2:
746 : case 6:
747 0 : if (unitCode == 1)
748 0 : searchCode += 1;
749 0 : if (unitCode == 2)
750 0 : searchCode += 2;
751 0 : else if (unitCode == 3)
752 0 : searchCode += 6;
753 0 : break;
754 0 : case 4:
755 : // FIXME? The following cond is not possible:
756 : // if( unitCode == 2 )
757 : // searchCode += 4;
758 0 : break;
759 : }
760 : }
761 0 : break;
762 : }
763 0 : }
764 : }
765 : else // Find state plane prj str by all inputs.
766 : {
767 3 : if (code < 0 || code > INT_MAX / 10)
768 0 : return OGRERR_FAILURE;
769 :
770 : // Need to have a special EPSG-ESRI zone code mapping first.
771 360 : for (int i = 0; statePlaneZoneMapping[i] != 0; i += 3)
772 : {
773 357 : if (code == statePlaneZoneMapping[i] &&
774 0 : (statePlaneZoneMapping[i + 1] == -1 ||
775 0 : pcsCode == statePlaneZoneMapping[i + 1]))
776 : {
777 0 : code = statePlaneZoneMapping[i + 2];
778 0 : break;
779 : }
780 : }
781 3 : searchCode = code * 10;
782 3 : if (!datumName)
783 : {
784 0 : CPLError(CE_Failure, CPLE_AppDefined, "datumName is NULL.");
785 0 : return OGRERR_FAILURE;
786 : }
787 3 : if (EQUAL(datumName, "HARN"))
788 : {
789 0 : if (EQUAL(unitsName, "international_feet"))
790 0 : searchCode += 5;
791 0 : else if (strstr(unitsName, "feet") || strstr(unitsName, "foot"))
792 0 : searchCode += 3;
793 : }
794 3 : else if (strstr(datumName, "NAD") && strstr(datumName, "83"))
795 : {
796 3 : if (EQUAL(unitsName, "meters"))
797 0 : searchCode += 1;
798 3 : else if (EQUAL(unitsName, "international_feet"))
799 0 : searchCode += 6;
800 3 : else if (strstr(unitsName, "feet") || strstr(unitsName, "foot"))
801 3 : searchCode += 2;
802 : }
803 0 : else if (strstr(datumName, "NAD") && strstr(datumName, "27") &&
804 0 : !EQUAL(unitsName, "meters"))
805 : {
806 0 : searchCode += 4;
807 : }
808 : else
809 0 : searchCode = -1;
810 : }
811 3 : if (searchCode > 0)
812 : {
813 3 : char codeS[20] = {};
814 3 : snprintf(codeS, sizeof(codeS), "%d", static_cast<int>(searchCode));
815 3 : return importFromDict("esri_StatePlane_extra.wkt", codeS);
816 : }
817 0 : return OGRERR_FAILURE;
818 : }
|
