Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: OpenGIS Simple Features Reference Implementation
4 : * Purpose: OGRSpatialReference translation from OziExplorer
5 : * georeferencing information.
6 : * Author: Andrey Kiselev, dron@ak4719.spb.edu
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 2009, Andrey Kiselev <dron@ak4719.spb.edu>
10 : * Copyright (c) 2009-2012, Even Rouault <even dot rouault at spatialys.com>
11 : *
12 : * SPDX-License-Identifier: MIT
13 : ****************************************************************************/
14 :
15 : #include "cpl_port.h"
16 : #include "ogr_spatialref.h"
17 :
18 : #include <cstdlib>
19 : #include <cstring>
20 :
21 : #include "cpl_conv.h"
22 : #include "cpl_csv.h"
23 : #include "cpl_error.h"
24 : #include "cpl_string.h"
25 : #include "ogr_core.h"
26 : #include "ogr_srs_api.h"
27 :
28 : /************************************************************************/
29 : /* OSRImportFromOzi() */
30 : /************************************************************************/
31 :
32 : /**
33 : * Import coordinate system from OziExplorer projection definition.
34 : *
35 : * This function will import projection definition in style, used by
36 : * OziExplorer software.
37 : *
38 : * @param hSRS spatial reference object.
39 : * @param papszLines Map file lines. This is an array of strings containing
40 : * the whole OziExplorer .MAP file. The array is terminated by a NULL pointer.
41 : *
42 : * @return OGRERR_NONE on success or an error code in case of failure.
43 : */
44 :
45 3 : OGRErr OSRImportFromOzi(OGRSpatialReferenceH hSRS,
46 : const char *const *papszLines)
47 :
48 : {
49 3 : VALIDATE_POINTER1(hSRS, "OSRImportFromOzi", OGRERR_FAILURE);
50 :
51 3 : return OGRSpatialReference::FromHandle(hSRS)->importFromOzi(papszLines);
52 : }
53 :
54 : /************************************************************************/
55 : /* importFromOzi() */
56 : /************************************************************************/
57 :
58 : /**
59 : * Import coordinate system from OziExplorer projection definition.
60 : *
61 : * This method will import projection definition in style, used by
62 : * OziExplorer software.
63 : *
64 : * @param papszLines Map file lines. This is an array of strings containing
65 : * the whole OziExplorer .MAP file. The array is terminated by a NULL pointer.
66 : *
67 : * @return OGRERR_NONE on success or an error code in case of failure.
68 : *
69 : */
70 :
71 3 : OGRErr OGRSpatialReference::importFromOzi(const char *const *papszLines)
72 : {
73 : const char *pszDatum;
74 3 : const char *pszProj = nullptr;
75 3 : const char *pszProjParams = nullptr;
76 :
77 3 : Clear();
78 :
79 3 : const int nLines = CSLCount(papszLines);
80 3 : if (nLines < 5)
81 0 : return OGRERR_NOT_ENOUGH_DATA;
82 :
83 3 : pszDatum = papszLines[4];
84 :
85 9 : for (int iLine = 5; iLine < nLines; iLine++)
86 : {
87 6 : if (STARTS_WITH_CI(papszLines[iLine], "Map Projection"))
88 : {
89 3 : pszProj = papszLines[iLine];
90 : }
91 3 : else if (STARTS_WITH_CI(papszLines[iLine], "Projection Setup"))
92 : {
93 3 : pszProjParams = papszLines[iLine];
94 : }
95 : }
96 :
97 3 : if (!(pszDatum && pszProj && pszProjParams))
98 0 : return OGRERR_NOT_ENOUGH_DATA;
99 :
100 : /* -------------------------------------------------------------------- */
101 : /* Operate on the basis of the projection name. */
102 : /* -------------------------------------------------------------------- */
103 3 : char **papszProj = CSLTokenizeStringComplex(pszProj, ",", TRUE, TRUE);
104 : char **papszProjParams =
105 3 : CSLTokenizeStringComplex(pszProjParams, ",", TRUE, TRUE);
106 3 : char **papszDatum = nullptr;
107 :
108 3 : if (CSLCount(papszProj) < 2)
109 : {
110 0 : goto not_enough_data;
111 : }
112 :
113 3 : if (STARTS_WITH_CI(papszProj[1], "Latitude/Longitude"))
114 : {
115 : // Do nothing.
116 : }
117 1 : else if (STARTS_WITH_CI(papszProj[1], "Mercator"))
118 : {
119 0 : if (CSLCount(papszProjParams) < 6)
120 0 : goto not_enough_data;
121 0 : double dfScale = CPLAtof(papszProjParams[3]);
122 : // If unset, default to scale = 1.
123 0 : if (papszProjParams[3][0] == 0)
124 0 : dfScale = 1;
125 0 : SetMercator(CPLAtof(papszProjParams[1]), CPLAtof(papszProjParams[2]),
126 0 : dfScale, CPLAtof(papszProjParams[4]),
127 0 : CPLAtof(papszProjParams[5]));
128 : }
129 1 : else if (STARTS_WITH_CI(papszProj[1], "Transverse Mercator"))
130 : {
131 0 : if (CSLCount(papszProjParams) < 6)
132 0 : goto not_enough_data;
133 0 : SetTM(CPLAtof(papszProjParams[1]), CPLAtof(papszProjParams[2]),
134 0 : CPLAtof(papszProjParams[3]), CPLAtof(papszProjParams[4]),
135 0 : CPLAtof(papszProjParams[5]));
136 : }
137 1 : else if (STARTS_WITH_CI(papszProj[1], "Lambert Conformal Conic"))
138 : {
139 1 : if (CSLCount(papszProjParams) < 8)
140 0 : goto not_enough_data;
141 1 : SetLCC(CPLAtof(papszProjParams[6]), CPLAtof(papszProjParams[7]),
142 1 : CPLAtof(papszProjParams[1]), CPLAtof(papszProjParams[2]),
143 1 : CPLAtof(papszProjParams[4]), CPLAtof(papszProjParams[5]));
144 : }
145 0 : else if (STARTS_WITH_CI(papszProj[1], "Sinusoidal"))
146 : {
147 0 : if (CSLCount(papszProjParams) < 6)
148 0 : goto not_enough_data;
149 0 : SetSinusoidal(CPLAtof(papszProjParams[2]), CPLAtof(papszProjParams[4]),
150 0 : CPLAtof(papszProjParams[5]));
151 : }
152 0 : else if (STARTS_WITH_CI(papszProj[1], "Albers Equal Area"))
153 : {
154 0 : if (CSLCount(papszProjParams) < 8)
155 0 : goto not_enough_data;
156 0 : SetACEA(CPLAtof(papszProjParams[6]), CPLAtof(papszProjParams[7]),
157 0 : CPLAtof(papszProjParams[1]), CPLAtof(papszProjParams[2]),
158 0 : CPLAtof(papszProjParams[4]), CPLAtof(papszProjParams[5]));
159 : }
160 0 : else if (STARTS_WITH_CI(papszProj[1],
161 0 : "(UTM) Universal Transverse Mercator") &&
162 : nLines > 5)
163 : {
164 : // Look for the UTM zone in the calibration point data.
165 0 : int iLine = 5; // Used after for.
166 0 : for (; iLine < nLines; iLine++)
167 : {
168 0 : if (STARTS_WITH_CI(papszLines[iLine], "Point"))
169 : {
170 0 : char **papszTok = CSLTokenizeString2(papszLines[iLine], ",",
171 : CSLT_ALLOWEMPTYTOKENS |
172 : CSLT_STRIPLEADSPACES |
173 : CSLT_STRIPENDSPACES);
174 0 : if (CSLCount(papszTok) < 17 || EQUAL(papszTok[2], "") ||
175 0 : EQUAL(papszTok[13], "") || EQUAL(papszTok[14], "") ||
176 0 : EQUAL(papszTok[15], "") || EQUAL(papszTok[16], ""))
177 : {
178 0 : CSLDestroy(papszTok);
179 0 : continue;
180 : }
181 0 : SetUTM(atoi(papszTok[13]), EQUAL(papszTok[16], "N"));
182 0 : CSLDestroy(papszTok);
183 0 : break;
184 : }
185 : }
186 0 : if (iLine == nLines) // Try to guess the UTM zone.
187 : {
188 0 : float fMinLongitude = 1000.0f;
189 0 : float fMaxLongitude = -1000.0f;
190 0 : float fMinLatitude = 1000.0f;
191 0 : float fMaxLatitude = -1000.0f;
192 0 : bool bFoundMMPLL = false;
193 0 : for (iLine = 5; iLine < nLines; iLine++)
194 : {
195 0 : if (STARTS_WITH_CI(papszLines[iLine], "MMPLL"))
196 : {
197 0 : char **papszTok = CSLTokenizeString2(
198 0 : papszLines[iLine], ",",
199 : CSLT_ALLOWEMPTYTOKENS | CSLT_STRIPLEADSPACES |
200 : CSLT_STRIPENDSPACES);
201 0 : if (CSLCount(papszTok) < 4)
202 : {
203 0 : CSLDestroy(papszTok);
204 0 : continue;
205 : }
206 : const float fLongitude =
207 0 : static_cast<float>(CPLAtofM(papszTok[2]));
208 : const float fLatitude =
209 0 : static_cast<float>(CPLAtofM(papszTok[3]));
210 0 : CSLDestroy(papszTok);
211 :
212 0 : bFoundMMPLL = true;
213 :
214 0 : if (fMinLongitude > fLongitude)
215 0 : fMinLongitude = fLongitude;
216 0 : if (fMaxLongitude < fLongitude)
217 0 : fMaxLongitude = fLongitude;
218 0 : if (fMinLatitude > fLatitude)
219 0 : fMinLatitude = fLatitude;
220 0 : if (fMaxLatitude < fLatitude)
221 0 : fMaxLatitude = fLatitude;
222 : }
223 : }
224 0 : const float fMedianLatitude = (fMinLatitude + fMaxLatitude) / 2;
225 0 : const float fMedianLongitude = (fMinLongitude + fMaxLongitude) / 2;
226 0 : if (bFoundMMPLL && fMaxLatitude <= 90)
227 : {
228 0 : int nUtmZone = 0;
229 0 : if (fMedianLatitude >= 56 && fMedianLatitude <= 64 &&
230 0 : fMedianLongitude >= 3 && fMedianLongitude <= 12)
231 0 : nUtmZone = 32; // Norway exception.
232 0 : else if (fMedianLatitude >= 72 && fMedianLatitude <= 84 &&
233 0 : fMedianLongitude >= 0 && fMedianLongitude <= 42)
234 : // Svalbard exception.
235 0 : nUtmZone =
236 0 : static_cast<int>((fMedianLongitude + 3) / 12) * 2 + 31;
237 : else
238 0 : nUtmZone =
239 0 : static_cast<int>((fMedianLongitude + 180) / 6) + 1;
240 0 : SetUTM(nUtmZone, fMedianLatitude >= 0);
241 : }
242 : else
243 : {
244 0 : CPLDebug("OSR_Ozi", "UTM Zone not found");
245 : }
246 0 : }
247 : }
248 0 : else if (STARTS_WITH_CI(papszProj[1], "(I) France Zone I"))
249 : {
250 0 : SetLCC1SP(49.5, 2.337229167, 0.99987734, 600000, 1200000);
251 : }
252 0 : else if (STARTS_WITH_CI(papszProj[1], "(II) France Zone II"))
253 : {
254 0 : SetLCC1SP(46.8, 2.337229167, 0.99987742, 600000, 2200000);
255 : }
256 0 : else if (STARTS_WITH_CI(papszProj[1], "(III) France Zone III"))
257 : {
258 0 : SetLCC1SP(44.1, 2.337229167, 0.99987750, 600000, 3200000);
259 : }
260 0 : else if (STARTS_WITH_CI(papszProj[1], "(IV) France Zone IV"))
261 : {
262 0 : SetLCC1SP(42.165, 2.337229167, 0.99994471, 234.358, 4185861.369);
263 : }
264 :
265 : /*
266 : * Note: The following projections have not been implemented yet
267 : *
268 : */
269 :
270 : /*
271 : else if( STARTS_WITH_CI(papszProj[1], "(BNG) British National Grid") )
272 : {
273 : }
274 : else if( STARTS_WITH_CI(papszProj[1], "(IG) Irish Grid") )
275 : {
276 : }
277 :
278 : else if( STARTS_WITH_CI(papszProj[1], "(NZG) New Zealand Grid") )
279 : {
280 : }
281 : else if( STARTS_WITH_CI(papszProj[1], "(NZTM2) New Zealand TM 2000") )
282 : {
283 : }
284 : else if( STARTS_WITH_CI(papszProj[1], "(SG) Swedish Grid") )
285 : {
286 : }
287 : else if( STARTS_WITH_CI(papszProj[1], "(SUI) Swiss Grid") )
288 : {
289 : }
290 : else if( STARTS_WITH_CI(papszProj[1], "(A)Lambert Azimuthual Equal
291 : Area") )
292 : {
293 : }
294 : else if( STARTS_WITH_CI(papszProj[1], "(EQC) Equidistant Conic") )
295 : {
296 : }
297 : else if( STARTS_WITH_CI(papszProj[1], "Polyconic (American)") )
298 : {
299 : }
300 : else if( STARTS_WITH_CI(papszProj[1], "Van Der Grinten") )
301 : {
302 : }
303 : else if( STARTS_WITH_CI(papszProj[1], "Vertical Near-Sided Perspective")
304 : )
305 : {
306 : }
307 : else if( STARTS_WITH_CI(papszProj[1], "(WIV) Wagner IV") )
308 : {
309 : }
310 : else if( STARTS_WITH_CI(papszProj[1], "Bonne") )
311 : {
312 : }
313 : else if( STARTS_WITH_CI(papszProj[1],
314 : "(MT0) Montana State Plane Zone 2500") )
315 : {
316 : }
317 : else if( STARTS_WITH_CI(papszProj[1], "ITA1) Italy Grid Zone 1") )
318 : {
319 : }
320 : else if( STARTS_WITH_CI(papszProj[1], "ITA2) Italy Grid Zone 2") )
321 : {
322 : }
323 : else if( STARTS_WITH_CI(papszProj[1],
324 : "(VICMAP-TM) Victoria Aust.(pseudo AMG)") )
325 : {
326 : }
327 : else if( STARTS_WITH_CI(papszProj[1], "VICGRID) Victoria Australia") )
328 : {
329 : }
330 : else if( STARTS_WITH_CI(papszProj[1],
331 : "(VG94) VICGRID94 Victoria Australia") )
332 : {
333 : }
334 : else if( STARTS_WITH_CI(papszProj[1], "Gnomonic") )
335 : {
336 : }
337 : */
338 : else
339 : {
340 0 : CPLDebug("OSR_Ozi", "Unsupported projection: \"%s\"", papszProj[1]);
341 0 : SetLocalCS(
342 0 : CPLString().Printf(R"("Ozi" projection "%s")", papszProj[1]));
343 : }
344 :
345 : /* -------------------------------------------------------------------- */
346 : /* Try to translate the datum/spheroid. */
347 : /* -------------------------------------------------------------------- */
348 3 : papszDatum = CSLTokenizeString2(
349 : pszDatum, ",",
350 : CSLT_ALLOWEMPTYTOKENS | CSLT_STRIPLEADSPACES | CSLT_STRIPENDSPACES);
351 3 : if (papszDatum == nullptr)
352 0 : goto not_enough_data;
353 :
354 3 : if (!IsLocal())
355 : {
356 : /* --------------------------------------------------------------------
357 : */
358 : /* Verify that we can find the CSV file containing the datums */
359 : /* --------------------------------------------------------------------
360 : */
361 3 : if (CSVScanFileByName(CSVFilename("ozi_datum.csv"), "EPSG_DATUM_CODE",
362 3 : "4326", CC_Integer) == nullptr)
363 : {
364 0 : CPLError(CE_Failure, CPLE_OpenFailed,
365 : "Unable to open OZI support file %s. "
366 : "Try setting the GDAL_DATA environment variable to point "
367 : "to the directory containing OZI csv files.",
368 : CSVFilename("ozi_datum.csv"));
369 0 : goto other_error;
370 : }
371 :
372 : /* --------------------------------------------------------------------
373 : */
374 : /* Search for matching datum */
375 : /* --------------------------------------------------------------------
376 : */
377 3 : const char *pszOziDatum = CSVFilename("ozi_datum.csv");
378 : CPLString osDName = CSVGetField(pszOziDatum, "NAME", papszDatum[0],
379 3 : CC_ApproxString, "NAME");
380 3 : if (osDName.empty())
381 : {
382 0 : CPLError(CE_Failure, CPLE_AppDefined,
383 : "Failed to find datum %s in ozi_datum.csv.",
384 : papszDatum[0]);
385 0 : goto other_error;
386 : }
387 :
388 : const int nDatumCode =
389 3 : atoi(CSVGetField(pszOziDatum, "NAME", papszDatum[0],
390 : CC_ApproxString, "EPSG_DATUM_CODE"));
391 :
392 3 : if (nDatumCode > 0) // There is a matching EPSG code
393 : {
394 4 : OGRSpatialReference oGCS;
395 2 : oGCS.importFromEPSG(nDatumCode);
396 2 : CopyGeogCSFrom(&oGCS);
397 : }
398 : else // We use the parameters from the CSV files
399 : {
400 : CPLString osEllipseCode =
401 : CSVGetField(pszOziDatum, "NAME", papszDatum[0], CC_ApproxString,
402 1 : "ELLIPSOID_CODE");
403 1 : const double dfDeltaX = CPLAtof(CSVGetField(
404 : pszOziDatum, "NAME", papszDatum[0], CC_ApproxString, "DELTAX"));
405 1 : const double dfDeltaY = CPLAtof(CSVGetField(
406 : pszOziDatum, "NAME", papszDatum[0], CC_ApproxString, "DELTAY"));
407 1 : const double dfDeltaZ = CPLAtof(CSVGetField(
408 : pszOziDatum, "NAME", papszDatum[0], CC_ApproxString, "DELTAZ"));
409 :
410 : /* --------------------------------------------------------------------
411 : */
412 : /* Verify that we can find the CSV file containing the
413 : * ellipsoids. */
414 : /* --------------------------------------------------------------------
415 : */
416 1 : if (CSVScanFileByName(CSVFilename("ozi_ellips.csv"),
417 : "ELLIPSOID_CODE", "20",
418 1 : CC_Integer) == nullptr)
419 : {
420 0 : CPLError(
421 : CE_Failure, CPLE_OpenFailed,
422 : "Unable to open OZI support file %s. "
423 : "Try setting the GDAL_DATA environment variable to point "
424 : "to the directory containing OZI csv files.",
425 : CSVFilename("ozi_ellips.csv"));
426 0 : goto other_error;
427 : }
428 :
429 : /* --------------------------------------------------------------------
430 : */
431 : /* Lookup the ellipse code. */
432 : /* --------------------------------------------------------------------
433 : */
434 1 : const char *pszOziEllipse = CSVFilename("ozi_ellips.csv");
435 :
436 : CPLString osEName =
437 : CSVGetField(pszOziEllipse, "ELLIPSOID_CODE", osEllipseCode,
438 1 : CC_ApproxString, "NAME");
439 1 : if (osEName.empty())
440 : {
441 0 : CPLError(CE_Failure, CPLE_AppDefined,
442 : "Failed to find ellipsoid %s in ozi_ellips.csv.",
443 : osEllipseCode.c_str());
444 0 : goto other_error;
445 : }
446 :
447 : const double dfA =
448 1 : CPLAtof(CSVGetField(pszOziEllipse, "ELLIPSOID_CODE",
449 : osEllipseCode, CC_ApproxString, "A"));
450 : const double dfInvF =
451 1 : CPLAtof(CSVGetField(pszOziEllipse, "ELLIPSOID_CODE",
452 : osEllipseCode, CC_ApproxString, "INVF"));
453 :
454 : /* --------------------------------------------------------------------
455 : */
456 : /* Create geographic coordinate system. */
457 : /* --------------------------------------------------------------------
458 : */
459 1 : SetGeogCS(osDName, osDName, osEName, dfA, dfInvF);
460 1 : SetTOWGS84(dfDeltaX, dfDeltaY, dfDeltaZ);
461 : }
462 : }
463 :
464 : /* -------------------------------------------------------------------- */
465 : /* Grid units translation */
466 : /* -------------------------------------------------------------------- */
467 3 : if (IsLocal() || IsProjected())
468 1 : SetLinearUnits(SRS_UL_METER, 1.0);
469 :
470 3 : CSLDestroy(papszProj);
471 3 : CSLDestroy(papszProjParams);
472 3 : CSLDestroy(papszDatum);
473 :
474 3 : return OGRERR_NONE;
475 :
476 0 : not_enough_data:
477 :
478 0 : CSLDestroy(papszProj);
479 0 : CSLDestroy(papszProjParams);
480 0 : CSLDestroy(papszDatum);
481 :
482 0 : return OGRERR_NOT_ENOUGH_DATA;
483 :
484 0 : other_error:
485 :
486 0 : CSLDestroy(papszProj);
487 0 : CSLDestroy(papszProjParams);
488 0 : CSLDestroy(papszDatum);
489 :
490 0 : return OGRERR_FAILURE;
491 : }
|
