Coverage Report

Created: 2025-06-13 06:18

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