/src/gdal/frmts/gtiff/gt_wkt_srs.cpp
Line | Count | Source (jump to first uncovered line) |
1 | | /****************************************************************************** |
2 | | * |
3 | | * Project: GeoTIFF Driver |
4 | | * Purpose: Implements translation between GeoTIFF normalized projection |
5 | | * definitions and OpenGIS WKT SRS format. |
6 | | * Author: Frank Warmerdam, warmerdam@pobox.com |
7 | | * |
8 | | ****************************************************************************** |
9 | | * Copyright (c) 1999, Frank Warmerdam |
10 | | * Copyright (c) 2008-2014, Even Rouault <even dot rouault at spatialys.com> |
11 | | * |
12 | | * SPDX-License-Identifier: MIT |
13 | | ****************************************************************************/ |
14 | | |
15 | | #include "cpl_port.h" |
16 | | |
17 | | #include "gt_wkt_srs.h" |
18 | | |
19 | | #include <cmath> |
20 | | #include <cstdio> |
21 | | #include <cstdlib> |
22 | | #include <cstring> |
23 | | |
24 | | #include <algorithm> |
25 | | #include <mutex> |
26 | | |
27 | | #include "cpl_conv.h" |
28 | | #include "cpl_error.h" |
29 | | #include "cpl_string.h" |
30 | | #include "cpl_vsi.h" |
31 | | #include "gt_citation.h" |
32 | | #include "gt_wkt_srs_for_gdal.h" |
33 | | #include "gt_wkt_srs_priv.h" |
34 | | #include "gtiff.h" |
35 | | #include "gdal.h" |
36 | | #include "geokeys.h" |
37 | | #include "geovalues.h" |
38 | | #include "ogr_core.h" |
39 | | #include "ogr_spatialref.h" |
40 | | #include "ogr_srs_api.h" |
41 | | #include "ogr_proj_p.h" |
42 | | #include "tiff.h" |
43 | | #include "tiffio.h" |
44 | | #include "tifvsi.h" |
45 | | #include "xtiffio.h" |
46 | | |
47 | | #include "proj.h" |
48 | | |
49 | | static const geokey_t ProjLinearUnitsInterpCorrectGeoKey = |
50 | | static_cast<geokey_t>(3059); |
51 | | |
52 | | #ifndef CT_HotineObliqueMercatorAzimuthCenter |
53 | 0 | #define CT_HotineObliqueMercatorAzimuthCenter 9815 |
54 | | #endif |
55 | | |
56 | | #if !defined(GTIFAtof) |
57 | 0 | #define GTIFAtof CPLAtof |
58 | | #endif |
59 | | |
60 | | // To remind myself not to use CPLString in this file! |
61 | | #define CPLString Please_do_not_use_CPLString_in_this_file |
62 | | |
63 | | static const char *const papszDatumEquiv[] = { |
64 | | "Militar_Geographische_Institut", |
65 | | "Militar_Geographische_Institute", |
66 | | "World_Geodetic_System_1984", |
67 | | "WGS_1984", |
68 | | "WGS_72_Transit_Broadcast_Ephemeris", |
69 | | "WGS_1972_Transit_Broadcast_Ephemeris", |
70 | | "World_Geodetic_System_1972", |
71 | | "WGS_1972", |
72 | | "European_Terrestrial_Reference_System_89", |
73 | | "European_Reference_System_1989", |
74 | | "D_North_American_1927", |
75 | | "North_American_Datum_1927", // #6863 |
76 | | nullptr}; |
77 | | |
78 | | // Older libgeotiff's won't list this. |
79 | | #ifndef CT_CylindricalEqualArea |
80 | 0 | #define CT_CylindricalEqualArea 28 |
81 | | #endif |
82 | | |
83 | | #if LIBGEOTIFF_VERSION < 1700 |
84 | | constexpr geokey_t CoordinateEpochGeoKey = static_cast<geokey_t>(5120); |
85 | | #endif |
86 | | |
87 | | // Exists since 8.0.1 |
88 | | #ifndef PROJ_AT_LEAST_VERSION |
89 | | #define PROJ_COMPUTE_VERSION(maj, min, patch) \ |
90 | | ((maj)*10000 + (min)*100 + (patch)) |
91 | | #define PROJ_VERSION_NUMBER \ |
92 | | PROJ_COMPUTE_VERSION(PROJ_VERSION_MAJOR, PROJ_VERSION_MINOR, \ |
93 | | PROJ_VERSION_PATCH) |
94 | | #define PROJ_AT_LEAST_VERSION(maj, min, patch) \ |
95 | | (PROJ_VERSION_NUMBER >= PROJ_COMPUTE_VERSION(maj, min, patch)) |
96 | | #endif |
97 | | |
98 | | /************************************************************************/ |
99 | | /* LibgeotiffOneTimeInit() */ |
100 | | /************************************************************************/ |
101 | | |
102 | | static std::mutex oDeleteMutex; |
103 | | |
104 | | void LibgeotiffOneTimeInit() |
105 | 0 | { |
106 | 0 | std::lock_guard<std::mutex> oLock(oDeleteMutex); |
107 | |
|
108 | 0 | static bool bOneTimeInitDone = false; |
109 | |
|
110 | 0 | if (bOneTimeInitDone) |
111 | 0 | return; |
112 | | |
113 | 0 | bOneTimeInitDone = true; |
114 | | |
115 | | // This isn't thread-safe, so better do it now |
116 | 0 | XTIFFInitialize(); |
117 | 0 | } |
118 | | |
119 | | /************************************************************************/ |
120 | | /* GTIFToCPLRecyleString() */ |
121 | | /* */ |
122 | | /* This changes a string from the libgeotiff heap to the GDAL */ |
123 | | /* heap. */ |
124 | | /************************************************************************/ |
125 | | |
126 | | static void GTIFToCPLRecycleString(char **ppszTarget) |
127 | | |
128 | 0 | { |
129 | 0 | if (*ppszTarget == nullptr) |
130 | 0 | return; |
131 | | |
132 | 0 | char *pszTempString = CPLStrdup(*ppszTarget); |
133 | 0 | GTIFFreeMemory(*ppszTarget); |
134 | 0 | *ppszTarget = pszTempString; |
135 | 0 | } |
136 | | |
137 | | /************************************************************************/ |
138 | | /* WKTMassageDatum() */ |
139 | | /* */ |
140 | | /* Massage an EPSG datum name into WMT format. Also transform */ |
141 | | /* specific exception cases into WKT versions. */ |
142 | | /************************************************************************/ |
143 | | |
144 | | static void WKTMassageDatum(char **ppszDatum) |
145 | | |
146 | 0 | { |
147 | 0 | char *pszDatum = *ppszDatum; |
148 | 0 | if (!pszDatum || pszDatum[0] == '\0') |
149 | 0 | return; |
150 | | |
151 | | /* -------------------------------------------------------------------- */ |
152 | | /* Translate non-alphanumeric values to underscores. */ |
153 | | /* -------------------------------------------------------------------- */ |
154 | 0 | for (int i = 0; pszDatum[i] != '\0'; i++) |
155 | 0 | { |
156 | 0 | if (pszDatum[i] != '+' && !(pszDatum[i] >= 'A' && pszDatum[i] <= 'Z') && |
157 | 0 | !(pszDatum[i] >= 'a' && pszDatum[i] <= 'z') && |
158 | 0 | !(pszDatum[i] >= '0' && pszDatum[i] <= '9')) |
159 | 0 | { |
160 | 0 | pszDatum[i] = '_'; |
161 | 0 | } |
162 | 0 | } |
163 | | |
164 | | /* -------------------------------------------------------------------- */ |
165 | | /* Remove repeated and trailing underscores. */ |
166 | | /* -------------------------------------------------------------------- */ |
167 | 0 | int j = 0; // Used after for. |
168 | 0 | for (int i = 1; pszDatum[i] != '\0'; i++) |
169 | 0 | { |
170 | 0 | if (pszDatum[j] == '_' && pszDatum[i] == '_') |
171 | 0 | continue; |
172 | | |
173 | 0 | pszDatum[++j] = pszDatum[i]; |
174 | 0 | } |
175 | 0 | if (pszDatum[j] == '_') |
176 | 0 | pszDatum[j] = '\0'; |
177 | 0 | else |
178 | 0 | pszDatum[j + 1] = '\0'; |
179 | | |
180 | | /* -------------------------------------------------------------------- */ |
181 | | /* Search for datum equivalences. Specific massaged names get */ |
182 | | /* mapped to OpenGIS specified names. */ |
183 | | /* -------------------------------------------------------------------- */ |
184 | 0 | for (int i = 0; papszDatumEquiv[i] != nullptr; i += 2) |
185 | 0 | { |
186 | 0 | if (EQUAL(*ppszDatum, papszDatumEquiv[i])) |
187 | 0 | { |
188 | 0 | CPLFree(*ppszDatum); |
189 | 0 | *ppszDatum = CPLStrdup(papszDatumEquiv[i + 1]); |
190 | 0 | return; |
191 | 0 | } |
192 | 0 | } |
193 | 0 | } |
194 | | |
195 | | /************************************************************************/ |
196 | | /* GTIFCleanupImageineNames() */ |
197 | | /* */ |
198 | | /* Erdas Imagine sometimes emits big copyright messages, and */ |
199 | | /* other stuff into citations. These can be pretty messy when */ |
200 | | /* turned into WKT, so we try to trim and clean the strings */ |
201 | | /* somewhat. */ |
202 | | /************************************************************************/ |
203 | | |
204 | | /* For example: |
205 | | GTCitationGeoKey (Ascii,215): "IMAGINE GeoTIFF Support\nCopyright 1991 - 2001 |
206 | | by ERDAS, Inc. All Rights Reserved\n@(#)$RCSfile$ $Revision: 34309 $ $Date: |
207 | | 2016-05-29 11:29:40 -0700 (Sun, 29 May 2016) $\nProjection Name = UTM\nUnits |
208 | | = meters\nGeoTIFF Units = meters" |
209 | | |
210 | | GeogCitationGeoKey (Ascii,267): "IMAGINE GeoTIFF Support\nCopyright 1991 - |
211 | | 2001 by ERDAS, Inc. All Rights Reserved\n@(#)$RCSfile$ $Revision: 34309 $ |
212 | | $Date: 2016-05-29 11:29:40 -0700 (Sun, 29 May 2016) $\nUnable to match |
213 | | Ellipsoid (Datum) to a GeographicTypeGeoKey value\nEllipsoid = Clarke |
214 | | 1866\nDatum = NAD27 (CONUS)" |
215 | | |
216 | | PCSCitationGeoKey (Ascii,214): "IMAGINE GeoTIFF Support\nCopyright 1991 - |
217 | | 2001 by ERDAS, Inc. All Rights Reserved\n@(#)$RCSfile$ $Revision: 34309 $ |
218 | | $Date: 2016-05-29 11:29:40 -0700 (Sun, 29 May 2016) $\nUTM Zone |
219 | | 10N\nEllipsoid = Clarke 1866\nDatum = NAD27 (CONUS)" |
220 | | */ |
221 | | |
222 | | static void GTIFCleanupImagineNames(char *pszCitation) |
223 | | |
224 | 0 | { |
225 | 0 | if (strstr(pszCitation, "IMAGINE GeoTIFF") == nullptr) |
226 | 0 | return; |
227 | | |
228 | | /* -------------------------------------------------------------------- */ |
229 | | /* First, we skip past all the copyright, and RCS stuff. We */ |
230 | | /* assume that this will have a "$" at the end of it all. */ |
231 | | /* -------------------------------------------------------------------- */ |
232 | 0 | char *pszSkip = pszCitation + strlen(pszCitation) - 1; |
233 | |
|
234 | 0 | for (; pszSkip != pszCitation && *pszSkip != '$'; pszSkip--) |
235 | 0 | { |
236 | 0 | } |
237 | |
|
238 | 0 | if (*pszSkip == '$') |
239 | 0 | pszSkip++; |
240 | 0 | if (*pszSkip == '\n') |
241 | 0 | pszSkip++; |
242 | |
|
243 | 0 | memmove(pszCitation, pszSkip, strlen(pszSkip) + 1); |
244 | | |
245 | | /* -------------------------------------------------------------------- */ |
246 | | /* Convert any newlines into spaces, they really gum up the */ |
247 | | /* WKT. */ |
248 | | /* -------------------------------------------------------------------- */ |
249 | 0 | for (int i = 0; pszCitation[i] != '\0'; i++) |
250 | 0 | { |
251 | 0 | if (pszCitation[i] == '\n') |
252 | 0 | pszCitation[i] = ' '; |
253 | 0 | } |
254 | 0 | } |
255 | | |
256 | | #if LIBGEOTIFF_VERSION < 1600 |
257 | | |
258 | | /************************************************************************/ |
259 | | /* GDALGTIFKeyGet() */ |
260 | | /************************************************************************/ |
261 | | |
262 | | static int GDALGTIFKeyGet(GTIF *hGTIF, geokey_t key, void *pData, int nIndex, |
263 | | int nCount, tagtype_t expected_tagtype) |
264 | | { |
265 | | tagtype_t tagtype = TYPE_UNKNOWN; |
266 | | if (!GTIFKeyInfo(hGTIF, key, nullptr, &tagtype)) |
267 | | return 0; |
268 | | if (tagtype != expected_tagtype) |
269 | | { |
270 | | CPLError(CE_Warning, CPLE_AppDefined, |
271 | | "Expected key %s to be of type %s. Got %s", GTIFKeyName(key), |
272 | | GTIFTypeName(expected_tagtype), GTIFTypeName(tagtype)); |
273 | | return 0; |
274 | | } |
275 | | return GTIFKeyGet(hGTIF, key, pData, nIndex, nCount); |
276 | | } |
277 | | |
278 | | /************************************************************************/ |
279 | | /* GDALGTIFKeyGetASCII() */ |
280 | | /************************************************************************/ |
281 | | |
282 | | int GDALGTIFKeyGetASCII(GTIF *hGTIF, geokey_t key, char *szStr, int szStrMaxLen) |
283 | | { |
284 | | return GDALGTIFKeyGet(hGTIF, key, szStr, 0, szStrMaxLen, TYPE_ASCII); |
285 | | } |
286 | | |
287 | | /************************************************************************/ |
288 | | /* GDALGTIFKeyGetSHORT() */ |
289 | | /************************************************************************/ |
290 | | |
291 | | int GDALGTIFKeyGetSHORT(GTIF *hGTIF, geokey_t key, unsigned short *pnVal, |
292 | | int nIndex, int nCount) |
293 | | { |
294 | | return GDALGTIFKeyGet(hGTIF, key, pnVal, nIndex, nCount, TYPE_SHORT); |
295 | | } |
296 | | |
297 | | /************************************************************************/ |
298 | | /* GDALGTIFKeyGetDOUBLE() */ |
299 | | /************************************************************************/ |
300 | | |
301 | | int GDALGTIFKeyGetDOUBLE(GTIF *hGTIF, geokey_t key, double *pdfVal, int nIndex, |
302 | | int nCount) |
303 | | { |
304 | | return GDALGTIFKeyGet(hGTIF, key, pdfVal, nIndex, nCount, TYPE_DOUBLE); |
305 | | } |
306 | | |
307 | | #endif |
308 | | |
309 | | /************************************************************************/ |
310 | | /* FillCompoundCRSWithManualVertCS() */ |
311 | | /************************************************************************/ |
312 | | |
313 | | static void FillCompoundCRSWithManualVertCS(GTIF *hGTIF, |
314 | | OGRSpatialReference &oSRS, |
315 | | const char *pszVertCSName, |
316 | | int verticalDatum, |
317 | | int verticalUnits) |
318 | 0 | { |
319 | | /* -------------------------------------------------------------------- */ |
320 | | /* Setup VERT_CS with citation if present. */ |
321 | | /* -------------------------------------------------------------------- */ |
322 | 0 | oSRS.SetNode("COMPD_CS|VERT_CS", pszVertCSName); |
323 | | |
324 | | /* -------------------------------------------------------------------- */ |
325 | | /* Setup the vertical datum. */ |
326 | | /* -------------------------------------------------------------------- */ |
327 | 0 | std::string osVDatumName = "unknown"; |
328 | 0 | const char *pszVDatumType = "2005"; // CS_VD_GeoidModelDerived |
329 | 0 | std::string osVDatumAuthName; |
330 | 0 | int nVDatumCode = 0; |
331 | |
|
332 | 0 | if (verticalDatum > 0 && verticalDatum != KvUserDefined) |
333 | 0 | { |
334 | 0 | osVDatumAuthName = "EPSG"; |
335 | 0 | nVDatumCode = verticalDatum; |
336 | |
|
337 | 0 | char szCode[12]; |
338 | 0 | snprintf(szCode, sizeof(szCode), "%d", verticalDatum); |
339 | 0 | auto ctx = |
340 | 0 | static_cast<PJ_CONTEXT *>(GTIFGetPROJContext(hGTIF, true, nullptr)); |
341 | 0 | auto datum = proj_create_from_database(ctx, "EPSG", szCode, |
342 | 0 | PJ_CATEGORY_DATUM, 0, nullptr); |
343 | 0 | if (datum) |
344 | 0 | { |
345 | 0 | const char *pszName = proj_get_name(datum); |
346 | 0 | if (pszName) |
347 | 0 | { |
348 | 0 | osVDatumName = pszName; |
349 | 0 | } |
350 | 0 | proj_destroy(datum); |
351 | 0 | } |
352 | 0 | } |
353 | 0 | else if (verticalDatum == KvUserDefined) |
354 | 0 | { |
355 | | // If the vertical datum is unknown, try to find the vertical CRS |
356 | | // from the database, and extra the datum information from it. |
357 | 0 | auto ctx = |
358 | 0 | static_cast<PJ_CONTEXT *>(GTIFGetPROJContext(hGTIF, true, nullptr)); |
359 | 0 | const auto type = PJ_TYPE_VERTICAL_CRS; |
360 | 0 | auto list = proj_create_from_name(ctx, nullptr, pszVertCSName, &type, 1, |
361 | 0 | true, // exact match |
362 | 0 | -1, // result set limit size, |
363 | 0 | nullptr); |
364 | 0 | if (list) |
365 | 0 | { |
366 | | // If we have several matches, check they all refer to the |
367 | | // same datum |
368 | 0 | bool bGoOn = true; |
369 | 0 | int ncount = proj_list_get_count(list); |
370 | 0 | for (int i = 0; bGoOn && i < ncount; ++i) |
371 | 0 | { |
372 | 0 | auto crs = proj_list_get(ctx, list, i); |
373 | 0 | if (crs) |
374 | 0 | { |
375 | 0 | auto datum = proj_crs_get_datum(ctx, crs); |
376 | 0 | if (datum) |
377 | 0 | { |
378 | 0 | osVDatumName = proj_get_name(datum); |
379 | 0 | const char *pszAuthName = |
380 | 0 | proj_get_id_auth_name(datum, 0); |
381 | 0 | const char *pszCode = proj_get_id_code(datum, 0); |
382 | 0 | if (pszCode && atoi(pszCode) && pszAuthName) |
383 | 0 | { |
384 | 0 | if (osVDatumAuthName.empty()) |
385 | 0 | { |
386 | 0 | osVDatumAuthName = pszAuthName; |
387 | 0 | nVDatumCode = atoi(pszCode); |
388 | 0 | } |
389 | 0 | else if (osVDatumAuthName != pszAuthName || |
390 | 0 | nVDatumCode != atoi(pszCode)) |
391 | 0 | { |
392 | 0 | osVDatumAuthName.clear(); |
393 | 0 | nVDatumCode = 0; |
394 | 0 | bGoOn = false; |
395 | 0 | } |
396 | 0 | } |
397 | 0 | proj_destroy(datum); |
398 | 0 | } |
399 | 0 | proj_destroy(crs); |
400 | 0 | } |
401 | 0 | } |
402 | 0 | } |
403 | 0 | proj_list_destroy(list); |
404 | 0 | } |
405 | |
|
406 | 0 | oSRS.SetNode("COMPD_CS|VERT_CS|VERT_DATUM", osVDatumName.c_str()); |
407 | 0 | oSRS.GetAttrNode("COMPD_CS|VERT_CS|VERT_DATUM") |
408 | 0 | ->AddChild(new OGR_SRSNode(pszVDatumType)); |
409 | 0 | if (!osVDatumAuthName.empty()) |
410 | 0 | oSRS.SetAuthority("COMPD_CS|VERT_CS|VERT_DATUM", |
411 | 0 | osVDatumAuthName.c_str(), nVDatumCode); |
412 | | |
413 | | /* -------------------------------------------------------------------- */ |
414 | | /* Set the vertical units. */ |
415 | | /* -------------------------------------------------------------------- */ |
416 | 0 | if (verticalUnits > 0 && verticalUnits != KvUserDefined && |
417 | 0 | verticalUnits != 9001) |
418 | 0 | { |
419 | 0 | char szCode[12]; |
420 | 0 | snprintf(szCode, sizeof(szCode), "%d", verticalUnits); |
421 | 0 | auto ctx = |
422 | 0 | static_cast<PJ_CONTEXT *>(GTIFGetPROJContext(hGTIF, true, nullptr)); |
423 | 0 | const char *pszName = nullptr; |
424 | 0 | double dfInMeters = 0.0; |
425 | 0 | if (proj_uom_get_info_from_database(ctx, "EPSG", szCode, &pszName, |
426 | 0 | &dfInMeters, nullptr)) |
427 | 0 | { |
428 | 0 | if (pszName) |
429 | 0 | oSRS.SetNode("COMPD_CS|VERT_CS|UNIT", pszName); |
430 | |
|
431 | 0 | char szInMeters[128] = {}; |
432 | 0 | CPLsnprintf(szInMeters, sizeof(szInMeters), "%.16g", dfInMeters); |
433 | 0 | oSRS.GetAttrNode("COMPD_CS|VERT_CS|UNIT") |
434 | 0 | ->AddChild(new OGR_SRSNode(szInMeters)); |
435 | 0 | } |
436 | |
|
437 | 0 | oSRS.SetAuthority("COMPD_CS|VERT_CS|UNIT", "EPSG", verticalUnits); |
438 | 0 | } |
439 | 0 | else |
440 | 0 | { |
441 | 0 | oSRS.SetNode("COMPD_CS|VERT_CS|UNIT", "metre"); |
442 | 0 | oSRS.GetAttrNode("COMPD_CS|VERT_CS|UNIT") |
443 | 0 | ->AddChild(new OGR_SRSNode("1.0")); |
444 | 0 | oSRS.SetAuthority("COMPD_CS|VERT_CS|UNIT", "EPSG", 9001); |
445 | 0 | } |
446 | | |
447 | | /* -------------------------------------------------------------------- */ |
448 | | /* Set the axis and VERT_CS authority. */ |
449 | | /* -------------------------------------------------------------------- */ |
450 | 0 | oSRS.SetNode("COMPD_CS|VERT_CS|AXIS", "Up"); |
451 | 0 | oSRS.GetAttrNode("COMPD_CS|VERT_CS|AXIS")->AddChild(new OGR_SRSNode("UP")); |
452 | 0 | } |
453 | | |
454 | | /************************************************************************/ |
455 | | /* GTIFGetEPSGOfficialName() */ |
456 | | /************************************************************************/ |
457 | | |
458 | | static char *GTIFGetEPSGOfficialName(GTIF *hGTIF, PJ_TYPE searchType, |
459 | | const char *pszName) |
460 | 0 | { |
461 | 0 | char *pszRet = nullptr; |
462 | | /* Search in database the corresponding EPSG 'official' name */ |
463 | 0 | auto ctx = |
464 | 0 | static_cast<PJ_CONTEXT *>(GTIFGetPROJContext(hGTIF, true, nullptr)); |
465 | 0 | auto list = proj_create_from_name(ctx, "EPSG", pszName, &searchType, 1, |
466 | 0 | false, /* approximate match */ |
467 | 0 | 1, nullptr); |
468 | 0 | if (list) |
469 | 0 | { |
470 | 0 | const auto listSize = proj_list_get_count(list); |
471 | 0 | if (listSize == 1) |
472 | 0 | { |
473 | 0 | auto obj = proj_list_get(ctx, list, 0); |
474 | 0 | if (obj) |
475 | 0 | { |
476 | 0 | const char *pszOfficialName = proj_get_name(obj); |
477 | 0 | if (pszOfficialName) |
478 | 0 | { |
479 | 0 | pszRet = CPLStrdup(pszOfficialName); |
480 | 0 | } |
481 | 0 | } |
482 | 0 | proj_destroy(obj); |
483 | 0 | } |
484 | 0 | proj_list_destroy(list); |
485 | 0 | } |
486 | 0 | return pszRet; |
487 | 0 | } |
488 | | |
489 | | /************************************************************************/ |
490 | | /* GTIFGetOGISDefnAsOSR() */ |
491 | | /************************************************************************/ |
492 | | |
493 | | OGRSpatialReferenceH GTIFGetOGISDefnAsOSR(GTIF *hGTIF, GTIFDefn *psDefn) |
494 | | |
495 | 0 | { |
496 | 0 | OGRSpatialReference oSRS; |
497 | |
|
498 | 0 | LibgeotiffOneTimeInit(); |
499 | |
|
500 | 0 | #if LIBGEOTIFF_VERSION >= 1600 |
501 | 0 | void *projContext = GTIFGetPROJContext(hGTIF, FALSE, nullptr); |
502 | 0 | #endif |
503 | | |
504 | | /* -------------------------------------------------------------------- */ |
505 | | /* Handle non-standard coordinate systems where GTModelTypeGeoKey */ |
506 | | /* is not defined, but ProjectedCSTypeGeoKey is defined (ticket #3019) */ |
507 | | /* -------------------------------------------------------------------- */ |
508 | 0 | if (psDefn->Model == KvUserDefined && psDefn->PCS != KvUserDefined) |
509 | 0 | { |
510 | 0 | psDefn->Model = ModelTypeProjected; |
511 | 0 | } |
512 | | |
513 | | /* ==================================================================== */ |
514 | | /* Read keys related to vertical component. */ |
515 | | /* ==================================================================== */ |
516 | 0 | unsigned short verticalCSType = 0; |
517 | 0 | unsigned short verticalDatum = 0; |
518 | 0 | unsigned short verticalUnits = 0; |
519 | |
|
520 | 0 | GDALGTIFKeyGetSHORT(hGTIF, VerticalCSTypeGeoKey, &verticalCSType, 0, 1); |
521 | 0 | GDALGTIFKeyGetSHORT(hGTIF, VerticalDatumGeoKey, &verticalDatum, 0, 1); |
522 | 0 | GDALGTIFKeyGetSHORT(hGTIF, VerticalUnitsGeoKey, &verticalUnits, 0, 1); |
523 | |
|
524 | 0 | if (verticalCSType != 0 || verticalDatum != 0 || verticalUnits != 0) |
525 | 0 | { |
526 | 0 | int versions[3]; |
527 | 0 | GTIFDirectoryInfo(hGTIF, versions, nullptr); |
528 | | // GeoTIFF 1.0 |
529 | 0 | if (versions[0] == 1 && versions[1] == 1 && versions[2] == 0) |
530 | 0 | { |
531 | | /* -------------------------------------------------------------------- |
532 | | */ |
533 | | /* The original geotiff specification appears to have */ |
534 | | /* misconstrued the EPSG codes 5101 to 5106 to be vertical */ |
535 | | /* coordinate system codes, when in fact they are vertical */ |
536 | | /* datum codes. So if these are found in the */ |
537 | | /* VerticalCSTypeGeoKey move them to the VerticalDatumGeoKey */ |
538 | | /* and insert the "normal" corresponding VerticalCSTypeGeoKey |
539 | | */ |
540 | | /* value. */ |
541 | | /* -------------------------------------------------------------------- |
542 | | */ |
543 | 0 | if ((verticalCSType >= 5101 && verticalCSType <= 5112) && |
544 | 0 | verticalDatum == 0) |
545 | 0 | { |
546 | 0 | verticalDatum = verticalCSType; |
547 | 0 | verticalCSType = verticalDatum + 600; |
548 | 0 | } |
549 | | |
550 | | /* -------------------------------------------------------------------- |
551 | | */ |
552 | | /* This addresses another case where the EGM96 Vertical Datum |
553 | | * code */ |
554 | | /* is misused as a Vertical CS code (#4922). */ |
555 | | /* -------------------------------------------------------------------- |
556 | | */ |
557 | 0 | if (verticalCSType == 5171) |
558 | 0 | { |
559 | 0 | verticalDatum = 5171; |
560 | 0 | verticalCSType = 5773; |
561 | 0 | } |
562 | 0 | } |
563 | | |
564 | | /* -------------------------------------------------------------------- |
565 | | */ |
566 | | /* Somewhat similarly, codes 5001 to 5033 were treated as */ |
567 | | /* vertical coordinate systems based on ellipsoidal heights. */ |
568 | | /* We use the corresponding geodetic datum as the vertical */ |
569 | | /* datum and clear the vertical coordinate system code since */ |
570 | | /* there isn't one in EPSG. */ |
571 | | /* -------------------------------------------------------------------- |
572 | | */ |
573 | 0 | if ((verticalCSType >= 5001 && verticalCSType <= 5033) && |
574 | 0 | verticalDatum == 0) |
575 | 0 | { |
576 | 0 | verticalDatum = verticalCSType + 1000; |
577 | 0 | verticalCSType = 0; |
578 | 0 | } |
579 | 0 | } |
580 | | |
581 | | /* -------------------------------------------------------------------- */ |
582 | | /* Handle non-standard coordinate systems as LOCAL_CS. */ |
583 | | /* -------------------------------------------------------------------- */ |
584 | 0 | if (psDefn->Model != ModelTypeProjected && |
585 | 0 | psDefn->Model != ModelTypeGeographic && |
586 | 0 | psDefn->Model != ModelTypeGeocentric) |
587 | 0 | { |
588 | 0 | char szPeStr[2400] = {'\0'}; |
589 | | |
590 | | /** check if there is a pe string citation key **/ |
591 | 0 | if (GDALGTIFKeyGetASCII(hGTIF, PCSCitationGeoKey, szPeStr, |
592 | 0 | sizeof(szPeStr)) && |
593 | 0 | strstr(szPeStr, "ESRI PE String = ")) |
594 | 0 | { |
595 | 0 | const char *pszWKT = szPeStr + strlen("ESRI PE String = "); |
596 | 0 | oSRS.importFromWkt(pszWKT); |
597 | |
|
598 | 0 | if (strstr(pszWKT, |
599 | 0 | "PROJCS[\"WGS_1984_Web_Mercator_Auxiliary_Sphere\"")) |
600 | 0 | { |
601 | 0 | oSRS.SetExtension( |
602 | 0 | "PROJCS", "PROJ4", |
603 | 0 | "+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 " |
604 | 0 | "+x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null " |
605 | 0 | "+wktext +no_defs"); // TODO(schwehr): Why 2 spaces? |
606 | 0 | } |
607 | |
|
608 | 0 | return OGRSpatialReference::ToHandle(oSRS.Clone()); |
609 | 0 | } |
610 | 0 | else |
611 | 0 | { |
612 | 0 | char *pszUnitsName = nullptr; |
613 | 0 | char szPCSName[300] = {'\0'}; |
614 | 0 | int nKeyCount = 0; |
615 | 0 | int anVersion[3] = {0}; |
616 | |
|
617 | 0 | GTIFDirectoryInfo(hGTIF, anVersion, &nKeyCount); |
618 | |
|
619 | 0 | if (nKeyCount > 0) // Use LOCAL_CS if we have any geokeys at all. |
620 | 0 | { |
621 | | // Handle citation. |
622 | 0 | strcpy(szPCSName, "unnamed"); |
623 | 0 | if (!GDALGTIFKeyGetASCII(hGTIF, GTCitationGeoKey, szPCSName, |
624 | 0 | sizeof(szPCSName))) |
625 | 0 | GDALGTIFKeyGetASCII(hGTIF, GeogCitationGeoKey, szPCSName, |
626 | 0 | sizeof(szPCSName)); |
627 | |
|
628 | 0 | GTIFCleanupImagineNames(szPCSName); |
629 | 0 | oSRS.SetLocalCS(szPCSName); |
630 | | |
631 | | // Handle units |
632 | 0 | if (psDefn->UOMLength != KvUserDefined) |
633 | 0 | { |
634 | 0 | #if LIBGEOTIFF_VERSION >= 1600 |
635 | 0 | GTIFGetUOMLengthInfoEx(projContext, |
636 | | #else |
637 | | GTIFGetUOMLengthInfo( |
638 | | #endif |
639 | 0 | psDefn->UOMLength, &pszUnitsName, |
640 | 0 | nullptr); |
641 | 0 | } |
642 | |
|
643 | 0 | if (pszUnitsName != nullptr) |
644 | 0 | { |
645 | 0 | char szUOMLength[12]; |
646 | 0 | snprintf(szUOMLength, sizeof(szUOMLength), "%d", |
647 | 0 | psDefn->UOMLength); |
648 | 0 | oSRS.SetTargetLinearUnits(nullptr, pszUnitsName, |
649 | 0 | psDefn->UOMLengthInMeters, "EPSG", |
650 | 0 | szUOMLength); |
651 | 0 | } |
652 | 0 | else |
653 | 0 | oSRS.SetLinearUnits("unknown", psDefn->UOMLengthInMeters); |
654 | |
|
655 | 0 | if (verticalUnits != 0) |
656 | 0 | { |
657 | 0 | char szVertCSCitation[2048] = {0}; |
658 | 0 | if (GDALGTIFKeyGetASCII(hGTIF, VerticalCitationGeoKey, |
659 | 0 | szVertCSCitation, |
660 | 0 | sizeof(szVertCSCitation))) |
661 | 0 | { |
662 | 0 | if (STARTS_WITH_CI(szVertCSCitation, "VCS Name = ")) |
663 | 0 | { |
664 | 0 | memmove(szVertCSCitation, |
665 | 0 | szVertCSCitation + strlen("VCS Name = "), |
666 | 0 | strlen(szVertCSCitation + |
667 | 0 | strlen("VCS Name = ")) + |
668 | 0 | 1); |
669 | 0 | char *pszPipeChar = strchr(szVertCSCitation, '|'); |
670 | 0 | if (pszPipeChar) |
671 | 0 | *pszPipeChar = '\0'; |
672 | 0 | } |
673 | 0 | } |
674 | 0 | else |
675 | 0 | { |
676 | 0 | strcpy(szVertCSCitation, "unknown"); |
677 | 0 | } |
678 | |
|
679 | 0 | const char *pszHorizontalName = oSRS.GetName(); |
680 | 0 | const std::string osHorizontalName( |
681 | 0 | pszHorizontalName ? pszHorizontalName : "unnamed"); |
682 | | /* -------------------------------------------------------------------- |
683 | | */ |
684 | | /* Promote to being a compound coordinate system. */ |
685 | | /* -------------------------------------------------------------------- |
686 | | */ |
687 | 0 | OGR_SRSNode *poOldRoot = oSRS.GetRoot()->Clone(); |
688 | |
|
689 | 0 | oSRS.Clear(); |
690 | | |
691 | | /* -------------------------------------------------------------------- |
692 | | */ |
693 | | /* Set COMPD_CS name. */ |
694 | | /* -------------------------------------------------------------------- |
695 | | */ |
696 | 0 | char szCTString[512]; |
697 | 0 | szCTString[0] = '\0'; |
698 | 0 | if (GDALGTIFKeyGetASCII(hGTIF, GTCitationGeoKey, szCTString, |
699 | 0 | sizeof(szCTString)) && |
700 | 0 | strstr(szCTString, " = ") == nullptr) |
701 | 0 | { |
702 | 0 | oSRS.SetNode("COMPD_CS", szCTString); |
703 | 0 | } |
704 | 0 | else |
705 | 0 | { |
706 | 0 | oSRS.SetNode("COMPD_CS", (osHorizontalName + " + " + |
707 | 0 | szVertCSCitation) |
708 | 0 | .c_str()); |
709 | 0 | } |
710 | |
|
711 | 0 | oSRS.GetRoot()->AddChild(poOldRoot); |
712 | |
|
713 | 0 | FillCompoundCRSWithManualVertCS( |
714 | 0 | hGTIF, oSRS, szVertCSCitation, verticalDatum, |
715 | 0 | verticalUnits); |
716 | 0 | } |
717 | |
|
718 | 0 | GTIFFreeMemory(pszUnitsName); |
719 | 0 | } |
720 | 0 | return OGRSpatialReference::ToHandle(oSRS.Clone()); |
721 | 0 | } |
722 | 0 | } |
723 | | |
724 | | /* -------------------------------------------------------------------- */ |
725 | | /* Handle Geocentric coordinate systems. */ |
726 | | /* -------------------------------------------------------------------- */ |
727 | 0 | if (psDefn->Model == ModelTypeGeocentric) |
728 | 0 | { |
729 | 0 | char szName[300] = {'\0'}; |
730 | |
|
731 | 0 | strcpy(szName, "unnamed"); |
732 | 0 | if (!GDALGTIFKeyGetASCII(hGTIF, GTCitationGeoKey, szName, |
733 | 0 | sizeof(szName))) |
734 | 0 | GDALGTIFKeyGetASCII(hGTIF, GeogCitationGeoKey, szName, |
735 | 0 | sizeof(szName)); |
736 | |
|
737 | 0 | oSRS.SetGeocCS(szName); |
738 | |
|
739 | 0 | char *pszUnitsName = nullptr; |
740 | |
|
741 | 0 | if (psDefn->UOMLength != KvUserDefined) |
742 | 0 | { |
743 | 0 | #if LIBGEOTIFF_VERSION >= 1600 |
744 | 0 | GTIFGetUOMLengthInfoEx(projContext, |
745 | | #else |
746 | | GTIFGetUOMLengthInfo( |
747 | | #endif |
748 | 0 | psDefn->UOMLength, &pszUnitsName, nullptr); |
749 | 0 | } |
750 | |
|
751 | 0 | if (pszUnitsName != nullptr) |
752 | 0 | { |
753 | 0 | char szUOMLength[12]; |
754 | 0 | snprintf(szUOMLength, sizeof(szUOMLength), "%d", psDefn->UOMLength); |
755 | 0 | oSRS.SetTargetLinearUnits(nullptr, pszUnitsName, |
756 | 0 | psDefn->UOMLengthInMeters, "EPSG", |
757 | 0 | szUOMLength); |
758 | 0 | } |
759 | 0 | else |
760 | 0 | oSRS.SetLinearUnits("unknown", psDefn->UOMLengthInMeters); |
761 | |
|
762 | 0 | GTIFFreeMemory(pszUnitsName); |
763 | 0 | } |
764 | | |
765 | | /* -------------------------------------------------------------------- */ |
766 | | /* #3901: In libgeotiff 1.3.0 and earlier we incorrectly */ |
767 | | /* interpreted linear projection parameter geokeys (false */ |
768 | | /* easting/northing) as being in meters instead of the */ |
769 | | /* coordinate system of the file. The following code attempts */ |
770 | | /* to provide mechanisms for fixing the issue if we are linked */ |
771 | | /* with an older version of libgeotiff. */ |
772 | | /* -------------------------------------------------------------------- */ |
773 | 0 | const char *pszLinearUnits = |
774 | 0 | CPLGetConfigOption("GTIFF_LINEAR_UNITS", "DEFAULT"); |
775 | | |
776 | | /* -------------------------------------------------------------------- */ |
777 | | /* #3901: If folks have broken GeoTIFF files generated with */ |
778 | | /* older versions of GDAL+libgeotiff, then they may need a */ |
779 | | /* hack to allow them to be read properly. This is that */ |
780 | | /* hack. We basically try to undue the conversion applied by */ |
781 | | /* libgeotiff to meters (or above) to simulate the old */ |
782 | | /* behavior. */ |
783 | | /* -------------------------------------------------------------------- */ |
784 | 0 | unsigned short bLinearUnitsMarkedCorrect = FALSE; |
785 | |
|
786 | 0 | GDALGTIFKeyGetSHORT(hGTIF, ProjLinearUnitsInterpCorrectGeoKey, |
787 | 0 | &bLinearUnitsMarkedCorrect, 0, 1); |
788 | |
|
789 | 0 | if (EQUAL(pszLinearUnits, "BROKEN") && |
790 | 0 | psDefn->Projection == KvUserDefined && !bLinearUnitsMarkedCorrect) |
791 | 0 | { |
792 | 0 | for (int iParam = 0; iParam < psDefn->nParms; iParam++) |
793 | 0 | { |
794 | 0 | switch (psDefn->ProjParmId[iParam]) |
795 | 0 | { |
796 | 0 | case ProjFalseEastingGeoKey: |
797 | 0 | case ProjFalseNorthingGeoKey: |
798 | 0 | case ProjFalseOriginEastingGeoKey: |
799 | 0 | case ProjFalseOriginNorthingGeoKey: |
800 | 0 | case ProjCenterEastingGeoKey: |
801 | 0 | case ProjCenterNorthingGeoKey: |
802 | 0 | if (psDefn->UOMLengthInMeters != 0 && |
803 | 0 | psDefn->UOMLengthInMeters != 1.0) |
804 | 0 | { |
805 | 0 | psDefn->ProjParm[iParam] /= psDefn->UOMLengthInMeters; |
806 | 0 | CPLDebug( |
807 | 0 | "GTIFF", |
808 | 0 | "Converting geokey to accommodate old broken file " |
809 | 0 | "due to GTIFF_LINEAR_UNITS=BROKEN setting."); |
810 | 0 | } |
811 | 0 | break; |
812 | | |
813 | 0 | default: |
814 | 0 | break; |
815 | 0 | } |
816 | 0 | } |
817 | 0 | } |
818 | | |
819 | | /* -------------------------------------------------------------------- */ |
820 | | /* If this is a projected SRS we set the PROJCS keyword first */ |
821 | | /* to ensure that the GEOGCS will be a child. */ |
822 | | /* -------------------------------------------------------------------- */ |
823 | 0 | OGRBoolean linearUnitIsSet = FALSE; |
824 | 0 | if (psDefn->Model == ModelTypeProjected) |
825 | 0 | { |
826 | 0 | char szCTString[512] = {'\0'}; |
827 | 0 | if (psDefn->PCS != KvUserDefined) |
828 | 0 | { |
829 | 0 | char *pszPCSName = nullptr; |
830 | |
|
831 | 0 | #if LIBGEOTIFF_VERSION >= 1600 |
832 | 0 | GTIFGetPCSInfoEx(projContext, |
833 | | #else |
834 | | GTIFGetPCSInfo( |
835 | | #endif |
836 | 0 | psDefn->PCS, &pszPCSName, nullptr, nullptr, |
837 | 0 | nullptr); |
838 | |
|
839 | 0 | oSRS.SetProjCS(pszPCSName ? pszPCSName : "unnamed"); |
840 | 0 | if (pszPCSName) |
841 | 0 | GTIFFreeMemory(pszPCSName); |
842 | |
|
843 | 0 | oSRS.SetLinearUnits("unknown", 1.0); |
844 | 0 | } |
845 | 0 | else |
846 | 0 | { |
847 | 0 | bool bTryGTCitationGeoKey = true; |
848 | 0 | if (GDALGTIFKeyGetASCII(hGTIF, PCSCitationGeoKey, szCTString, |
849 | 0 | sizeof(szCTString))) |
850 | 0 | { |
851 | 0 | bTryGTCitationGeoKey = false; |
852 | 0 | if (!SetCitationToSRS(hGTIF, szCTString, sizeof(szCTString), |
853 | 0 | PCSCitationGeoKey, &oSRS, |
854 | 0 | &linearUnitIsSet)) |
855 | 0 | { |
856 | 0 | if (!STARTS_WITH_CI(szCTString, "LUnits = ")) |
857 | 0 | { |
858 | 0 | oSRS.SetProjCS(szCTString); |
859 | 0 | oSRS.SetLinearUnits("unknown", 1.0); |
860 | 0 | } |
861 | 0 | else |
862 | 0 | { |
863 | 0 | bTryGTCitationGeoKey = true; |
864 | 0 | } |
865 | 0 | } |
866 | 0 | } |
867 | |
|
868 | 0 | if (bTryGTCitationGeoKey) |
869 | 0 | { |
870 | 0 | if (GDALGTIFKeyGetASCII(hGTIF, GTCitationGeoKey, szCTString, |
871 | 0 | sizeof(szCTString)) && |
872 | 0 | !SetCitationToSRS(hGTIF, szCTString, sizeof(szCTString), |
873 | 0 | GTCitationGeoKey, &oSRS, |
874 | 0 | &linearUnitIsSet)) |
875 | 0 | { |
876 | 0 | oSRS.SetNode("PROJCS", szCTString); |
877 | 0 | oSRS.SetLinearUnits("unknown", 1.0); |
878 | 0 | } |
879 | 0 | else |
880 | 0 | { |
881 | 0 | oSRS.SetNode("PROJCS", "unnamed"); |
882 | 0 | oSRS.SetLinearUnits("unknown", 1.0); |
883 | 0 | } |
884 | 0 | } |
885 | 0 | } |
886 | | |
887 | | /* Handle ESRI/Erdas style state plane and UTM in citation key */ |
888 | 0 | if (CheckCitationKeyForStatePlaneUTM(hGTIF, psDefn, &oSRS, |
889 | 0 | &linearUnitIsSet)) |
890 | 0 | { |
891 | 0 | return OGRSpatialReference::ToHandle(oSRS.Clone()); |
892 | 0 | } |
893 | | |
894 | | /* Handle ESRI PE string in citation */ |
895 | 0 | szCTString[0] = '\0'; |
896 | 0 | if (GDALGTIFKeyGetASCII(hGTIF, GTCitationGeoKey, szCTString, |
897 | 0 | sizeof(szCTString))) |
898 | 0 | SetCitationToSRS(hGTIF, szCTString, sizeof(szCTString), |
899 | 0 | GTCitationGeoKey, &oSRS, &linearUnitIsSet); |
900 | 0 | } |
901 | | |
902 | | /* ==================================================================== */ |
903 | | /* Setup the GeogCS */ |
904 | | /* ==================================================================== */ |
905 | 0 | char *pszGeogName = nullptr; |
906 | 0 | char *pszDatumName = nullptr; |
907 | 0 | char *pszPMName = nullptr; |
908 | 0 | char *pszSpheroidName = nullptr; |
909 | 0 | char *pszAngularUnits = nullptr; |
910 | 0 | char szGCSName[512] = {'\0'}; |
911 | |
|
912 | 0 | if (! |
913 | 0 | #if LIBGEOTIFF_VERSION >= 1600 |
914 | 0 | GTIFGetGCSInfoEx(projContext, |
915 | | #else |
916 | | GTIFGetGCSInfo( |
917 | | #endif |
918 | 0 | psDefn->GCS, &pszGeogName, nullptr, nullptr, |
919 | 0 | nullptr) && |
920 | 0 | GDALGTIFKeyGetASCII(hGTIF, GeogCitationGeoKey, szGCSName, |
921 | 0 | sizeof(szGCSName))) |
922 | 0 | { |
923 | 0 | GetGeogCSFromCitation(szGCSName, sizeof(szGCSName), GeogCitationGeoKey, |
924 | 0 | &pszGeogName, &pszDatumName, &pszPMName, |
925 | 0 | &pszSpheroidName, &pszAngularUnits); |
926 | 0 | } |
927 | 0 | else |
928 | 0 | { |
929 | 0 | GTIFToCPLRecycleString(&pszGeogName); |
930 | 0 | } |
931 | |
|
932 | 0 | if (pszGeogName && STARTS_WITH(pszGeogName, "GCS_")) |
933 | 0 | { |
934 | | // Morph from potential ESRI name |
935 | 0 | char *pszOfficialName = GTIFGetEPSGOfficialName( |
936 | 0 | hGTIF, PJ_TYPE_GEOGRAPHIC_2D_CRS, pszGeogName); |
937 | 0 | if (pszOfficialName) |
938 | 0 | { |
939 | 0 | CPLFree(pszGeogName); |
940 | 0 | pszGeogName = pszOfficialName; |
941 | 0 | } |
942 | 0 | } |
943 | |
|
944 | 0 | if (pszDatumName && strchr(pszDatumName, '_')) |
945 | 0 | { |
946 | | // Morph from potential ESRI name |
947 | 0 | char *pszOfficialName = GTIFGetEPSGOfficialName( |
948 | 0 | hGTIF, PJ_TYPE_GEODETIC_REFERENCE_FRAME, pszDatumName); |
949 | 0 | if (pszOfficialName) |
950 | 0 | { |
951 | 0 | CPLFree(pszDatumName); |
952 | 0 | pszDatumName = pszOfficialName; |
953 | 0 | } |
954 | 0 | } |
955 | |
|
956 | 0 | if (pszSpheroidName && strchr(pszSpheroidName, '_')) |
957 | 0 | { |
958 | | // Morph from potential ESRI name |
959 | 0 | char *pszOfficialName = |
960 | 0 | GTIFGetEPSGOfficialName(hGTIF, PJ_TYPE_ELLIPSOID, pszSpheroidName); |
961 | 0 | if (pszOfficialName) |
962 | 0 | { |
963 | 0 | CPLFree(pszSpheroidName); |
964 | 0 | pszSpheroidName = pszOfficialName; |
965 | 0 | } |
966 | 0 | } |
967 | |
|
968 | 0 | if (!pszDatumName) |
969 | 0 | { |
970 | 0 | #if LIBGEOTIFF_VERSION >= 1600 |
971 | 0 | GTIFGetDatumInfoEx(projContext, |
972 | | #else |
973 | | GTIFGetDatumInfo( |
974 | | #endif |
975 | 0 | psDefn->Datum, &pszDatumName, nullptr); |
976 | 0 | GTIFToCPLRecycleString(&pszDatumName); |
977 | 0 | } |
978 | |
|
979 | 0 | double dfSemiMajor = 0.0; |
980 | 0 | double dfInvFlattening = 0.0; |
981 | 0 | if (!pszSpheroidName) |
982 | 0 | { |
983 | 0 | #if LIBGEOTIFF_VERSION >= 1600 |
984 | 0 | GTIFGetEllipsoidInfoEx(projContext, |
985 | | #else |
986 | | GTIFGetEllipsoidInfo( |
987 | | #endif |
988 | 0 | psDefn->Ellipsoid, &pszSpheroidName, nullptr, |
989 | 0 | nullptr); |
990 | 0 | GTIFToCPLRecycleString(&pszSpheroidName); |
991 | 0 | } |
992 | 0 | else |
993 | 0 | { |
994 | 0 | CPL_IGNORE_RET_VAL(GDALGTIFKeyGetDOUBLE(hGTIF, GeogSemiMajorAxisGeoKey, |
995 | 0 | &(psDefn->SemiMajor), 0, 1)); |
996 | 0 | CPL_IGNORE_RET_VAL(GDALGTIFKeyGetDOUBLE(hGTIF, GeogInvFlatteningGeoKey, |
997 | 0 | &dfInvFlattening, 0, 1)); |
998 | 0 | if (std::isinf(dfInvFlattening)) |
999 | 0 | { |
1000 | | // Deal with the non-nominal case of |
1001 | | // https://github.com/OSGeo/PROJ/issues/2317 |
1002 | 0 | dfInvFlattening = 0; |
1003 | 0 | } |
1004 | 0 | } |
1005 | 0 | if (!pszPMName) |
1006 | 0 | { |
1007 | 0 | #if LIBGEOTIFF_VERSION >= 1600 |
1008 | 0 | GTIFGetPMInfoEx(projContext, |
1009 | | #else |
1010 | | GTIFGetPMInfo( |
1011 | | #endif |
1012 | 0 | psDefn->PM, &pszPMName, nullptr); |
1013 | 0 | GTIFToCPLRecycleString(&pszPMName); |
1014 | 0 | } |
1015 | 0 | else |
1016 | 0 | { |
1017 | 0 | CPL_IGNORE_RET_VAL( |
1018 | 0 | GDALGTIFKeyGetDOUBLE(hGTIF, GeogPrimeMeridianLongGeoKey, |
1019 | 0 | &(psDefn->PMLongToGreenwich), 0, 1)); |
1020 | 0 | } |
1021 | |
|
1022 | 0 | if (!pszAngularUnits) |
1023 | 0 | { |
1024 | 0 | #if LIBGEOTIFF_VERSION >= 1600 |
1025 | 0 | GTIFGetUOMAngleInfoEx(projContext, |
1026 | | #else |
1027 | | GTIFGetUOMAngleInfo( |
1028 | | #endif |
1029 | 0 | psDefn->UOMAngle, &pszAngularUnits, |
1030 | 0 | &psDefn->UOMAngleInDegrees); |
1031 | 0 | if (pszAngularUnits == nullptr) |
1032 | 0 | pszAngularUnits = CPLStrdup("unknown"); |
1033 | 0 | else |
1034 | 0 | GTIFToCPLRecycleString(&pszAngularUnits); |
1035 | 0 | } |
1036 | 0 | else |
1037 | 0 | { |
1038 | 0 | double dfRadians = 0.0; |
1039 | 0 | if (GDALGTIFKeyGetDOUBLE(hGTIF, GeogAngularUnitSizeGeoKey, &dfRadians, |
1040 | 0 | 0, 1)) |
1041 | 0 | { |
1042 | 0 | psDefn->UOMAngleInDegrees = dfRadians / CPLAtof(SRS_UA_DEGREE_CONV); |
1043 | 0 | } |
1044 | 0 | } |
1045 | | |
1046 | | // Avoid later division by zero. |
1047 | 0 | if (psDefn->UOMAngleInDegrees == 0) |
1048 | 0 | { |
1049 | 0 | CPLError(CE_Warning, CPLE_AppDefined, |
1050 | 0 | "Invalid value for GeogAngularUnitSizeGeoKey."); |
1051 | 0 | psDefn->UOMAngleInDegrees = 1; |
1052 | 0 | } |
1053 | |
|
1054 | 0 | dfSemiMajor = psDefn->SemiMajor; |
1055 | 0 | if (dfSemiMajor == 0.0) |
1056 | 0 | { |
1057 | 0 | CPLFree(pszSpheroidName); |
1058 | 0 | pszSpheroidName = CPLStrdup("unretrievable - using WGS84"); |
1059 | 0 | dfSemiMajor = SRS_WGS84_SEMIMAJOR; |
1060 | 0 | dfInvFlattening = SRS_WGS84_INVFLATTENING; |
1061 | 0 | } |
1062 | 0 | else if (dfInvFlattening == 0.0 && |
1063 | 0 | ((psDefn->SemiMinor / psDefn->SemiMajor) < 0.99999999999999999 || |
1064 | 0 | (psDefn->SemiMinor / psDefn->SemiMajor) > 1.00000000000000001)) |
1065 | 0 | { |
1066 | 0 | dfInvFlattening = |
1067 | 0 | OSRCalcInvFlattening(psDefn->SemiMajor, psDefn->SemiMinor); |
1068 | | |
1069 | | /* Take official inverse flattening definition in the WGS84 case */ |
1070 | 0 | if (fabs(dfSemiMajor - SRS_WGS84_SEMIMAJOR) < 1e-10 && |
1071 | 0 | fabs(dfInvFlattening - SRS_WGS84_INVFLATTENING) < 1e-10) |
1072 | 0 | dfInvFlattening = SRS_WGS84_INVFLATTENING; |
1073 | 0 | } |
1074 | 0 | if (!pszGeogName || strlen(pszGeogName) == 0) |
1075 | 0 | { |
1076 | 0 | CPLFree(pszGeogName); |
1077 | 0 | pszGeogName = CPLStrdup(pszDatumName ? pszDatumName : "unknown"); |
1078 | 0 | } |
1079 | |
|
1080 | 0 | oSRS.SetGeogCS(pszGeogName, pszDatumName, pszSpheroidName, dfSemiMajor, |
1081 | 0 | dfInvFlattening, pszPMName, psDefn->PMLongToGreenwich, |
1082 | 0 | pszAngularUnits, |
1083 | 0 | psDefn->UOMAngleInDegrees * CPLAtof(SRS_UA_DEGREE_CONV)); |
1084 | |
|
1085 | 0 | bool bGeog3DCRS = false; |
1086 | 0 | bool bSetDatumEllipsoidCode = true; |
1087 | 0 | bool bHasWarnedInconsistentGeogCRSEPSG = false; |
1088 | 0 | { |
1089 | 0 | const int nGCS = psDefn->GCS; |
1090 | 0 | if (nGCS != KvUserDefined && nGCS > 0 && |
1091 | 0 | psDefn->Model != ModelTypeGeocentric) |
1092 | 0 | { |
1093 | 0 | OGRSpatialReference oSRSGeog; |
1094 | 0 | const bool bGCSCodeValid = |
1095 | 0 | oSRSGeog.importFromEPSG(nGCS) == OGRERR_NONE; |
1096 | |
|
1097 | 0 | const std::string osGTiffSRSSource = |
1098 | 0 | CPLGetConfigOption("GTIFF_SRS_SOURCE", ""); |
1099 | | |
1100 | | // GeoTIFF 1.0 might put a Geographic 3D code in GeodeticCRSGeoKey |
1101 | 0 | bool bTryCompareToEPSG = oSRSGeog.GetAxesCount() == 2; |
1102 | |
|
1103 | 0 | if (psDefn->Datum != KvUserDefined) |
1104 | 0 | { |
1105 | 0 | char szCode[12]; |
1106 | 0 | snprintf(szCode, sizeof(szCode), "%d", psDefn->Datum); |
1107 | 0 | auto ctx = static_cast<PJ_CONTEXT *>( |
1108 | 0 | GTIFGetPROJContext(hGTIF, true, nullptr)); |
1109 | 0 | auto datum = proj_create_from_database( |
1110 | 0 | ctx, "EPSG", szCode, PJ_CATEGORY_DATUM, 0, nullptr); |
1111 | 0 | if (datum) |
1112 | 0 | { |
1113 | 0 | if (proj_get_type(datum) == |
1114 | 0 | PJ_TYPE_DYNAMIC_GEODETIC_REFERENCE_FRAME) |
1115 | 0 | { |
1116 | | // Current PROJ versions will not manage to |
1117 | | // consider a CRS with a regular datum and another one |
1118 | | // with a dynamic datum as being equivalent. |
1119 | 0 | bTryCompareToEPSG = false; |
1120 | 0 | } |
1121 | 0 | proj_destroy(datum); |
1122 | 0 | } |
1123 | 0 | } |
1124 | |
|
1125 | 0 | if (bTryCompareToEPSG && !oSRSGeog.IsSameGeogCS(&oSRS) && |
1126 | 0 | osGTiffSRSSource.empty()) |
1127 | 0 | { |
1128 | | // See https://github.com/OSGeo/gdal/issues/5399 |
1129 | | // where a file has inconsistent GeogSemiMinorAxisGeoKey / |
1130 | | // GeogInvFlatteningGeoKey values, which cause its datum to be |
1131 | | // considered as non-equivalent to the EPSG one. |
1132 | 0 | CPLError( |
1133 | 0 | CE_Warning, CPLE_AppDefined, |
1134 | 0 | "The definition of geographic CRS EPSG:%d got from GeoTIFF " |
1135 | 0 | "keys " |
1136 | 0 | "is not the same as the one from the EPSG registry, " |
1137 | 0 | "which may cause issues during reprojection operations. " |
1138 | 0 | "Set GTIFF_SRS_SOURCE configuration option to EPSG to " |
1139 | 0 | "use official parameters (overriding the ones from GeoTIFF " |
1140 | 0 | "keys), " |
1141 | 0 | "or to GEOKEYS to use custom values from GeoTIFF keys " |
1142 | 0 | "and drop the EPSG code.", |
1143 | 0 | nGCS); |
1144 | 0 | bHasWarnedInconsistentGeogCRSEPSG = true; |
1145 | 0 | } |
1146 | 0 | if (EQUAL(osGTiffSRSSource.c_str(), "EPSG")) |
1147 | 0 | { |
1148 | 0 | oSRS.CopyGeogCSFrom(&oSRSGeog); |
1149 | 0 | } |
1150 | 0 | else if (osGTiffSRSSource.empty() && oSRSGeog.IsDynamic() && |
1151 | 0 | psDefn->Model == ModelTypeGeographic) |
1152 | 0 | { |
1153 | | // We should perhaps be more careful and detect overrides |
1154 | | // of geokeys... |
1155 | 0 | oSRS = oSRSGeog; |
1156 | 0 | bSetDatumEllipsoidCode = false; |
1157 | 0 | } |
1158 | 0 | else if (bGCSCodeValid && osGTiffSRSSource.empty()) |
1159 | 0 | { |
1160 | 0 | oSRS.SetAuthority("GEOGCS", "EPSG", nGCS); |
1161 | 0 | } |
1162 | 0 | else |
1163 | 0 | { |
1164 | 0 | bSetDatumEllipsoidCode = false; |
1165 | 0 | } |
1166 | |
|
1167 | 0 | int nVertSRSCode = verticalCSType; |
1168 | 0 | if (verticalDatum == 6030 && nGCS == 4326) // DatumE_WGS84 |
1169 | 0 | { |
1170 | 0 | nVertSRSCode = 4979; |
1171 | 0 | } |
1172 | | |
1173 | | // Try to reconstruct a Geographic3D CRS from the |
1174 | | // GeodeticCRSGeoKey and the VerticalGeoKey, when they are |
1175 | | // consistent |
1176 | 0 | if (nVertSRSCode > 0 && nVertSRSCode != KvUserDefined) |
1177 | 0 | { |
1178 | 0 | OGRSpatialReference oTmpVertSRS; |
1179 | 0 | if (oSRSGeog.IsGeographic() && oSRSGeog.GetAxesCount() == 2 && |
1180 | 0 | oTmpVertSRS.importFromEPSG(nVertSRSCode) == OGRERR_NONE && |
1181 | 0 | oTmpVertSRS.IsGeographic() && |
1182 | 0 | oTmpVertSRS.GetAxesCount() == 3) |
1183 | 0 | { |
1184 | 0 | const char *pszTmpCode = |
1185 | 0 | oSRSGeog.GetAuthorityCode("GEOGCS|DATUM"); |
1186 | 0 | const char *pszTmpVertCode = |
1187 | 0 | oTmpVertSRS.GetAuthorityCode("GEOGCS|DATUM"); |
1188 | 0 | if (pszTmpCode && pszTmpVertCode && |
1189 | 0 | atoi(pszTmpCode) == atoi(pszTmpVertCode)) |
1190 | 0 | { |
1191 | 0 | verticalCSType = 0; |
1192 | 0 | verticalDatum = 0; |
1193 | 0 | verticalUnits = 0; |
1194 | 0 | oSRS.CopyGeogCSFrom(&oTmpVertSRS); |
1195 | 0 | bSetDatumEllipsoidCode = false; |
1196 | 0 | bGeog3DCRS = true; |
1197 | 0 | } |
1198 | 0 | } |
1199 | 0 | } |
1200 | 0 | } |
1201 | 0 | } |
1202 | 0 | if (bSetDatumEllipsoidCode) |
1203 | 0 | { |
1204 | 0 | if (psDefn->Datum != KvUserDefined) |
1205 | 0 | oSRS.SetAuthority("DATUM", "EPSG", psDefn->Datum); |
1206 | |
|
1207 | 0 | if (psDefn->Ellipsoid != KvUserDefined) |
1208 | 0 | oSRS.SetAuthority("SPHEROID", "EPSG", psDefn->Ellipsoid); |
1209 | 0 | } |
1210 | |
|
1211 | 0 | CPLFree(pszGeogName); |
1212 | 0 | CPLFree(pszDatumName); |
1213 | 0 | CPLFree(pszSpheroidName); |
1214 | 0 | CPLFree(pszPMName); |
1215 | 0 | CPLFree(pszAngularUnits); |
1216 | | |
1217 | | /* -------------------------------------------------------------------- */ |
1218 | | /* Set projection units if not yet done */ |
1219 | | /* -------------------------------------------------------------------- */ |
1220 | 0 | if (psDefn->Model == ModelTypeProjected && !linearUnitIsSet) |
1221 | 0 | { |
1222 | 0 | char *pszUnitsName = nullptr; |
1223 | |
|
1224 | 0 | if (psDefn->UOMLength != KvUserDefined) |
1225 | 0 | { |
1226 | 0 | #if LIBGEOTIFF_VERSION >= 1600 |
1227 | 0 | GTIFGetUOMLengthInfoEx(projContext, |
1228 | | #else |
1229 | | GTIFGetUOMLengthInfo( |
1230 | | #endif |
1231 | 0 | psDefn->UOMLength, &pszUnitsName, nullptr); |
1232 | 0 | } |
1233 | |
|
1234 | 0 | if (pszUnitsName != nullptr) |
1235 | 0 | { |
1236 | 0 | char szUOMLength[12]; |
1237 | 0 | snprintf(szUOMLength, sizeof(szUOMLength), "%d", psDefn->UOMLength); |
1238 | 0 | oSRS.SetTargetLinearUnits(nullptr, pszUnitsName, |
1239 | 0 | psDefn->UOMLengthInMeters, "EPSG", |
1240 | 0 | szUOMLength); |
1241 | 0 | } |
1242 | 0 | else |
1243 | 0 | oSRS.SetLinearUnits("unknown", psDefn->UOMLengthInMeters); |
1244 | |
|
1245 | 0 | GTIFFreeMemory(pszUnitsName); |
1246 | 0 | } |
1247 | | |
1248 | | /* ==================================================================== */ |
1249 | | /* Try to import PROJCS from ProjectedCSTypeGeoKey if we */ |
1250 | | /* have essentially only it. We could relax a bit the constraints */ |
1251 | | /* but that should do for now. This may mask shortcomings in the */ |
1252 | | /* libgeotiff GTIFGetDefn() function. */ |
1253 | | /* ==================================================================== */ |
1254 | 0 | unsigned short tmp = 0; |
1255 | 0 | bool bGotFromEPSG = false; |
1256 | 0 | if (psDefn->Model == ModelTypeProjected && psDefn->PCS != KvUserDefined && |
1257 | 0 | GDALGTIFKeyGetSHORT(hGTIF, ProjectionGeoKey, &tmp, 0, 1) == 0 && |
1258 | 0 | GDALGTIFKeyGetSHORT(hGTIF, ProjCoordTransGeoKey, &tmp, 0, 1) == 0 && |
1259 | 0 | GDALGTIFKeyGetSHORT(hGTIF, GeographicTypeGeoKey, &tmp, 0, 1) == 0 && |
1260 | 0 | GDALGTIFKeyGetSHORT(hGTIF, GeogGeodeticDatumGeoKey, &tmp, 0, 1) == 0 && |
1261 | 0 | GDALGTIFKeyGetSHORT(hGTIF, GeogEllipsoidGeoKey, &tmp, 0, 1) == 0 && |
1262 | 0 | CPLTestBool(CPLGetConfigOption("GTIFF_IMPORT_FROM_EPSG", "YES"))) |
1263 | 0 | { |
1264 | | // Save error state as importFromEPSGA() will call CPLReset() |
1265 | 0 | CPLErrorNum errNo = CPLGetLastErrorNo(); |
1266 | 0 | CPLErr eErr = CPLGetLastErrorType(); |
1267 | 0 | const char *pszTmp = CPLGetLastErrorMsg(); |
1268 | 0 | char *pszLastErrorMsg = CPLStrdup(pszTmp ? pszTmp : ""); |
1269 | 0 | CPLPushErrorHandler(CPLQuietErrorHandler); |
1270 | 0 | OGRSpatialReference oSRSTmp; |
1271 | 0 | OGRErr eImportErr = oSRSTmp.importFromEPSG(psDefn->PCS); |
1272 | 0 | CPLPopErrorHandler(); |
1273 | | // Restore error state |
1274 | 0 | CPLErrorSetState(eErr, errNo, pszLastErrorMsg); |
1275 | 0 | CPLFree(pszLastErrorMsg); |
1276 | 0 | bGotFromEPSG = eImportErr == OGRERR_NONE; |
1277 | |
|
1278 | 0 | if (bGotFromEPSG) |
1279 | 0 | { |
1280 | | // See #6210. In case there's an overridden linear units, take it |
1281 | | // into account |
1282 | 0 | const char *pszUnitsName = nullptr; |
1283 | 0 | double dfUOMLengthInMeters = oSRS.GetLinearUnits(&pszUnitsName); |
1284 | | // Non exact comparison, as there's a slight difference between |
1285 | | // the evaluation of US Survey foot hardcoded in geo_normalize.c to |
1286 | | // 12.0 / 39.37, and the corresponding value returned by |
1287 | | // PROJ >= 6.0.0 and <= 7.0.0 for EPSG:9003 |
1288 | 0 | if (fabs(dfUOMLengthInMeters - oSRSTmp.GetLinearUnits(nullptr)) > |
1289 | 0 | 1e-15 * dfUOMLengthInMeters) |
1290 | 0 | { |
1291 | 0 | CPLDebug("GTiff", "Modify EPSG:%d to have %s linear units...", |
1292 | 0 | psDefn->PCS, pszUnitsName ? pszUnitsName : "unknown"); |
1293 | |
|
1294 | 0 | const char *pszUnitAuthorityCode = |
1295 | 0 | oSRS.GetAuthorityCode("PROJCS|UNIT"); |
1296 | 0 | const char *pszUnitAuthorityName = |
1297 | 0 | oSRS.GetAuthorityName("PROJCS|UNIT"); |
1298 | |
|
1299 | 0 | if (pszUnitsName) |
1300 | 0 | oSRSTmp.SetLinearUnitsAndUpdateParameters( |
1301 | 0 | pszUnitsName, dfUOMLengthInMeters, pszUnitAuthorityCode, |
1302 | 0 | pszUnitAuthorityName); |
1303 | 0 | } |
1304 | |
|
1305 | 0 | if (bGeog3DCRS) |
1306 | 0 | { |
1307 | 0 | oSRSTmp.CopyGeogCSFrom(&oSRS); |
1308 | 0 | oSRSTmp.UpdateCoordinateSystemFromGeogCRS(); |
1309 | 0 | } |
1310 | 0 | oSRS = std::move(oSRSTmp); |
1311 | 0 | } |
1312 | 0 | } |
1313 | |
|
1314 | 0 | #if !defined(GEO_NORMALIZE_DISABLE_TOWGS84) |
1315 | 0 | if (psDefn->TOWGS84Count > 0 && bGotFromEPSG && |
1316 | 0 | CPLTestBool(CPLGetConfigOption("OSR_STRIP_TOWGS84", "YES"))) |
1317 | 0 | { |
1318 | 0 | CPLDebug("OSR", "TOWGS84 information has been removed. " |
1319 | 0 | "It can be kept by setting the OSR_STRIP_TOWGS84 " |
1320 | 0 | "configuration option to NO"); |
1321 | 0 | } |
1322 | 0 | else if (psDefn->TOWGS84Count > 0 && |
1323 | 0 | (!bGotFromEPSG || |
1324 | 0 | !CPLTestBool(CPLGetConfigOption("OSR_STRIP_TOWGS84", "YES")))) |
1325 | 0 | { |
1326 | 0 | if (bGotFromEPSG) |
1327 | 0 | { |
1328 | 0 | double adfTOWGS84[7] = {0.0}; |
1329 | 0 | CPL_IGNORE_RET_VAL(oSRS.GetTOWGS84(adfTOWGS84)); |
1330 | 0 | bool bSame = true; |
1331 | 0 | for (int i = 0; i < 7; i++) |
1332 | 0 | { |
1333 | 0 | if (fabs(adfTOWGS84[i] - psDefn->TOWGS84[i]) > 1e-5) |
1334 | 0 | { |
1335 | 0 | bSame = false; |
1336 | 0 | break; |
1337 | 0 | } |
1338 | 0 | } |
1339 | 0 | if (!bSame) |
1340 | 0 | { |
1341 | 0 | CPLDebug("GTiff", |
1342 | 0 | "Modify EPSG:%d to have " |
1343 | 0 | "TOWGS84=%f,%f,%f,%f,%f,%f,%f " |
1344 | 0 | "coming from GeogTOWGS84GeoKey, instead of " |
1345 | 0 | "%f,%f,%f,%f,%f,%f,%f coming from EPSG", |
1346 | 0 | psDefn->PCS, psDefn->TOWGS84[0], psDefn->TOWGS84[1], |
1347 | 0 | psDefn->TOWGS84[2], psDefn->TOWGS84[3], |
1348 | 0 | psDefn->TOWGS84[4], psDefn->TOWGS84[5], |
1349 | 0 | psDefn->TOWGS84[6], adfTOWGS84[0], adfTOWGS84[1], |
1350 | 0 | adfTOWGS84[2], adfTOWGS84[3], adfTOWGS84[4], |
1351 | 0 | adfTOWGS84[5], adfTOWGS84[6]); |
1352 | 0 | } |
1353 | 0 | } |
1354 | |
|
1355 | 0 | oSRS.SetTOWGS84(psDefn->TOWGS84[0], psDefn->TOWGS84[1], |
1356 | 0 | psDefn->TOWGS84[2], psDefn->TOWGS84[3], |
1357 | 0 | psDefn->TOWGS84[4], psDefn->TOWGS84[5], |
1358 | 0 | psDefn->TOWGS84[6]); |
1359 | 0 | } |
1360 | 0 | #endif |
1361 | | |
1362 | | /* ==================================================================== */ |
1363 | | /* Handle projection parameters. */ |
1364 | | /* ==================================================================== */ |
1365 | 0 | if (psDefn->Model == ModelTypeProjected && !bGotFromEPSG) |
1366 | 0 | { |
1367 | | /* -------------------------------------------------------------------- |
1368 | | */ |
1369 | | /* Make a local copy of params, and convert back into the */ |
1370 | | /* angular units of the GEOGCS and the linear units of the */ |
1371 | | /* projection. */ |
1372 | | /* -------------------------------------------------------------------- |
1373 | | */ |
1374 | 0 | double adfParam[10] = {0.0}; |
1375 | 0 | int i = 0; // Used after for. |
1376 | |
|
1377 | 0 | for (; i < std::min(10, psDefn->nParms); i++) |
1378 | 0 | adfParam[i] = psDefn->ProjParm[i]; |
1379 | |
|
1380 | 0 | for (; i < 10; i++) |
1381 | 0 | adfParam[i] = 0.0; |
1382 | |
|
1383 | | #if LIBGEOTIFF_VERSION <= 1730 |
1384 | | // libgeotiff <= 1.7.3 is unfortunately inconsistent. When it synthetizes the |
1385 | | // projection parameters from the EPSG ProjectedCRS code, it returns |
1386 | | // them normalized in degrees. But when it gets them from |
1387 | | // ProjCoordTransGeoKey and other Proj....GeoKey's it return them in |
1388 | | // a raw way, that is in the units of GeogAngularUnitSizeGeoKey |
1389 | | // The below oSRS.SetXXXXX() methods assume the angular projection |
1390 | | // parameters to be in degrees, so convert them to degrees in that later case. |
1391 | | // From GDAL 3.0 to 3.9.0, we didn't do that conversion... |
1392 | | // And all versions of GDAL <= 3.9.0 when writing those geokeys, wrote |
1393 | | // them as degrees, hence this GTIFF_READ_ANGULAR_PARAMS_IN_DEGREE |
1394 | | // config option that can be set to YES to avoid that conversion and |
1395 | | // assume that the angular parameters have been written as degree. |
1396 | | if (GDALGTIFKeyGetSHORT(hGTIF, ProjCoordTransGeoKey, &tmp, 0, 1) && |
1397 | | !CPLTestBool(CPLGetConfigOption( |
1398 | | "GTIFF_READ_ANGULAR_PARAMS_IN_DEGREE", "NO"))) |
1399 | | { |
1400 | | adfParam[0] *= psDefn->UOMAngleInDegrees; |
1401 | | adfParam[1] *= psDefn->UOMAngleInDegrees; |
1402 | | adfParam[2] *= psDefn->UOMAngleInDegrees; |
1403 | | adfParam[3] *= psDefn->UOMAngleInDegrees; |
1404 | | } |
1405 | | #else |
1406 | | // If GTIFF_READ_ANGULAR_PARAMS_IN_DEGREE=YES (non-nominal case), undo |
1407 | | // the conversion to degrees, that has been done by libgeotiff > 1.7.3 |
1408 | 0 | if (GDALGTIFKeyGetSHORT(hGTIF, ProjCoordTransGeoKey, &tmp, 0, 1) && |
1409 | 0 | psDefn->UOMAngleInDegrees != 0 && psDefn->UOMAngleInDegrees != 1 && |
1410 | 0 | CPLTestBool(CPLGetConfigOption( |
1411 | 0 | "GTIFF_READ_ANGULAR_PARAMS_IN_DEGREE", "NO"))) |
1412 | 0 | { |
1413 | 0 | adfParam[0] /= psDefn->UOMAngleInDegrees; |
1414 | 0 | adfParam[1] /= psDefn->UOMAngleInDegrees; |
1415 | 0 | adfParam[2] /= psDefn->UOMAngleInDegrees; |
1416 | 0 | adfParam[3] /= psDefn->UOMAngleInDegrees; |
1417 | 0 | } |
1418 | 0 | #endif |
1419 | | |
1420 | | /* -------------------------------------------------------------------- |
1421 | | */ |
1422 | | /* Translation the fundamental projection. */ |
1423 | | /* -------------------------------------------------------------------- |
1424 | | */ |
1425 | 0 | switch (psDefn->CTProjection) |
1426 | 0 | { |
1427 | 0 | case CT_TransverseMercator: |
1428 | 0 | oSRS.SetTM(adfParam[0], adfParam[1], adfParam[4], adfParam[5], |
1429 | 0 | adfParam[6]); |
1430 | 0 | break; |
1431 | | |
1432 | 0 | case CT_TransvMercator_SouthOriented: |
1433 | 0 | oSRS.SetTMSO(adfParam[0], adfParam[1], adfParam[4], adfParam[5], |
1434 | 0 | adfParam[6]); |
1435 | 0 | break; |
1436 | | |
1437 | 0 | case CT_Mercator: |
1438 | | // If a lat_ts was specified use 2SP, otherwise use 1SP. |
1439 | 0 | if (psDefn->ProjParmId[2] == ProjStdParallel1GeoKey) |
1440 | 0 | { |
1441 | 0 | if (psDefn->ProjParmId[4] == ProjScaleAtNatOriginGeoKey) |
1442 | 0 | CPLError(CE_Warning, CPLE_AppDefined, |
1443 | 0 | "Mercator projection should not define " |
1444 | 0 | "both StdParallel1 and ScaleAtNatOrigin. " |
1445 | 0 | "Using StdParallel1 and ignoring " |
1446 | 0 | "ScaleAtNatOrigin."); |
1447 | 0 | oSRS.SetMercator2SP(adfParam[2], adfParam[0], adfParam[1], |
1448 | 0 | adfParam[5], adfParam[6]); |
1449 | 0 | } |
1450 | 0 | else |
1451 | 0 | oSRS.SetMercator(adfParam[0], adfParam[1], adfParam[4], |
1452 | 0 | adfParam[5], adfParam[6]); |
1453 | | |
1454 | | // Override hack for google mercator. |
1455 | 0 | if (psDefn->Projection == 1024 || psDefn->Projection == 9841) |
1456 | 0 | { |
1457 | 0 | oSRS.SetExtension( |
1458 | 0 | "PROJCS", "PROJ4", |
1459 | 0 | "+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 " |
1460 | 0 | "+lon_0=0.0 " |
1461 | 0 | "+x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null " |
1462 | 0 | "+wktext +no_defs"); // TODO(schwehr): Why 2 spaces? |
1463 | 0 | } |
1464 | 0 | break; |
1465 | | |
1466 | 0 | case CT_ObliqueStereographic: |
1467 | 0 | oSRS.SetOS(adfParam[0], adfParam[1], adfParam[4], adfParam[5], |
1468 | 0 | adfParam[6]); |
1469 | 0 | break; |
1470 | | |
1471 | 0 | case CT_Stereographic: |
1472 | 0 | oSRS.SetStereographic(adfParam[0], adfParam[1], adfParam[4], |
1473 | 0 | adfParam[5], adfParam[6]); |
1474 | 0 | break; |
1475 | | |
1476 | 0 | case CT_ObliqueMercator: // Hotine. |
1477 | 0 | oSRS.SetHOM(adfParam[0], adfParam[1], adfParam[2], adfParam[3], |
1478 | 0 | adfParam[4], adfParam[5], adfParam[6]); |
1479 | 0 | break; |
1480 | | |
1481 | 0 | case CT_HotineObliqueMercatorAzimuthCenter: |
1482 | 0 | oSRS.SetHOMAC(adfParam[0], adfParam[1], adfParam[2], |
1483 | 0 | adfParam[3], adfParam[4], adfParam[5], |
1484 | 0 | adfParam[6]); |
1485 | 0 | break; |
1486 | | |
1487 | 0 | case CT_ObliqueMercator_Laborde: |
1488 | 0 | oSRS.SetLOM(adfParam[0], adfParam[1], adfParam[2], adfParam[4], |
1489 | 0 | adfParam[5], adfParam[6]); |
1490 | 0 | break; |
1491 | | |
1492 | 0 | case CT_EquidistantConic: |
1493 | 0 | oSRS.SetEC(adfParam[0], adfParam[1], adfParam[2], adfParam[3], |
1494 | 0 | adfParam[5], adfParam[6]); |
1495 | 0 | break; |
1496 | | |
1497 | 0 | case CT_CassiniSoldner: |
1498 | 0 | oSRS.SetCS(adfParam[0], adfParam[1], adfParam[5], adfParam[6]); |
1499 | 0 | break; |
1500 | | |
1501 | 0 | case CT_Polyconic: |
1502 | 0 | oSRS.SetPolyconic(adfParam[0], adfParam[1], adfParam[5], |
1503 | 0 | adfParam[6]); |
1504 | 0 | break; |
1505 | | |
1506 | 0 | case CT_AzimuthalEquidistant: |
1507 | 0 | oSRS.SetAE(adfParam[0], adfParam[1], adfParam[5], adfParam[6]); |
1508 | 0 | break; |
1509 | | |
1510 | 0 | case CT_MillerCylindrical: |
1511 | 0 | oSRS.SetMC(adfParam[0], adfParam[1], adfParam[5], adfParam[6]); |
1512 | 0 | break; |
1513 | | |
1514 | 0 | case CT_Equirectangular: |
1515 | 0 | oSRS.SetEquirectangular2(adfParam[0], adfParam[1], adfParam[2], |
1516 | 0 | adfParam[5], adfParam[6]); |
1517 | 0 | break; |
1518 | | |
1519 | 0 | case CT_Gnomonic: |
1520 | 0 | oSRS.SetGnomonic(adfParam[0], adfParam[1], adfParam[5], |
1521 | 0 | adfParam[6]); |
1522 | 0 | break; |
1523 | | |
1524 | 0 | case CT_LambertAzimEqualArea: |
1525 | 0 | oSRS.SetLAEA(adfParam[0], adfParam[1], adfParam[5], |
1526 | 0 | adfParam[6]); |
1527 | 0 | break; |
1528 | | |
1529 | 0 | case CT_Orthographic: |
1530 | 0 | oSRS.SetOrthographic(adfParam[0], adfParam[1], adfParam[5], |
1531 | 0 | adfParam[6]); |
1532 | 0 | break; |
1533 | | |
1534 | 0 | case CT_Robinson: |
1535 | 0 | oSRS.SetRobinson(adfParam[1], adfParam[5], adfParam[6]); |
1536 | 0 | break; |
1537 | | |
1538 | 0 | case CT_Sinusoidal: |
1539 | 0 | oSRS.SetSinusoidal(adfParam[1], adfParam[5], adfParam[6]); |
1540 | 0 | break; |
1541 | | |
1542 | 0 | case CT_VanDerGrinten: |
1543 | 0 | oSRS.SetVDG(adfParam[1], adfParam[5], adfParam[6]); |
1544 | 0 | break; |
1545 | | |
1546 | 0 | case CT_PolarStereographic: |
1547 | 0 | oSRS.SetPS(adfParam[0], adfParam[1], adfParam[4], adfParam[5], |
1548 | 0 | adfParam[6]); |
1549 | 0 | break; |
1550 | | |
1551 | 0 | case CT_LambertConfConic_2SP: |
1552 | 0 | oSRS.SetLCC(adfParam[2], adfParam[3], adfParam[0], adfParam[1], |
1553 | 0 | adfParam[5], adfParam[6]); |
1554 | 0 | break; |
1555 | | |
1556 | 0 | case CT_LambertConfConic_1SP: |
1557 | 0 | oSRS.SetLCC1SP(adfParam[0], adfParam[1], adfParam[4], |
1558 | 0 | adfParam[5], adfParam[6]); |
1559 | 0 | break; |
1560 | | |
1561 | 0 | case CT_AlbersEqualArea: |
1562 | 0 | oSRS.SetACEA(adfParam[0], adfParam[1], adfParam[2], adfParam[3], |
1563 | 0 | adfParam[5], adfParam[6]); |
1564 | 0 | break; |
1565 | | |
1566 | 0 | case CT_NewZealandMapGrid: |
1567 | 0 | oSRS.SetNZMG(adfParam[0], adfParam[1], adfParam[5], |
1568 | 0 | adfParam[6]); |
1569 | 0 | break; |
1570 | | |
1571 | 0 | case CT_CylindricalEqualArea: |
1572 | 0 | oSRS.SetCEA(adfParam[0], adfParam[1], adfParam[5], adfParam[6]); |
1573 | 0 | break; |
1574 | 0 | default: |
1575 | 0 | if (oSRS.IsProjected()) |
1576 | 0 | { |
1577 | 0 | const char *pszName = oSRS.GetName(); |
1578 | 0 | std::string osName(pszName ? pszName : "unnamed"); |
1579 | 0 | oSRS.Clear(); |
1580 | 0 | oSRS.SetLocalCS(osName.c_str()); |
1581 | 0 | } |
1582 | 0 | break; |
1583 | 0 | } |
1584 | 0 | } |
1585 | | |
1586 | 0 | if (psDefn->Model == ModelTypeProjected && psDefn->PCS != KvUserDefined && |
1587 | 0 | !bGotFromEPSG) |
1588 | 0 | { |
1589 | 0 | OGRSpatialReference oSRSTest(oSRS); |
1590 | 0 | OGRSpatialReference oSRSTmp; |
1591 | |
|
1592 | 0 | const bool bPCSCodeValid = |
1593 | 0 | oSRSTmp.importFromEPSG(psDefn->PCS) == OGRERR_NONE; |
1594 | 0 | oSRSTmp.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER); |
1595 | | |
1596 | | // Force axis to avoid issues with SRS with northing, easting order |
1597 | 0 | oSRSTest.SetAxes(nullptr, "X", OAO_East, "Y", OAO_North); |
1598 | 0 | oSRSTmp.SetAxes(nullptr, "X", OAO_East, "Y", OAO_North); |
1599 | |
|
1600 | 0 | const std::string osGTiffSRSSource = |
1601 | 0 | CPLGetConfigOption("GTIFF_SRS_SOURCE", ""); |
1602 | 0 | const char *const apszOptions[] = { |
1603 | 0 | "IGNORE_DATA_AXIS_TO_SRS_AXIS_MAPPING=YES", nullptr}; |
1604 | 0 | if (!bHasWarnedInconsistentGeogCRSEPSG && |
1605 | 0 | !oSRSTmp.IsSame(&oSRS, apszOptions) && |
1606 | 0 | EQUAL(osGTiffSRSSource.c_str(), "")) |
1607 | 0 | { |
1608 | | // See https://github.com/OSGeo/gdal/issues/5399 |
1609 | | // where a file has inconsistent GeogSemiMinorAxisGeoKey / |
1610 | | // GeogInvFlatteningGeoKey values, which cause its datum to be |
1611 | | // considered as non-equivalent to the EPSG one. |
1612 | 0 | CPLError( |
1613 | 0 | CE_Warning, CPLE_AppDefined, |
1614 | 0 | "The definition of projected CRS EPSG:%d got from GeoTIFF keys " |
1615 | 0 | "is not the same as the one from the EPSG registry, " |
1616 | 0 | "which may cause issues during reprojection operations. " |
1617 | 0 | "Set GTIFF_SRS_SOURCE configuration option to EPSG to " |
1618 | 0 | "use official parameters (overriding the ones from GeoTIFF " |
1619 | 0 | "keys), " |
1620 | 0 | "or to GEOKEYS to use custom values from GeoTIFF keys " |
1621 | 0 | "and drop the EPSG code.", |
1622 | 0 | psDefn->PCS); |
1623 | 0 | } |
1624 | 0 | if (EQUAL(osGTiffSRSSource.c_str(), "EPSG")) |
1625 | 0 | { |
1626 | 0 | oSRS = std::move(oSRSTmp); |
1627 | 0 | } |
1628 | 0 | else if (bPCSCodeValid && EQUAL(osGTiffSRSSource.c_str(), "")) |
1629 | 0 | { |
1630 | 0 | oSRS.SetAuthority(nullptr, "EPSG", psDefn->PCS); |
1631 | 0 | } |
1632 | 0 | } |
1633 | |
|
1634 | 0 | if (oSRS.IsProjected() && oSRS.GetAxesCount() == 2) |
1635 | 0 | { |
1636 | 0 | const char *pszProjCRSName = oSRS.GetAttrValue("PROJCS"); |
1637 | 0 | if (pszProjCRSName) |
1638 | 0 | { |
1639 | | // Hack to be able to read properly what we have written for |
1640 | | // ESRI:102113 (ESRI ancient WebMercator). |
1641 | 0 | if (EQUAL(pszProjCRSName, "WGS_1984_Web_Mercator")) |
1642 | 0 | oSRS.SetFromUserInput("ESRI:102113"); |
1643 | | // And for EPSG:900913 |
1644 | 0 | else if (EQUAL(pszProjCRSName, "Google Maps Global Mercator")) |
1645 | 0 | oSRS.importFromEPSG(900913); |
1646 | 0 | else if (strchr(pszProjCRSName, '_')) |
1647 | 0 | { |
1648 | | // Morph from potential ESRI name |
1649 | 0 | char *pszOfficialName = GTIFGetEPSGOfficialName( |
1650 | 0 | hGTIF, PJ_TYPE_PROJECTED_CRS, pszProjCRSName); |
1651 | 0 | if (pszOfficialName) |
1652 | 0 | { |
1653 | 0 | oSRS.SetProjCS(pszOfficialName); |
1654 | 0 | CPLFree(pszOfficialName); |
1655 | 0 | } |
1656 | 0 | } |
1657 | 0 | } |
1658 | 0 | } |
1659 | | |
1660 | | /* ==================================================================== */ |
1661 | | /* Handle vertical coordinate system information if we have it. */ |
1662 | | /* ==================================================================== */ |
1663 | 0 | bool bNeedManualVertCS = false; |
1664 | 0 | char citation[2048] = {'\0'}; |
1665 | | |
1666 | | // See https://github.com/OSGeo/gdal/pull/4197 |
1667 | 0 | if (verticalCSType > KvUserDefined || verticalDatum > KvUserDefined || |
1668 | 0 | verticalUnits > KvUserDefined) |
1669 | 0 | { |
1670 | 0 | CPLError(CE_Warning, CPLE_AppDefined, |
1671 | 0 | "At least one of VerticalCSTypeGeoKey, VerticalDatumGeoKey or " |
1672 | 0 | "VerticalUnitsGeoKey has a value in the private user range. " |
1673 | 0 | "Ignoring vertical information."); |
1674 | 0 | verticalCSType = 0; |
1675 | 0 | verticalDatum = 0; |
1676 | 0 | verticalUnits = 0; |
1677 | 0 | } |
1678 | |
|
1679 | 0 | if ((verticalCSType != 0 || verticalDatum != 0 || verticalUnits != 0) && |
1680 | 0 | (oSRS.IsGeographic() || oSRS.IsProjected() || oSRS.IsLocal())) |
1681 | 0 | { |
1682 | 0 | std::string osVertCRSName; |
1683 | 0 | if (GDALGTIFKeyGetASCII(hGTIF, VerticalCitationGeoKey, citation, |
1684 | 0 | sizeof(citation))) |
1685 | 0 | { |
1686 | 0 | if (STARTS_WITH_CI(citation, "VCS Name = ")) |
1687 | 0 | { |
1688 | 0 | memmove(citation, citation + strlen("VCS Name = "), |
1689 | 0 | strlen(citation + strlen("VCS Name = ")) + 1); |
1690 | 0 | char *pszPipeChar = strchr(citation, '|'); |
1691 | 0 | if (pszPipeChar) |
1692 | 0 | *pszPipeChar = '\0'; |
1693 | 0 | osVertCRSName = citation; |
1694 | 0 | } |
1695 | 0 | } |
1696 | |
|
1697 | 0 | OGRSpatialReference oVertSRS; |
1698 | 0 | bool bCanBuildCompoundCRS = oSRS.GetRoot() != nullptr; |
1699 | 0 | if (verticalCSType != KvUserDefined && verticalCSType > 0) |
1700 | 0 | { |
1701 | 0 | if (!(oVertSRS.importFromEPSG(verticalCSType) == OGRERR_NONE && |
1702 | 0 | oVertSRS.IsVertical())) |
1703 | 0 | { |
1704 | 0 | bCanBuildCompoundCRS = false; |
1705 | 0 | } |
1706 | 0 | else |
1707 | 0 | { |
1708 | 0 | osVertCRSName = oVertSRS.GetName(); |
1709 | 0 | } |
1710 | 0 | } |
1711 | 0 | if (osVertCRSName.empty()) |
1712 | 0 | osVertCRSName = "unknown"; |
1713 | |
|
1714 | 0 | if (bCanBuildCompoundCRS) |
1715 | 0 | { |
1716 | 0 | const bool bHorizontalHasCode = |
1717 | 0 | oSRS.GetAuthorityCode(nullptr) != nullptr; |
1718 | 0 | const char *pszHorizontalName = oSRS.GetName(); |
1719 | 0 | const std::string osHorizontalName( |
1720 | 0 | pszHorizontalName ? pszHorizontalName : "unnamed"); |
1721 | | /* -------------------------------------------------------------------- |
1722 | | */ |
1723 | | /* Promote to being a compound coordinate system. */ |
1724 | | /* -------------------------------------------------------------------- |
1725 | | */ |
1726 | 0 | OGR_SRSNode *poOldRoot = oSRS.GetRoot()->Clone(); |
1727 | |
|
1728 | 0 | oSRS.Clear(); |
1729 | | |
1730 | | /* -------------------------------------------------------------------- |
1731 | | */ |
1732 | | /* Set COMPD_CS name. */ |
1733 | | /* -------------------------------------------------------------------- |
1734 | | */ |
1735 | 0 | char szCTString[512]; |
1736 | 0 | szCTString[0] = '\0'; |
1737 | 0 | if (GDALGTIFKeyGetASCII(hGTIF, GTCitationGeoKey, szCTString, |
1738 | 0 | sizeof(szCTString)) && |
1739 | 0 | strstr(szCTString, " = ") == nullptr) |
1740 | 0 | { |
1741 | 0 | oSRS.SetNode("COMPD_CS", szCTString); |
1742 | 0 | } |
1743 | 0 | else |
1744 | 0 | { |
1745 | 0 | oSRS.SetNode( |
1746 | 0 | "COMPD_CS", |
1747 | 0 | (osHorizontalName + " + " + osVertCRSName).c_str()); |
1748 | 0 | } |
1749 | |
|
1750 | 0 | oSRS.GetRoot()->AddChild(poOldRoot); |
1751 | | |
1752 | | /* -------------------------------------------------------------------- |
1753 | | */ |
1754 | | /* If we have the vertical cs, try to look it up, and use the |
1755 | | */ |
1756 | | /* definition provided by that. */ |
1757 | | /* -------------------------------------------------------------------- |
1758 | | */ |
1759 | 0 | bNeedManualVertCS = true; |
1760 | |
|
1761 | 0 | if (!oVertSRS.IsEmpty()) |
1762 | 0 | { |
1763 | 0 | oSRS.GetRoot()->AddChild(oVertSRS.GetRoot()->Clone()); |
1764 | 0 | bNeedManualVertCS = false; |
1765 | | |
1766 | | // GeoTIFF doesn't store EPSG code of CompoundCRS, so |
1767 | | // if we have an EPSG code for both the horizontal and vertical |
1768 | | // parts, check if there's a known CompoundCRS associating |
1769 | | // both |
1770 | 0 | if (bHorizontalHasCode && verticalCSType != KvUserDefined && |
1771 | 0 | verticalCSType > 0) |
1772 | 0 | { |
1773 | 0 | const auto *poSRSMatch = oSRS.FindBestMatch(100); |
1774 | 0 | if (poSRSMatch) |
1775 | 0 | oSRS = *poSRSMatch; |
1776 | 0 | delete poSRSMatch; |
1777 | 0 | } |
1778 | 0 | } |
1779 | 0 | } |
1780 | 0 | } |
1781 | | |
1782 | | /* -------------------------------------------------------------------- */ |
1783 | | /* Collect some information from the VerticalCS if not provided */ |
1784 | | /* via geokeys. */ |
1785 | | /* -------------------------------------------------------------------- */ |
1786 | 0 | if (bNeedManualVertCS) |
1787 | 0 | { |
1788 | 0 | FillCompoundCRSWithManualVertCS(hGTIF, oSRS, citation, verticalDatum, |
1789 | 0 | verticalUnits); |
1790 | 0 | } |
1791 | | |
1792 | | // Hack for tiff_read.py:test_tiff_grads so as to normalize angular |
1793 | | // parameters to grad |
1794 | 0 | if (psDefn->UOMAngleInDegrees != 1.0) |
1795 | 0 | { |
1796 | 0 | char *pszWKT = nullptr; |
1797 | 0 | const char *const apszOptions[] = { |
1798 | 0 | "FORMAT=WKT1", "ADD_TOWGS84_ON_EXPORT_TO_WKT1=NO", nullptr}; |
1799 | 0 | if (oSRS.exportToWkt(&pszWKT, apszOptions) == OGRERR_NONE) |
1800 | 0 | { |
1801 | 0 | const char *const apszOptionsImport[] = { |
1802 | 0 | #if PROJ_AT_LEAST_VERSION(9, 1, 0) |
1803 | 0 | "UNSET_IDENTIFIERS_IF_INCOMPATIBLE_DEF=NO", |
1804 | 0 | #endif |
1805 | 0 | nullptr |
1806 | 0 | }; |
1807 | 0 | oSRS.importFromWkt(pszWKT, apszOptionsImport); |
1808 | 0 | } |
1809 | 0 | CPLFree(pszWKT); |
1810 | 0 | } |
1811 | |
|
1812 | 0 | oSRS.StripTOWGS84IfKnownDatumAndAllowed(); |
1813 | |
|
1814 | 0 | double dfCoordinateEpoch = 0.0; |
1815 | 0 | if (GDALGTIFKeyGetDOUBLE(hGTIF, CoordinateEpochGeoKey, &dfCoordinateEpoch, |
1816 | 0 | 0, 1)) |
1817 | 0 | { |
1818 | 0 | oSRS.SetCoordinateEpoch(dfCoordinateEpoch); |
1819 | 0 | } |
1820 | |
|
1821 | 0 | return OGRSpatialReference::ToHandle(oSRS.Clone()); |
1822 | 0 | } |
1823 | | |
1824 | | /************************************************************************/ |
1825 | | /* GTIFGetOGISDefn() */ |
1826 | | /************************************************************************/ |
1827 | | |
1828 | | char *GTIFGetOGISDefn(GTIF *hGTIF, GTIFDefn *psDefn) |
1829 | 0 | { |
1830 | 0 | OGRSpatialReferenceH hSRS = GTIFGetOGISDefnAsOSR(hGTIF, psDefn); |
1831 | |
|
1832 | 0 | char *pszWKT = nullptr; |
1833 | 0 | if (hSRS && OGRSpatialReference::FromHandle(hSRS)->exportToWkt(&pszWKT) == |
1834 | 0 | OGRERR_NONE) |
1835 | 0 | { |
1836 | 0 | OSRDestroySpatialReference(hSRS); |
1837 | 0 | return pszWKT; |
1838 | 0 | } |
1839 | 0 | CPLFree(pszWKT); |
1840 | 0 | OSRDestroySpatialReference(hSRS); |
1841 | |
|
1842 | 0 | return nullptr; |
1843 | 0 | } |
1844 | | |
1845 | | /************************************************************************/ |
1846 | | /* OGCDatumName2EPSGDatumCode() */ |
1847 | | /************************************************************************/ |
1848 | | |
1849 | | static int OGCDatumName2EPSGDatumCode(GTIF *psGTIF, const char *pszOGCName) |
1850 | | |
1851 | 0 | { |
1852 | 0 | int nReturn = KvUserDefined; |
1853 | | |
1854 | | /* -------------------------------------------------------------------- */ |
1855 | | /* Do we know it as a built in? */ |
1856 | | /* -------------------------------------------------------------------- */ |
1857 | 0 | if (EQUAL(pszOGCName, "NAD27") || |
1858 | 0 | EQUAL(pszOGCName, "North_American_Datum_1927")) |
1859 | 0 | return Datum_North_American_Datum_1927; |
1860 | 0 | else if (EQUAL(pszOGCName, "NAD83") || |
1861 | 0 | EQUAL(pszOGCName, "North_American_Datum_1983")) |
1862 | 0 | return Datum_North_American_Datum_1983; |
1863 | 0 | else if (EQUAL(pszOGCName, "WGS84") || EQUAL(pszOGCName, "WGS_1984") || |
1864 | 0 | EQUAL(pszOGCName, "WGS 84")) |
1865 | 0 | return Datum_WGS84; |
1866 | 0 | else if (EQUAL(pszOGCName, "WGS72") || EQUAL(pszOGCName, "WGS_1972")) |
1867 | 0 | return Datum_WGS72; |
1868 | | |
1869 | | /* Search in database */ |
1870 | 0 | auto ctx = |
1871 | 0 | static_cast<PJ_CONTEXT *>(GTIFGetPROJContext(psGTIF, true, nullptr)); |
1872 | 0 | const PJ_TYPE searchType = PJ_TYPE_GEODETIC_REFERENCE_FRAME; |
1873 | 0 | auto list = proj_create_from_name(ctx, "EPSG", pszOGCName, &searchType, 1, |
1874 | 0 | true, /* approximate match */ |
1875 | 0 | 10, nullptr); |
1876 | 0 | if (list) |
1877 | 0 | { |
1878 | 0 | const auto listSize = proj_list_get_count(list); |
1879 | 0 | for (int i = 0; nReturn == KvUserDefined && i < listSize; i++) |
1880 | 0 | { |
1881 | 0 | auto datum = proj_list_get(ctx, list, i); |
1882 | 0 | if (datum) |
1883 | 0 | { |
1884 | 0 | const char *pszDatumName = proj_get_name(datum); |
1885 | 0 | if (pszDatumName) |
1886 | 0 | { |
1887 | 0 | char *pszTmp = CPLStrdup(pszDatumName); |
1888 | 0 | WKTMassageDatum(&pszTmp); |
1889 | 0 | if (EQUAL(pszTmp, pszOGCName)) |
1890 | 0 | { |
1891 | 0 | const char *pszCode = proj_get_id_code(datum, 0); |
1892 | 0 | if (pszCode) |
1893 | 0 | { |
1894 | 0 | nReturn = atoi(pszCode); |
1895 | 0 | } |
1896 | 0 | } |
1897 | 0 | CPLFree(pszTmp); |
1898 | 0 | } |
1899 | 0 | } |
1900 | 0 | proj_destroy(datum); |
1901 | 0 | } |
1902 | 0 | proj_list_destroy(list); |
1903 | 0 | } |
1904 | |
|
1905 | 0 | return nReturn; |
1906 | 0 | } |
1907 | | |
1908 | | /************************************************************************/ |
1909 | | /* GTIFSetFromOGISDefn() */ |
1910 | | /* */ |
1911 | | /* Write GeoTIFF projection tags from an OGC WKT definition. */ |
1912 | | /************************************************************************/ |
1913 | | |
1914 | | int GTIFSetFromOGISDefn(GTIF *psGTIF, const char *pszOGCWKT) |
1915 | | |
1916 | 0 | { |
1917 | | /* -------------------------------------------------------------------- */ |
1918 | | /* Create an OGRSpatialReference object corresponding to the */ |
1919 | | /* string. */ |
1920 | | /* -------------------------------------------------------------------- */ |
1921 | |
|
1922 | 0 | OGRSpatialReference oSRS; |
1923 | 0 | if (oSRS.importFromWkt(pszOGCWKT) != OGRERR_NONE) |
1924 | 0 | { |
1925 | 0 | return FALSE; |
1926 | 0 | } |
1927 | 0 | return GTIFSetFromOGISDefnEx(psGTIF, OGRSpatialReference::ToHandle(&oSRS), |
1928 | 0 | GEOTIFF_KEYS_STANDARD, GEOTIFF_VERSION_1_0); |
1929 | 0 | } |
1930 | | |
1931 | | int GTIFSetFromOGISDefnEx(GTIF *psGTIF, OGRSpatialReferenceH hSRS, |
1932 | | GTIFFKeysFlavorEnum eFlavor, |
1933 | | GeoTIFFVersionEnum eVersion) |
1934 | 0 | { |
1935 | 0 | std::map<geokey_t, std::string> oMapAsciiKeys; |
1936 | |
|
1937 | 0 | GTIFKeySet(psGTIF, GTRasterTypeGeoKey, TYPE_SHORT, 1, RasterPixelIsArea); |
1938 | |
|
1939 | 0 | const OGRSpatialReference *poSRS = OGRSpatialReference::FromHandle(hSRS); |
1940 | | |
1941 | | /* -------------------------------------------------------------------- */ |
1942 | | /* Set version number. */ |
1943 | | /* -------------------------------------------------------------------- */ |
1944 | 0 | if (eVersion == GEOTIFF_VERSION_AUTO) |
1945 | 0 | { |
1946 | 0 | if (poSRS->IsCompound() || |
1947 | 0 | (poSRS->IsGeographic() && poSRS->GetAxesCount() == 3)) |
1948 | 0 | { |
1949 | 0 | eVersion = GEOTIFF_VERSION_1_1; |
1950 | 0 | } |
1951 | 0 | else |
1952 | 0 | { |
1953 | 0 | eVersion = GEOTIFF_VERSION_1_0; |
1954 | 0 | } |
1955 | 0 | } |
1956 | 0 | CPLAssert(eVersion == GEOTIFF_VERSION_1_0 || |
1957 | 0 | eVersion == GEOTIFF_VERSION_1_1); |
1958 | 0 | if (eVersion >= GEOTIFF_VERSION_1_1) |
1959 | 0 | { |
1960 | 0 | #if LIBGEOTIFF_VERSION >= 1600 |
1961 | 0 | GTIFSetVersionNumbers(psGTIF, GEOTIFF_SPEC_1_1_VERSION, |
1962 | 0 | GEOTIFF_SPEC_1_1_KEY_REVISION, |
1963 | 0 | GEOTIFF_SPEC_1_1_MINOR_REVISION); |
1964 | | #else |
1965 | | CPLError(CE_Warning, CPLE_AppDefined, |
1966 | | "Setting GeoTIFF 1.1 requires libgeotiff >= 1.6. Key values " |
1967 | | "will be written as GeoTIFF 1.1, but the version number " |
1968 | | "will be seen as 1.0, which might confuse GeoTIFF readers"); |
1969 | | #endif |
1970 | 0 | } |
1971 | | |
1972 | | /* -------------------------------------------------------------------- */ |
1973 | | /* Get the ellipsoid definition. */ |
1974 | | /* -------------------------------------------------------------------- */ |
1975 | 0 | short nSpheroid = KvUserDefined; |
1976 | |
|
1977 | 0 | if (poSRS->GetAuthorityName("PROJCS|GEOGCS|DATUM|SPHEROID") != nullptr && |
1978 | 0 | EQUAL(poSRS->GetAuthorityName("PROJCS|GEOGCS|DATUM|SPHEROID"), "EPSG")) |
1979 | 0 | { |
1980 | 0 | nSpheroid = static_cast<short>( |
1981 | 0 | atoi(poSRS->GetAuthorityCode("PROJCS|GEOGCS|DATUM|SPHEROID"))); |
1982 | 0 | } |
1983 | 0 | else if (poSRS->GetAuthorityName("GEOGCS|DATUM|SPHEROID") != nullptr && |
1984 | 0 | EQUAL(poSRS->GetAuthorityName("GEOGCS|DATUM|SPHEROID"), "EPSG")) |
1985 | 0 | { |
1986 | 0 | nSpheroid = static_cast<short>( |
1987 | 0 | atoi(poSRS->GetAuthorityCode("GEOGCS|DATUM|SPHEROID"))); |
1988 | 0 | } |
1989 | |
|
1990 | 0 | OGRErr eErr = OGRERR_NONE; |
1991 | 0 | double dfSemiMajor = 0; |
1992 | 0 | double dfInvFlattening = 0; |
1993 | 0 | bool bHasEllipsoid = false; |
1994 | 0 | if (!poSRS->IsLocal()) |
1995 | 0 | { |
1996 | 0 | bHasEllipsoid = true; |
1997 | 0 | if (poSRS->IsCompound()) |
1998 | 0 | { |
1999 | 0 | OGRSpatialReference oSRSTmp(*poSRS); |
2000 | 0 | oSRSTmp.StripVertical(); |
2001 | 0 | bHasEllipsoid = CPL_TO_BOOL(!oSRSTmp.IsLocal()); |
2002 | 0 | } |
2003 | 0 | if (bHasEllipsoid) |
2004 | 0 | { |
2005 | 0 | dfSemiMajor = poSRS->GetSemiMajor(&eErr); |
2006 | 0 | dfInvFlattening = poSRS->GetInvFlattening(&eErr); |
2007 | 0 | if (eErr != OGRERR_NONE) |
2008 | 0 | { |
2009 | 0 | dfSemiMajor = 0.0; |
2010 | 0 | dfInvFlattening = 0.0; |
2011 | 0 | } |
2012 | 0 | } |
2013 | 0 | } |
2014 | | |
2015 | | /* -------------------------------------------------------------------- */ |
2016 | | /* Get the Datum so we can special case a few PCS codes. */ |
2017 | | /* -------------------------------------------------------------------- */ |
2018 | 0 | int nDatum = KvUserDefined; |
2019 | |
|
2020 | 0 | if (poSRS->GetAuthorityName("PROJCS|GEOGCS|DATUM") != nullptr && |
2021 | 0 | EQUAL(poSRS->GetAuthorityName("PROJCS|GEOGCS|DATUM"), "EPSG")) |
2022 | 0 | nDatum = atoi(poSRS->GetAuthorityCode("PROJCS|GEOGCS|DATUM")); |
2023 | 0 | else if (poSRS->GetAuthorityName("GEOGCS|DATUM") != nullptr && |
2024 | 0 | EQUAL(poSRS->GetAuthorityName("GEOGCS|DATUM"), "EPSG")) |
2025 | 0 | nDatum = atoi(poSRS->GetAuthorityCode("GEOGCS|DATUM")); |
2026 | 0 | else if (poSRS->GetAttrValue("DATUM") != nullptr) |
2027 | 0 | nDatum = |
2028 | 0 | OGCDatumName2EPSGDatumCode(psGTIF, poSRS->GetAttrValue("DATUM")); |
2029 | | |
2030 | | /* -------------------------------------------------------------------- */ |
2031 | | /* Get the GCS if possible. */ |
2032 | | /* -------------------------------------------------------------------- */ |
2033 | 0 | int nGCS = KvUserDefined; |
2034 | |
|
2035 | 0 | if (poSRS->GetAuthorityName("PROJCS|GEOGCS") != nullptr && |
2036 | 0 | EQUAL(poSRS->GetAuthorityName("PROJCS|GEOGCS"), "EPSG")) |
2037 | 0 | nGCS = atoi(poSRS->GetAuthorityCode("PROJCS|GEOGCS")); |
2038 | 0 | else if (poSRS->GetAuthorityName("GEOGCS") != nullptr && |
2039 | 0 | EQUAL(poSRS->GetAuthorityName("GEOGCS"), "EPSG")) |
2040 | 0 | nGCS = atoi(poSRS->GetAuthorityCode("GEOGCS")); |
2041 | |
|
2042 | 0 | int nVerticalCSKeyValue = 0; |
2043 | 0 | bool hasEllipsoidHeight = !poSRS->IsCompound() && poSRS->IsGeographic() && |
2044 | 0 | poSRS->GetAxesCount() == 3; |
2045 | 0 | if (nGCS == 4937 && eVersion >= GEOTIFF_VERSION_1_1) |
2046 | 0 | { |
2047 | | // Workaround a bug of PROJ 6.3.0 |
2048 | | // See https://github.com/OSGeo/PROJ/pull/1880 |
2049 | | // EPSG:4937 = ETRS89 3D |
2050 | 0 | hasEllipsoidHeight = true; |
2051 | 0 | nVerticalCSKeyValue = nGCS; |
2052 | 0 | nGCS = 4258; // ETRS89 2D |
2053 | 0 | } |
2054 | 0 | else if (nGCS != KvUserDefined) |
2055 | 0 | { |
2056 | 0 | OGRSpatialReference oGeogCRS; |
2057 | 0 | if (oGeogCRS.importFromEPSG(nGCS) == OGRERR_NONE && |
2058 | 0 | oGeogCRS.IsGeographic() && oGeogCRS.GetAxesCount() == 3) |
2059 | 0 | { |
2060 | 0 | hasEllipsoidHeight = true; |
2061 | 0 | if (eVersion >= GEOTIFF_VERSION_1_1) |
2062 | 0 | { |
2063 | 0 | const auto candidate_nVerticalCSKeyValue = nGCS; |
2064 | 0 | nGCS = KvUserDefined; |
2065 | | |
2066 | | // In case of a geographic 3D CRS, find the corresponding |
2067 | | // geographic 2D CRS |
2068 | 0 | auto ctx = static_cast<PJ_CONTEXT *>( |
2069 | 0 | GTIFGetPROJContext(psGTIF, true, nullptr)); |
2070 | 0 | const auto type = PJ_TYPE_GEOGRAPHIC_2D_CRS; |
2071 | 0 | auto list = proj_create_from_name(ctx, "EPSG", |
2072 | 0 | oGeogCRS.GetName(), &type, 1, |
2073 | 0 | false, // exact match |
2074 | 0 | 1, // result set limit size, |
2075 | 0 | nullptr); |
2076 | 0 | if (list && proj_list_get_count(list) == 1) |
2077 | 0 | { |
2078 | 0 | auto crs2D = proj_list_get(ctx, list, 0); |
2079 | 0 | if (crs2D) |
2080 | 0 | { |
2081 | 0 | const char *pszCode = proj_get_id_code(crs2D, 0); |
2082 | 0 | if (pszCode) |
2083 | 0 | { |
2084 | 0 | nVerticalCSKeyValue = candidate_nVerticalCSKeyValue; |
2085 | 0 | nGCS = atoi(pszCode); |
2086 | 0 | } |
2087 | 0 | proj_destroy(crs2D); |
2088 | 0 | } |
2089 | 0 | } |
2090 | 0 | proj_list_destroy(list); |
2091 | 0 | } |
2092 | 0 | } |
2093 | 0 | } |
2094 | | |
2095 | | // Deprecated way of encoding ellipsoidal height |
2096 | 0 | if (hasEllipsoidHeight && nVerticalCSKeyValue == 0) |
2097 | 0 | { |
2098 | 0 | if (nGCS == 4979 || nDatum == 6326 || nSpheroid == 7030) |
2099 | 0 | { |
2100 | 0 | nVerticalCSKeyValue = 5030; // WGS_84_ellipsoid |
2101 | 0 | if (nGCS == 4979 || nDatum == 6326) |
2102 | 0 | { |
2103 | 0 | nGCS = 4326; |
2104 | 0 | } |
2105 | 0 | } |
2106 | 0 | else if (nDatum >= 6001 && nDatum <= 6033) |
2107 | 0 | { |
2108 | 0 | nVerticalCSKeyValue = nDatum - 1000; |
2109 | 0 | } |
2110 | 0 | else if (nSpheroid >= 7001 && nSpheroid <= 7033) |
2111 | 0 | { |
2112 | 0 | nVerticalCSKeyValue = nSpheroid - 2000; |
2113 | 0 | } |
2114 | 0 | } |
2115 | |
|
2116 | 0 | if (nGCS > 32767) |
2117 | 0 | nGCS = KvUserDefined; |
2118 | | |
2119 | | /* -------------------------------------------------------------------- */ |
2120 | | /* Get the linear units. */ |
2121 | | /* -------------------------------------------------------------------- */ |
2122 | 0 | const char *pszLinearUOMNameTmp = nullptr; |
2123 | 0 | const double dfLinearUOM = poSRS->GetLinearUnits(&pszLinearUOMNameTmp); |
2124 | 0 | const std::string osLinearUOMName(pszLinearUOMNameTmp ? pszLinearUOMNameTmp |
2125 | 0 | : ""); |
2126 | 0 | int nUOMLengthCode = 9001; // Meters. |
2127 | |
|
2128 | 0 | if (poSRS->GetAuthorityName("PROJCS|UNIT") != nullptr && |
2129 | 0 | EQUAL(poSRS->GetAuthorityName("PROJCS|UNIT"), "EPSG") && |
2130 | 0 | poSRS->GetAttrNode("PROJCS|UNIT") != poSRS->GetAttrNode("GEOGCS|UNIT")) |
2131 | 0 | nUOMLengthCode = atoi(poSRS->GetAuthorityCode("PROJCS|UNIT")); |
2132 | 0 | else if (EQUAL(osLinearUOMName.c_str(), SRS_UL_FOOT) || |
2133 | 0 | fabs(dfLinearUOM - GTIFAtof(SRS_UL_FOOT_CONV)) < 0.0000001) |
2134 | 0 | nUOMLengthCode = 9002; // International foot. |
2135 | 0 | else if (EQUAL(osLinearUOMName.c_str(), SRS_UL_US_FOOT) || |
2136 | 0 | std::abs(dfLinearUOM - GTIFAtof(SRS_UL_US_FOOT_CONV)) < 0.0000001) |
2137 | 0 | nUOMLengthCode = 9003; // US survey foot. |
2138 | 0 | else if (fabs(dfLinearUOM - 1.0) > 0.00000001) |
2139 | 0 | nUOMLengthCode = KvUserDefined; |
2140 | | |
2141 | | /* -------------------------------------------------------------------- */ |
2142 | | /* Get some authority values. */ |
2143 | | /* -------------------------------------------------------------------- */ |
2144 | 0 | int nPCS = KvUserDefined; |
2145 | |
|
2146 | 0 | if (poSRS->GetAuthorityName("PROJCS") != nullptr && |
2147 | 0 | EQUAL(poSRS->GetAuthorityName("PROJCS"), "EPSG")) |
2148 | 0 | { |
2149 | 0 | nPCS = atoi(poSRS->GetAuthorityCode("PROJCS")); |
2150 | 0 | if (nPCS > 32767) |
2151 | 0 | nPCS = KvUserDefined; |
2152 | 0 | } |
2153 | | |
2154 | | /* -------------------------------------------------------------------- */ |
2155 | | /* Handle the projection transformation. */ |
2156 | | /* -------------------------------------------------------------------- */ |
2157 | 0 | bool bWritePEString = false; |
2158 | 0 | bool bUnknownProjection = false; |
2159 | |
|
2160 | 0 | const OGRSpatialReference *poSRSCompatibleOfWKT1 = poSRS; |
2161 | 0 | #if PROJ_AT_LEAST_VERSION(6, 3, 0) |
2162 | 0 | OGRSpatialReference oSRSCompatibleOfWKT1; |
2163 | 0 | if (poSRS->IsProjected() && !poSRS->IsCompound() && |
2164 | 0 | poSRS->GetAxesCount() == 3) |
2165 | 0 | { |
2166 | 0 | oSRSCompatibleOfWKT1 = *poSRS; |
2167 | 0 | oSRSCompatibleOfWKT1.DemoteTo2D(nullptr); |
2168 | 0 | poSRSCompatibleOfWKT1 = &oSRSCompatibleOfWKT1; |
2169 | 0 | bWritePEString = true; |
2170 | 0 | } |
2171 | 0 | #endif |
2172 | |
|
2173 | 0 | double dfAngUnitValue = 0; |
2174 | 0 | std::string osAngUnitName; |
2175 | 0 | if (bHasEllipsoid) |
2176 | 0 | { |
2177 | 0 | const char *angUnitNameTmp = ""; |
2178 | 0 | dfAngUnitValue = poSRS->GetAngularUnits(&angUnitNameTmp); |
2179 | 0 | osAngUnitName = angUnitNameTmp; |
2180 | 0 | } |
2181 | | |
2182 | | // Convert angular projection parameters from its normalized value in degree |
2183 | | // to the units of GeogAngularUnitsGeoKey. |
2184 | | // Note: for GDAL <= 3.9.0, we always have written them in degrees ! |
2185 | | // We can set GTIFF_WRITE_ANGULAR_PARAMS_IN_DEGREE=YES to get that |
2186 | | // non-conformant behavior... |
2187 | 0 | const auto ConvertAngularParam = [dfAngUnitValue](double dfValInDeg) |
2188 | 0 | { |
2189 | 0 | constexpr double DEG_TO_RAD = M_PI / 180.0; |
2190 | 0 | return dfAngUnitValue != 0 && |
2191 | 0 | std::fabs(dfAngUnitValue - DEG_TO_RAD) > 1e-10 && |
2192 | 0 | !CPLTestBool(CPLGetConfigOption( |
2193 | 0 | "GTIFF_WRITE_ANGULAR_PARAMS_IN_DEGREE", "NO")) |
2194 | 0 | ? dfValInDeg * DEG_TO_RAD / dfAngUnitValue |
2195 | 0 | : dfValInDeg; |
2196 | 0 | }; |
2197 | |
|
2198 | 0 | const char *pszProjection = |
2199 | 0 | poSRSCompatibleOfWKT1->GetAttrValue("PROJECTION"); |
2200 | 0 | if (nPCS != KvUserDefined) |
2201 | 0 | { |
2202 | | // If ESRI_PE flavor is explicitly required, then for EPSG:3857 |
2203 | | // we will have to write a completely non-standard definition |
2204 | | // that requires not setting GTModelTypeGeoKey to ProjectedCSTypeGeoKey. |
2205 | 0 | if (eFlavor == GEOTIFF_KEYS_ESRI_PE && nPCS == 3857) |
2206 | 0 | { |
2207 | 0 | bWritePEString = true; |
2208 | 0 | } |
2209 | 0 | else |
2210 | 0 | { |
2211 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2212 | 0 | ModelTypeProjected); |
2213 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, nPCS); |
2214 | 0 | } |
2215 | 0 | } |
2216 | 0 | else if (poSRSCompatibleOfWKT1->IsGeocentric()) |
2217 | 0 | { |
2218 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2219 | 0 | ModelTypeGeocentric); |
2220 | 0 | } |
2221 | 0 | else if (pszProjection == nullptr) |
2222 | 0 | { |
2223 | 0 | if (poSRSCompatibleOfWKT1->IsGeographic()) |
2224 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2225 | 0 | ModelTypeGeographic); |
2226 | | // Otherwise, presumably something like LOCAL_CS. |
2227 | 0 | } |
2228 | 0 | else if (EQUAL(pszProjection, SRS_PT_ALBERS_CONIC_EQUAL_AREA)) |
2229 | 0 | { |
2230 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2231 | 0 | ModelTypeProjected); |
2232 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2233 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2234 | |
|
2235 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2236 | 0 | CT_AlbersEqualArea); |
2237 | |
|
2238 | 0 | GTIFKeySet(psGTIF, ProjStdParallelGeoKey, TYPE_DOUBLE, 1, |
2239 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2240 | 0 | SRS_PP_STANDARD_PARALLEL_1, 0.0))); |
2241 | |
|
2242 | 0 | GTIFKeySet(psGTIF, ProjStdParallel2GeoKey, TYPE_DOUBLE, 1, |
2243 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2244 | 0 | SRS_PP_STANDARD_PARALLEL_2, 0.0))); |
2245 | |
|
2246 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1, |
2247 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2248 | 0 | SRS_PP_LATITUDE_OF_CENTER, 0.0))); |
2249 | |
|
2250 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1, |
2251 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2252 | 0 | SRS_PP_LONGITUDE_OF_CENTER, 0.0))); |
2253 | |
|
2254 | 0 | GTIFKeySet( |
2255 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2256 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2257 | |
|
2258 | 0 | GTIFKeySet( |
2259 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2260 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2261 | 0 | } |
2262 | | |
2263 | 0 | else if (poSRSCompatibleOfWKT1->GetUTMZone() != 0) |
2264 | 0 | { |
2265 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2266 | 0 | ModelTypeProjected); |
2267 | |
|
2268 | 0 | int bNorth = 0; |
2269 | 0 | const int nZone = poSRSCompatibleOfWKT1->GetUTMZone(&bNorth); |
2270 | |
|
2271 | 0 | if (nDatum == Datum_North_American_Datum_1983 && nZone >= 3 && |
2272 | 0 | nZone <= 22 && bNorth && nUOMLengthCode == 9001) |
2273 | 0 | { |
2274 | 0 | nPCS = 26900 + nZone; |
2275 | |
|
2276 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, nPCS); |
2277 | 0 | } |
2278 | 0 | else if (nDatum == Datum_North_American_Datum_1927 && nZone >= 3 && |
2279 | 0 | nZone <= 22 && bNorth && nUOMLengthCode == 9001) |
2280 | 0 | { |
2281 | 0 | nPCS = 26700 + nZone; |
2282 | |
|
2283 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, nPCS); |
2284 | 0 | } |
2285 | 0 | else if (nDatum == Datum_WGS84 && nUOMLengthCode == 9001) |
2286 | 0 | { |
2287 | 0 | if (bNorth) |
2288 | 0 | nPCS = 32600 + nZone; |
2289 | 0 | else |
2290 | 0 | nPCS = 32700 + nZone; |
2291 | |
|
2292 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, nPCS); |
2293 | 0 | } |
2294 | 0 | else |
2295 | 0 | { |
2296 | 0 | const int nProjection = nZone + (bNorth ? 16000 : 16100); |
2297 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, |
2298 | 0 | KvUserDefined); |
2299 | |
|
2300 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, nProjection); |
2301 | 0 | } |
2302 | 0 | } |
2303 | | |
2304 | 0 | else if (EQUAL(pszProjection, SRS_PT_TRANSVERSE_MERCATOR)) |
2305 | 0 | { |
2306 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2307 | 0 | ModelTypeProjected); |
2308 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2309 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2310 | |
|
2311 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2312 | 0 | CT_TransverseMercator); |
2313 | |
|
2314 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1, |
2315 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2316 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2317 | |
|
2318 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1, |
2319 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2320 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2321 | |
|
2322 | 0 | GTIFKeySet( |
2323 | 0 | psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1, |
2324 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_SCALE_FACTOR, 1.0)); |
2325 | |
|
2326 | 0 | GTIFKeySet( |
2327 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2328 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2329 | |
|
2330 | 0 | GTIFKeySet( |
2331 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2332 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2333 | 0 | } |
2334 | | |
2335 | 0 | else if (EQUAL(pszProjection, SRS_PT_TRANSVERSE_MERCATOR_SOUTH_ORIENTED)) |
2336 | 0 | { |
2337 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2338 | 0 | ModelTypeProjected); |
2339 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2340 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2341 | |
|
2342 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2343 | 0 | CT_TransvMercator_SouthOriented); |
2344 | |
|
2345 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1, |
2346 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2347 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2348 | |
|
2349 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1, |
2350 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2351 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2352 | |
|
2353 | 0 | GTIFKeySet( |
2354 | 0 | psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1, |
2355 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_SCALE_FACTOR, 1.0)); |
2356 | |
|
2357 | 0 | GTIFKeySet( |
2358 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2359 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2360 | |
|
2361 | 0 | GTIFKeySet( |
2362 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2363 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2364 | 0 | } |
2365 | | |
2366 | 0 | else if (EQUAL(pszProjection, SRS_PT_MERCATOR_2SP) || |
2367 | 0 | EQUAL(pszProjection, SRS_PT_MERCATOR_1SP)) |
2368 | | |
2369 | 0 | { |
2370 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2371 | 0 | ModelTypeProjected); |
2372 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2373 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2374 | |
|
2375 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, CT_Mercator); |
2376 | |
|
2377 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1, |
2378 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2379 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2380 | |
|
2381 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1, |
2382 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2383 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2384 | |
|
2385 | 0 | if (EQUAL(pszProjection, SRS_PT_MERCATOR_2SP)) |
2386 | 0 | GTIFKeySet( |
2387 | 0 | psGTIF, ProjStdParallel1GeoKey, TYPE_DOUBLE, 1, |
2388 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2389 | 0 | SRS_PP_STANDARD_PARALLEL_1, 0.0))); |
2390 | 0 | else |
2391 | 0 | GTIFKeySet(psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1, |
2392 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm( |
2393 | 0 | SRS_PP_SCALE_FACTOR, 1.0)); |
2394 | |
|
2395 | 0 | GTIFKeySet( |
2396 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2397 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2398 | |
|
2399 | 0 | GTIFKeySet( |
2400 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2401 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2402 | 0 | } |
2403 | | |
2404 | 0 | else if (EQUAL(pszProjection, SRS_PT_OBLIQUE_STEREOGRAPHIC)) |
2405 | 0 | { |
2406 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2407 | 0 | ModelTypeProjected); |
2408 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2409 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2410 | |
|
2411 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2412 | 0 | CT_ObliqueStereographic); |
2413 | |
|
2414 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1, |
2415 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2416 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2417 | |
|
2418 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1, |
2419 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2420 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2421 | |
|
2422 | 0 | GTIFKeySet( |
2423 | 0 | psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1, |
2424 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_SCALE_FACTOR, 1.0)); |
2425 | |
|
2426 | 0 | GTIFKeySet( |
2427 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2428 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2429 | |
|
2430 | 0 | GTIFKeySet( |
2431 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2432 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2433 | 0 | } |
2434 | | |
2435 | 0 | else if (EQUAL(pszProjection, SRS_PT_STEREOGRAPHIC)) |
2436 | 0 | { |
2437 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2438 | 0 | ModelTypeProjected); |
2439 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2440 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2441 | |
|
2442 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2443 | 0 | CT_Stereographic); |
2444 | |
|
2445 | 0 | GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1, |
2446 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2447 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2448 | |
|
2449 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2450 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2451 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2452 | |
|
2453 | 0 | GTIFKeySet( |
2454 | 0 | psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1, |
2455 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_SCALE_FACTOR, 1.0)); |
2456 | |
|
2457 | 0 | GTIFKeySet( |
2458 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2459 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2460 | |
|
2461 | 0 | GTIFKeySet( |
2462 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2463 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2464 | 0 | } |
2465 | | |
2466 | 0 | else if (EQUAL(pszProjection, SRS_PT_POLAR_STEREOGRAPHIC)) |
2467 | 0 | { |
2468 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2469 | 0 | ModelTypeProjected); |
2470 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2471 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2472 | |
|
2473 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2474 | 0 | CT_PolarStereographic); |
2475 | |
|
2476 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1, |
2477 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2478 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2479 | |
|
2480 | 0 | GTIFKeySet(psGTIF, ProjStraightVertPoleLongGeoKey, TYPE_DOUBLE, 1, |
2481 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2482 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2483 | |
|
2484 | 0 | GTIFKeySet( |
2485 | 0 | psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1, |
2486 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_SCALE_FACTOR, 1.0)); |
2487 | |
|
2488 | 0 | GTIFKeySet( |
2489 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2490 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2491 | |
|
2492 | 0 | GTIFKeySet( |
2493 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2494 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2495 | 0 | } |
2496 | | |
2497 | 0 | else if (EQUAL(pszProjection, SRS_PT_HOTINE_OBLIQUE_MERCATOR)) |
2498 | 0 | { |
2499 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2500 | 0 | ModelTypeProjected); |
2501 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2502 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2503 | |
|
2504 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2505 | 0 | CT_ObliqueMercator); |
2506 | |
|
2507 | 0 | GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1, |
2508 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2509 | 0 | SRS_PP_LATITUDE_OF_CENTER, 0.0))); |
2510 | |
|
2511 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2512 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2513 | 0 | SRS_PP_LONGITUDE_OF_CENTER, 0.0))); |
2514 | |
|
2515 | 0 | GTIFKeySet(psGTIF, ProjAzimuthAngleGeoKey, TYPE_DOUBLE, 1, |
2516 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_AZIMUTH, 0.0)); |
2517 | |
|
2518 | 0 | GTIFKeySet(psGTIF, ProjRectifiedGridAngleGeoKey, TYPE_DOUBLE, 1, |
2519 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2520 | 0 | SRS_PP_RECTIFIED_GRID_ANGLE, 0.0))); |
2521 | |
|
2522 | 0 | GTIFKeySet( |
2523 | 0 | psGTIF, ProjScaleAtCenterGeoKey, TYPE_DOUBLE, 1, |
2524 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_SCALE_FACTOR, 1.0)); |
2525 | |
|
2526 | 0 | GTIFKeySet( |
2527 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2528 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2529 | |
|
2530 | 0 | GTIFKeySet( |
2531 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2532 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2533 | 0 | } |
2534 | | |
2535 | 0 | else if (EQUAL(pszProjection, |
2536 | 0 | SRS_PT_HOTINE_OBLIQUE_MERCATOR_AZIMUTH_CENTER)) |
2537 | 0 | { |
2538 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2539 | 0 | ModelTypeProjected); |
2540 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2541 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2542 | |
|
2543 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2544 | 0 | CT_HotineObliqueMercatorAzimuthCenter); |
2545 | |
|
2546 | 0 | GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1, |
2547 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2548 | 0 | SRS_PP_LATITUDE_OF_CENTER, 0.0))); |
2549 | |
|
2550 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2551 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2552 | 0 | SRS_PP_LONGITUDE_OF_CENTER, 0.0))); |
2553 | |
|
2554 | 0 | GTIFKeySet(psGTIF, ProjAzimuthAngleGeoKey, TYPE_DOUBLE, 1, |
2555 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_AZIMUTH, 0.0)); |
2556 | |
|
2557 | 0 | GTIFKeySet(psGTIF, ProjRectifiedGridAngleGeoKey, TYPE_DOUBLE, 1, |
2558 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2559 | 0 | SRS_PP_RECTIFIED_GRID_ANGLE, 0.0))); |
2560 | |
|
2561 | 0 | GTIFKeySet( |
2562 | 0 | psGTIF, ProjScaleAtCenterGeoKey, TYPE_DOUBLE, 1, |
2563 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_SCALE_FACTOR, 1.0)); |
2564 | |
|
2565 | 0 | GTIFKeySet( |
2566 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2567 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2568 | |
|
2569 | 0 | GTIFKeySet( |
2570 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2571 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2572 | 0 | } |
2573 | | |
2574 | 0 | else if (EQUAL(pszProjection, "Laborde_Oblique_Mercator")) |
2575 | 0 | { |
2576 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2577 | 0 | ModelTypeProjected); |
2578 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2579 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2580 | |
|
2581 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2582 | 0 | CT_ObliqueMercator_Laborde); |
2583 | |
|
2584 | 0 | GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1, |
2585 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2586 | 0 | SRS_PP_LATITUDE_OF_CENTER, 0.0))); |
2587 | |
|
2588 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2589 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2590 | 0 | SRS_PP_LONGITUDE_OF_CENTER, 0.0))); |
2591 | |
|
2592 | 0 | GTIFKeySet(psGTIF, ProjAzimuthAngleGeoKey, TYPE_DOUBLE, 1, |
2593 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_AZIMUTH, 0.0)); |
2594 | |
|
2595 | 0 | GTIFKeySet( |
2596 | 0 | psGTIF, ProjScaleAtCenterGeoKey, TYPE_DOUBLE, 1, |
2597 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_SCALE_FACTOR, 1.0)); |
2598 | |
|
2599 | 0 | GTIFKeySet( |
2600 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2601 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2602 | |
|
2603 | 0 | GTIFKeySet( |
2604 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2605 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2606 | 0 | } |
2607 | | |
2608 | 0 | else if (EQUAL(pszProjection, SRS_PT_CASSINI_SOLDNER)) |
2609 | 0 | { |
2610 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2611 | 0 | ModelTypeProjected); |
2612 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2613 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2614 | |
|
2615 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2616 | 0 | CT_CassiniSoldner); |
2617 | |
|
2618 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1, |
2619 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2620 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2621 | |
|
2622 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1, |
2623 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2624 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2625 | |
|
2626 | 0 | GTIFKeySet( |
2627 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2628 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2629 | |
|
2630 | 0 | GTIFKeySet( |
2631 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2632 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2633 | 0 | } |
2634 | | |
2635 | 0 | else if (EQUAL(pszProjection, SRS_PT_EQUIDISTANT_CONIC)) |
2636 | 0 | { |
2637 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2638 | 0 | ModelTypeProjected); |
2639 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2640 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2641 | |
|
2642 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2643 | 0 | CT_EquidistantConic); |
2644 | |
|
2645 | 0 | GTIFKeySet(psGTIF, ProjStdParallel1GeoKey, TYPE_DOUBLE, 1, |
2646 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2647 | 0 | SRS_PP_STANDARD_PARALLEL_1, 0.0))); |
2648 | |
|
2649 | 0 | GTIFKeySet(psGTIF, ProjStdParallel2GeoKey, TYPE_DOUBLE, 1, |
2650 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2651 | 0 | SRS_PP_STANDARD_PARALLEL_2, 0.0))); |
2652 | |
|
2653 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1, |
2654 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2655 | 0 | SRS_PP_LATITUDE_OF_CENTER, 0.0))); |
2656 | |
|
2657 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1, |
2658 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2659 | 0 | SRS_PP_LONGITUDE_OF_CENTER, 0.0))); |
2660 | |
|
2661 | 0 | GTIFKeySet( |
2662 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2663 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2664 | |
|
2665 | 0 | GTIFKeySet( |
2666 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2667 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2668 | 0 | } |
2669 | | |
2670 | 0 | else if (EQUAL(pszProjection, SRS_PT_POLYCONIC)) |
2671 | 0 | { |
2672 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2673 | 0 | ModelTypeProjected); |
2674 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2675 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2676 | |
|
2677 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, CT_Polyconic); |
2678 | |
|
2679 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1, |
2680 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm( |
2681 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0)); |
2682 | |
|
2683 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1, |
2684 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm( |
2685 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0)); |
2686 | |
|
2687 | 0 | GTIFKeySet( |
2688 | 0 | psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1, |
2689 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_SCALE_FACTOR, 1.0)); |
2690 | |
|
2691 | 0 | GTIFKeySet( |
2692 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2693 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2694 | |
|
2695 | 0 | GTIFKeySet( |
2696 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2697 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2698 | 0 | } |
2699 | | |
2700 | 0 | else if (EQUAL(pszProjection, SRS_PT_AZIMUTHAL_EQUIDISTANT)) |
2701 | 0 | { |
2702 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2703 | 0 | ModelTypeProjected); |
2704 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2705 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2706 | |
|
2707 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2708 | 0 | CT_AzimuthalEquidistant); |
2709 | |
|
2710 | 0 | GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1, |
2711 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2712 | 0 | SRS_PP_LATITUDE_OF_CENTER, 0.0))); |
2713 | |
|
2714 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2715 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2716 | 0 | SRS_PP_LONGITUDE_OF_CENTER, 0.0))); |
2717 | |
|
2718 | 0 | GTIFKeySet( |
2719 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2720 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2721 | |
|
2722 | 0 | GTIFKeySet( |
2723 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2724 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2725 | 0 | } |
2726 | | |
2727 | 0 | else if (EQUAL(pszProjection, SRS_PT_MILLER_CYLINDRICAL)) |
2728 | 0 | { |
2729 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2730 | 0 | ModelTypeProjected); |
2731 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2732 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2733 | |
|
2734 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2735 | 0 | CT_MillerCylindrical); |
2736 | |
|
2737 | 0 | GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1, |
2738 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2739 | 0 | SRS_PP_LATITUDE_OF_CENTER, 0.0))); |
2740 | |
|
2741 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2742 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2743 | 0 | SRS_PP_LONGITUDE_OF_CENTER, 0.0))); |
2744 | |
|
2745 | 0 | GTIFKeySet( |
2746 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2747 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2748 | |
|
2749 | 0 | GTIFKeySet( |
2750 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2751 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2752 | 0 | } |
2753 | | |
2754 | 0 | else if (EQUAL(pszProjection, SRS_PT_EQUIRECTANGULAR)) |
2755 | 0 | { |
2756 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2757 | 0 | ModelTypeProjected); |
2758 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2759 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2760 | |
|
2761 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2762 | 0 | CT_Equirectangular); |
2763 | |
|
2764 | 0 | GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1, |
2765 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2766 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2767 | |
|
2768 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2769 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2770 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2771 | |
|
2772 | 0 | GTIFKeySet(psGTIF, ProjStdParallel1GeoKey, TYPE_DOUBLE, 1, |
2773 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2774 | 0 | SRS_PP_STANDARD_PARALLEL_1, 0.0))); |
2775 | |
|
2776 | 0 | GTIFKeySet( |
2777 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2778 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2779 | |
|
2780 | 0 | GTIFKeySet( |
2781 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2782 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2783 | 0 | } |
2784 | | |
2785 | 0 | else if (EQUAL(pszProjection, SRS_PT_GNOMONIC)) |
2786 | 0 | { |
2787 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2788 | 0 | ModelTypeProjected); |
2789 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2790 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2791 | |
|
2792 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, CT_Gnomonic); |
2793 | |
|
2794 | 0 | GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1, |
2795 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2796 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2797 | |
|
2798 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2799 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2800 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2801 | |
|
2802 | 0 | GTIFKeySet( |
2803 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2804 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2805 | |
|
2806 | 0 | GTIFKeySet( |
2807 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2808 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2809 | 0 | } |
2810 | | |
2811 | 0 | else if (EQUAL(pszProjection, SRS_PT_LAMBERT_AZIMUTHAL_EQUAL_AREA)) |
2812 | 0 | { |
2813 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2814 | 0 | ModelTypeProjected); |
2815 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2816 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2817 | |
|
2818 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2819 | 0 | CT_LambertAzimEqualArea); |
2820 | |
|
2821 | 0 | GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1, |
2822 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2823 | 0 | SRS_PP_LATITUDE_OF_CENTER, 0.0))); |
2824 | |
|
2825 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2826 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2827 | 0 | SRS_PP_LONGITUDE_OF_CENTER, 0.0))); |
2828 | |
|
2829 | 0 | GTIFKeySet( |
2830 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2831 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2832 | |
|
2833 | 0 | GTIFKeySet( |
2834 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2835 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2836 | 0 | } |
2837 | | |
2838 | 0 | else if (EQUAL(pszProjection, SRS_PT_ORTHOGRAPHIC)) |
2839 | 0 | { |
2840 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2841 | 0 | ModelTypeProjected); |
2842 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2843 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2844 | |
|
2845 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2846 | 0 | CT_Orthographic); |
2847 | |
|
2848 | 0 | GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1, |
2849 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2850 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2851 | |
|
2852 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2853 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2854 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2855 | |
|
2856 | 0 | GTIFKeySet( |
2857 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2858 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2859 | |
|
2860 | 0 | GTIFKeySet( |
2861 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2862 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2863 | 0 | } |
2864 | | |
2865 | 0 | else if (EQUAL(pszProjection, SRS_PT_NEW_ZEALAND_MAP_GRID)) |
2866 | 0 | { |
2867 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2868 | 0 | ModelTypeProjected); |
2869 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2870 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2871 | |
|
2872 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2873 | 0 | CT_NewZealandMapGrid); |
2874 | |
|
2875 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1, |
2876 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2877 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2878 | |
|
2879 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1, |
2880 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2881 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2882 | |
|
2883 | 0 | GTIFKeySet( |
2884 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2885 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2886 | |
|
2887 | 0 | GTIFKeySet( |
2888 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2889 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2890 | 0 | } |
2891 | | |
2892 | 0 | else if (EQUAL(pszProjection, SRS_PT_ROBINSON)) |
2893 | 0 | { |
2894 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2895 | 0 | ModelTypeProjected); |
2896 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2897 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2898 | |
|
2899 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, CT_Robinson); |
2900 | |
|
2901 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2902 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2903 | 0 | SRS_PP_LONGITUDE_OF_CENTER, 0.0))); |
2904 | |
|
2905 | 0 | GTIFKeySet( |
2906 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2907 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2908 | |
|
2909 | 0 | GTIFKeySet( |
2910 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2911 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2912 | 0 | } |
2913 | | |
2914 | 0 | else if (EQUAL(pszProjection, SRS_PT_SINUSOIDAL)) |
2915 | 0 | { |
2916 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2917 | 0 | ModelTypeProjected); |
2918 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2919 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2920 | |
|
2921 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, CT_Sinusoidal); |
2922 | |
|
2923 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2924 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2925 | 0 | SRS_PP_LONGITUDE_OF_CENTER, 0.0))); |
2926 | |
|
2927 | 0 | GTIFKeySet( |
2928 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2929 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2930 | |
|
2931 | 0 | GTIFKeySet( |
2932 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2933 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2934 | 0 | } |
2935 | | |
2936 | 0 | else if (EQUAL(pszProjection, SRS_PT_VANDERGRINTEN)) |
2937 | 0 | { |
2938 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2939 | 0 | ModelTypeProjected); |
2940 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2941 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2942 | |
|
2943 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2944 | 0 | CT_VanDerGrinten); |
2945 | |
|
2946 | 0 | GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1, |
2947 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2948 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2949 | |
|
2950 | 0 | GTIFKeySet( |
2951 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
2952 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2953 | |
|
2954 | 0 | GTIFKeySet( |
2955 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
2956 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2957 | 0 | } |
2958 | | |
2959 | 0 | else if (EQUAL(pszProjection, SRS_PT_LAMBERT_CONFORMAL_CONIC_2SP)) |
2960 | 0 | { |
2961 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2962 | 0 | ModelTypeProjected); |
2963 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2964 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2965 | |
|
2966 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
2967 | 0 | CT_LambertConfConic_2SP); |
2968 | |
|
2969 | 0 | GTIFKeySet(psGTIF, ProjFalseOriginLatGeoKey, TYPE_DOUBLE, 1, |
2970 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2971 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
2972 | |
|
2973 | 0 | GTIFKeySet(psGTIF, ProjFalseOriginLongGeoKey, TYPE_DOUBLE, 1, |
2974 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2975 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
2976 | |
|
2977 | 0 | GTIFKeySet(psGTIF, ProjStdParallel1GeoKey, TYPE_DOUBLE, 1, |
2978 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2979 | 0 | SRS_PP_STANDARD_PARALLEL_1, 0.0))); |
2980 | |
|
2981 | 0 | GTIFKeySet(psGTIF, ProjStdParallel2GeoKey, TYPE_DOUBLE, 1, |
2982 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
2983 | 0 | SRS_PP_STANDARD_PARALLEL_2, 0.0))); |
2984 | |
|
2985 | 0 | GTIFKeySet( |
2986 | 0 | psGTIF, ProjFalseOriginEastingGeoKey, TYPE_DOUBLE, 1, |
2987 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
2988 | |
|
2989 | 0 | GTIFKeySet( |
2990 | 0 | psGTIF, ProjFalseOriginNorthingGeoKey, TYPE_DOUBLE, 1, |
2991 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
2992 | 0 | } |
2993 | | |
2994 | 0 | else if (EQUAL(pszProjection, SRS_PT_LAMBERT_CONFORMAL_CONIC_1SP)) |
2995 | 0 | { |
2996 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
2997 | 0 | ModelTypeProjected); |
2998 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
2999 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
3000 | |
|
3001 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
3002 | 0 | CT_LambertConfConic_1SP); |
3003 | |
|
3004 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1, |
3005 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
3006 | 0 | SRS_PP_LATITUDE_OF_ORIGIN, 0.0))); |
3007 | |
|
3008 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1, |
3009 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
3010 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
3011 | |
|
3012 | 0 | GTIFKeySet( |
3013 | 0 | psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1, |
3014 | 0 | poSRSCompatibleOfWKT1->GetNormProjParm(SRS_PP_SCALE_FACTOR, 1.0)); |
3015 | |
|
3016 | 0 | GTIFKeySet( |
3017 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
3018 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
3019 | |
|
3020 | 0 | GTIFKeySet( |
3021 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
3022 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
3023 | 0 | } |
3024 | | |
3025 | 0 | else if (EQUAL(pszProjection, SRS_PT_CYLINDRICAL_EQUAL_AREA)) |
3026 | 0 | { |
3027 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, |
3028 | 0 | ModelTypeProjected); |
3029 | 0 | GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
3030 | 0 | GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, KvUserDefined); |
3031 | |
|
3032 | 0 | GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1, |
3033 | 0 | CT_CylindricalEqualArea); |
3034 | |
|
3035 | 0 | GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1, |
3036 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
3037 | 0 | SRS_PP_CENTRAL_MERIDIAN, 0.0))); |
3038 | |
|
3039 | 0 | GTIFKeySet(psGTIF, ProjStdParallel1GeoKey, TYPE_DOUBLE, 1, |
3040 | 0 | ConvertAngularParam(poSRSCompatibleOfWKT1->GetNormProjParm( |
3041 | 0 | SRS_PP_STANDARD_PARALLEL_1, 0.0))); |
3042 | |
|
3043 | 0 | GTIFKeySet( |
3044 | 0 | psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1, |
3045 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_EASTING, 0.0)); |
3046 | |
|
3047 | 0 | GTIFKeySet( |
3048 | 0 | psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1, |
3049 | 0 | poSRSCompatibleOfWKT1->GetProjParm(SRS_PP_FALSE_NORTHING, 0.0)); |
3050 | 0 | } |
3051 | 0 | else |
3052 | 0 | { |
3053 | 0 | bWritePEString = true; |
3054 | 0 | bUnknownProjection = true; |
3055 | 0 | } |
3056 | | |
3057 | | // Note that VERTCS is an ESRI "spelling" of VERT_CS so we assume if |
3058 | | // we find it that we should try to treat this as a PE string. |
3059 | 0 | if (eFlavor == GEOTIFF_KEYS_ESRI_PE || |
3060 | 0 | poSRS->GetAttrValue("VERTCS") != nullptr) |
3061 | 0 | { |
3062 | 0 | bWritePEString = true; |
3063 | 0 | } |
3064 | |
|
3065 | 0 | if (nPCS == KvUserDefined) |
3066 | 0 | { |
3067 | 0 | const char *pszPROJ4Ext = |
3068 | 0 | poSRS->GetExtension("PROJCS", "PROJ4", nullptr); |
3069 | 0 | if (pszPROJ4Ext && |
3070 | 0 | strstr(pszPROJ4Ext, "+proj=merc +a=6378137 +b=6378137")) |
3071 | 0 | { |
3072 | 0 | bWritePEString = true; |
3073 | 0 | } |
3074 | 0 | } |
3075 | |
|
3076 | 0 | bWritePEString &= |
3077 | 0 | CPLTestBool(CPLGetConfigOption("GTIFF_ESRI_CITATION", "YES")); |
3078 | |
|
3079 | 0 | bool peStrStored = false; |
3080 | |
|
3081 | 0 | if (bWritePEString) |
3082 | 0 | { |
3083 | | // Anything we can't map, store as an ESRI PE string with a citation |
3084 | | // key. |
3085 | 0 | char *pszPEString = nullptr; |
3086 | | // We cheat a bit, but if we have a custom_proj4, do not morph to ESRI |
3087 | | // so as to keep the EXTENSION PROJ4 node |
3088 | 0 | const char *const apszOptionsDefault[] = {nullptr}; |
3089 | 0 | const char *const apszOptionsEsri[] = {"FORMAT=WKT1_ESRI", nullptr}; |
3090 | 0 | const char *const *papszOptions = apszOptionsDefault; |
3091 | 0 | if (!(bUnknownProjection && |
3092 | 0 | poSRS->GetExtension("PROJCS", "PROJ4", nullptr) != nullptr) && |
3093 | 0 | !(poSRS->IsProjected() && !poSRS->IsCompound() && |
3094 | 0 | poSRS->GetAxesCount() == 3)) |
3095 | 0 | { |
3096 | 0 | papszOptions = apszOptionsEsri; |
3097 | 0 | } |
3098 | 0 | poSRS->exportToWkt(&pszPEString, papszOptions); |
3099 | 0 | const int peStrLen = static_cast<int>(strlen(pszPEString)); |
3100 | 0 | if (peStrLen > 0) |
3101 | 0 | { |
3102 | 0 | char *outPeStr = static_cast<char *>( |
3103 | 0 | CPLMalloc(peStrLen + strlen("ESRI PE String = ") + 1)); |
3104 | 0 | strcpy(outPeStr, "ESRI PE String = "); |
3105 | 0 | strcat(outPeStr, pszPEString); |
3106 | 0 | oMapAsciiKeys[PCSCitationGeoKey] = outPeStr; |
3107 | 0 | peStrStored = true; |
3108 | 0 | CPLFree(outPeStr); |
3109 | 0 | } |
3110 | 0 | CPLFree(pszPEString); |
3111 | 0 | GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1, KvUserDefined); |
3112 | | |
3113 | | // Not completely sure we really need to imitate ArcGIS to that point |
3114 | | // but that cannot hurt. |
3115 | 0 | if (nPCS == 3857) |
3116 | 0 | { |
3117 | 0 | oMapAsciiKeys[GTCitationGeoKey] = |
3118 | 0 | "PCS Name = WGS_1984_Web_Mercator_Auxiliary_Sphere"; |
3119 | 0 | GTIFKeySet(psGTIF, GeographicTypeGeoKey, TYPE_SHORT, 1, GCS_WGS_84); |
3120 | 0 | GTIFKeySet(psGTIF, GeogSemiMajorAxisGeoKey, TYPE_DOUBLE, 1, |
3121 | 0 | 6378137.0); |
3122 | 0 | GTIFKeySet(psGTIF, GeogInvFlatteningGeoKey, TYPE_DOUBLE, 1, |
3123 | 0 | 298.257223563); |
3124 | 0 | } |
3125 | 0 | } |
3126 | | |
3127 | | /* -------------------------------------------------------------------- */ |
3128 | | /* Is there a false easting/northing set? If so, write out a */ |
3129 | | /* special geokey tag to indicate that GDAL has written these */ |
3130 | | /* with the proper interpretation of the linear units. */ |
3131 | | /* -------------------------------------------------------------------- */ |
3132 | 0 | double dfFE = 0.0; |
3133 | 0 | double dfFN = 0.0; |
3134 | |
|
3135 | 0 | if (eVersion == GEOTIFF_VERSION_1_0 && |
3136 | 0 | (GDALGTIFKeyGetDOUBLE(psGTIF, ProjFalseEastingGeoKey, &dfFE, 0, 1) || |
3137 | 0 | GDALGTIFKeyGetDOUBLE(psGTIF, ProjFalseNorthingGeoKey, &dfFN, 0, 1) || |
3138 | 0 | GDALGTIFKeyGetDOUBLE(psGTIF, ProjFalseOriginEastingGeoKey, &dfFE, 0, |
3139 | 0 | 1) || |
3140 | 0 | GDALGTIFKeyGetDOUBLE(psGTIF, ProjFalseOriginNorthingGeoKey, &dfFN, 0, |
3141 | 0 | 1)) && |
3142 | 0 | (dfFE != 0.0 || dfFN != 0.0) && nUOMLengthCode != 9001) |
3143 | 0 | { |
3144 | 0 | GTIFKeySet(psGTIF, ProjLinearUnitsInterpCorrectGeoKey, TYPE_SHORT, 1, |
3145 | 0 | static_cast<short>(1)); |
3146 | 0 | } |
3147 | | |
3148 | | /* -------------------------------------------------------------------- */ |
3149 | | /* Write linear units information. */ |
3150 | | /* -------------------------------------------------------------------- */ |
3151 | 0 | if (poSRS->IsGeocentric()) |
3152 | 0 | { |
3153 | 0 | GTIFKeySet(psGTIF, GeogLinearUnitsGeoKey, TYPE_SHORT, 1, |
3154 | 0 | nUOMLengthCode); |
3155 | 0 | if (nUOMLengthCode == KvUserDefined) |
3156 | 0 | GTIFKeySet(psGTIF, GeogLinearUnitSizeGeoKey, TYPE_DOUBLE, 1, |
3157 | 0 | dfLinearUOM); |
3158 | 0 | } |
3159 | 0 | else if (!poSRS->IsGeographic() && |
3160 | 0 | (nPCS == KvUserDefined || eVersion == GEOTIFF_VERSION_1_0)) |
3161 | 0 | { |
3162 | 0 | GTIFKeySet(psGTIF, ProjLinearUnitsGeoKey, TYPE_SHORT, 1, |
3163 | 0 | nUOMLengthCode); |
3164 | 0 | if (nUOMLengthCode == KvUserDefined) |
3165 | 0 | GTIFKeySet(psGTIF, ProjLinearUnitSizeGeoKey, TYPE_DOUBLE, 1, |
3166 | 0 | dfLinearUOM); |
3167 | | |
3168 | | // If linear units name is available and user defined, store it as |
3169 | | // citation. |
3170 | 0 | if (!peStrStored && nUOMLengthCode == KvUserDefined && |
3171 | 0 | !osLinearUOMName.empty() && |
3172 | 0 | CPLTestBool(CPLGetConfigOption("GTIFF_ESRI_CITATION", "YES"))) |
3173 | 0 | { |
3174 | 0 | SetLinearUnitCitation(oMapAsciiKeys, osLinearUOMName.c_str()); |
3175 | 0 | } |
3176 | 0 | } |
3177 | | |
3178 | | /* -------------------------------------------------------------------- */ |
3179 | | /* Write angular units. */ |
3180 | | /* -------------------------------------------------------------------- */ |
3181 | |
|
3182 | 0 | if (bHasEllipsoid && |
3183 | 0 | (nGCS == KvUserDefined || eVersion == GEOTIFF_VERSION_1_0)) |
3184 | 0 | { |
3185 | 0 | if (EQUAL(osAngUnitName.c_str(), "Degree")) |
3186 | 0 | GTIFKeySet(psGTIF, GeogAngularUnitsGeoKey, TYPE_SHORT, 1, |
3187 | 0 | Angular_Degree); |
3188 | 0 | else if (EQUAL(osAngUnitName.c_str(), "arc-second")) |
3189 | 0 | GTIFKeySet(psGTIF, GeogAngularUnitsGeoKey, TYPE_SHORT, 1, |
3190 | 0 | Angular_Arc_Second); |
3191 | 0 | else if (EQUAL(osAngUnitName.c_str(), "arc-minute")) |
3192 | 0 | GTIFKeySet(psGTIF, GeogAngularUnitsGeoKey, TYPE_SHORT, 1, |
3193 | 0 | Angular_Arc_Minute); |
3194 | 0 | else if (EQUAL(osAngUnitName.c_str(), "grad")) |
3195 | 0 | GTIFKeySet(psGTIF, GeogAngularUnitsGeoKey, TYPE_SHORT, 1, |
3196 | 0 | Angular_Grad); |
3197 | 0 | else if (EQUAL(osAngUnitName.c_str(), "gon")) |
3198 | 0 | GTIFKeySet(psGTIF, GeogAngularUnitsGeoKey, TYPE_SHORT, 1, |
3199 | 0 | Angular_Gon); |
3200 | 0 | else if (EQUAL(osAngUnitName.c_str(), "radian")) |
3201 | 0 | GTIFKeySet(psGTIF, GeogAngularUnitsGeoKey, TYPE_SHORT, 1, |
3202 | 0 | Angular_Radian); |
3203 | | // else if (EQUAL(osAngUnitName.c_str(), "microradian")) |
3204 | | // GTIFKeySet(psGTIF, GeogAngularUnitsGeoKey, TYPE_SHORT, 1, |
3205 | | // 9109); |
3206 | 0 | else |
3207 | 0 | { |
3208 | | // GeogCitationGeoKey may be rewritten if the gcs is user defined. |
3209 | 0 | oMapAsciiKeys[GeogCitationGeoKey] = osAngUnitName; |
3210 | 0 | GTIFKeySet(psGTIF, GeogAngularUnitSizeGeoKey, TYPE_DOUBLE, 1, |
3211 | 0 | dfAngUnitValue); |
3212 | 0 | } |
3213 | 0 | } |
3214 | | |
3215 | | /* -------------------------------------------------------------------- */ |
3216 | | /* Try to write a citation from the main coordinate system */ |
3217 | | /* name. */ |
3218 | | /* -------------------------------------------------------------------- */ |
3219 | 0 | if (poSRS->GetName() != nullptr && |
3220 | 0 | ((poSRS->IsProjected() && |
3221 | 0 | (nPCS == KvUserDefined || eVersion == GEOTIFF_VERSION_1_0)) || |
3222 | 0 | poSRS->IsCompound() || poSRS->IsLocal() || |
3223 | 0 | (poSRS->IsGeocentric() && |
3224 | 0 | (nGCS == KvUserDefined || eVersion == GEOTIFF_VERSION_1_0)))) |
3225 | 0 | { |
3226 | 0 | if (!(bWritePEString && nPCS == 3857)) |
3227 | 0 | { |
3228 | 0 | oMapAsciiKeys[GTCitationGeoKey] = poSRS->GetName(); |
3229 | 0 | } |
3230 | 0 | } |
3231 | | |
3232 | | /* -------------------------------------------------------------------- */ |
3233 | | /* Try to write a GCS citation. */ |
3234 | | /* -------------------------------------------------------------------- */ |
3235 | 0 | if (nGCS == KvUserDefined || eVersion == GEOTIFF_VERSION_1_0) |
3236 | 0 | { |
3237 | 0 | const OGR_SRSNode *poGCS = poSRS->GetAttrNode("GEOGCS"); |
3238 | |
|
3239 | 0 | if (poGCS != nullptr && poGCS->GetChild(0) != nullptr) |
3240 | 0 | { |
3241 | 0 | oMapAsciiKeys[GeogCitationGeoKey] = poGCS->GetChild(0)->GetValue(); |
3242 | 0 | } |
3243 | 0 | } |
3244 | | |
3245 | | /* -------------------------------------------------------------------- */ |
3246 | | /* Try to identify the GCS/datum, scanning the EPSG datum file for */ |
3247 | | /* a match. */ |
3248 | | /* -------------------------------------------------------------------- */ |
3249 | 0 | if (nPCS == KvUserDefined) |
3250 | 0 | { |
3251 | 0 | if (nGCS == KvUserDefined && poSRS->IsGeographic() && |
3252 | 0 | std::fabs(poSRS->GetAngularUnits() - CPLAtof(SRS_UA_DEGREE_CONV)) < |
3253 | 0 | 1e-9) |
3254 | 0 | { |
3255 | 0 | if (nDatum == Datum_North_American_Datum_1927) |
3256 | 0 | nGCS = GCS_NAD27; |
3257 | 0 | else if (nDatum == Datum_North_American_Datum_1983) |
3258 | 0 | nGCS = GCS_NAD83; |
3259 | 0 | else if (nDatum == Datum_WGS84 || nDatum == DatumE_WGS84) |
3260 | 0 | nGCS = GCS_WGS_84; |
3261 | 0 | } |
3262 | |
|
3263 | 0 | if (nGCS != KvUserDefined) |
3264 | 0 | { |
3265 | 0 | GTIFKeySet(psGTIF, GeographicTypeGeoKey, TYPE_SHORT, 1, nGCS); |
3266 | 0 | } |
3267 | 0 | else if (nDatum != KvUserDefined) |
3268 | 0 | { |
3269 | 0 | GTIFKeySet(psGTIF, GeographicTypeGeoKey, TYPE_SHORT, 1, |
3270 | 0 | KvUserDefined); |
3271 | 0 | GTIFKeySet(psGTIF, GeogGeodeticDatumGeoKey, TYPE_SHORT, 1, nDatum); |
3272 | 0 | } |
3273 | 0 | else if (nSpheroid != KvUserDefined) |
3274 | 0 | { |
3275 | 0 | GTIFKeySet(psGTIF, GeographicTypeGeoKey, TYPE_SHORT, 1, |
3276 | 0 | KvUserDefined); |
3277 | 0 | GTIFKeySet(psGTIF, GeogGeodeticDatumGeoKey, TYPE_SHORT, 1, |
3278 | 0 | KvUserDefined); |
3279 | 0 | GTIFKeySet(psGTIF, GeogEllipsoidGeoKey, TYPE_SHORT, 1, nSpheroid); |
3280 | 0 | } |
3281 | 0 | else if (dfSemiMajor != 0.0) |
3282 | 0 | { |
3283 | 0 | GTIFKeySet(psGTIF, GeographicTypeGeoKey, TYPE_SHORT, 1, |
3284 | 0 | KvUserDefined); |
3285 | 0 | GTIFKeySet(psGTIF, GeogGeodeticDatumGeoKey, TYPE_SHORT, 1, |
3286 | 0 | KvUserDefined); |
3287 | 0 | GTIFKeySet(psGTIF, GeogEllipsoidGeoKey, TYPE_SHORT, 1, |
3288 | 0 | KvUserDefined); |
3289 | 0 | GTIFKeySet(psGTIF, GeogSemiMajorAxisGeoKey, TYPE_DOUBLE, 1, |
3290 | 0 | dfSemiMajor); |
3291 | 0 | if (dfInvFlattening == 0.0) |
3292 | 0 | GTIFKeySet(psGTIF, GeogSemiMinorAxisGeoKey, TYPE_DOUBLE, 1, |
3293 | 0 | dfSemiMajor); |
3294 | 0 | else |
3295 | 0 | GTIFKeySet(psGTIF, GeogInvFlatteningGeoKey, TYPE_DOUBLE, 1, |
3296 | 0 | dfInvFlattening); |
3297 | 0 | } |
3298 | 0 | else if (poSRS->GetAttrValue("DATUM") != nullptr && |
3299 | 0 | strstr(poSRS->GetAttrValue("DATUM"), "unknown") == nullptr && |
3300 | 0 | strstr(poSRS->GetAttrValue("DATUM"), "unnamed") == nullptr) |
3301 | | |
3302 | 0 | { |
3303 | 0 | CPLError(CE_Warning, CPLE_AppDefined, |
3304 | 0 | "Couldn't translate `%s' to a GeoTIFF datum.", |
3305 | 0 | poSRS->GetAttrValue("DATUM")); |
3306 | 0 | } |
3307 | |
|
3308 | 0 | if (nGCS == KvUserDefined || eVersion == GEOTIFF_VERSION_1_0) |
3309 | 0 | { |
3310 | | // Always set InvFlattening if it is available. |
3311 | | // So that it doesn't need to calculate from SemiMinor. |
3312 | 0 | if (dfInvFlattening != 0.0) |
3313 | 0 | GTIFKeySet(psGTIF, GeogInvFlatteningGeoKey, TYPE_DOUBLE, 1, |
3314 | 0 | dfInvFlattening); |
3315 | | // Always set SemiMajor to keep the precision and in case of |
3316 | | // editing. |
3317 | 0 | if (dfSemiMajor != 0.0) |
3318 | 0 | GTIFKeySet(psGTIF, GeogSemiMajorAxisGeoKey, TYPE_DOUBLE, 1, |
3319 | 0 | dfSemiMajor); |
3320 | |
|
3321 | 0 | if (nGCS == KvUserDefined && |
3322 | 0 | CPLTestBool(CPLGetConfigOption("GTIFF_ESRI_CITATION", "YES"))) |
3323 | 0 | { |
3324 | 0 | SetGeogCSCitation(psGTIF, oMapAsciiKeys, poSRS, |
3325 | 0 | osAngUnitName.c_str(), nDatum, nSpheroid); |
3326 | 0 | } |
3327 | 0 | } |
3328 | 0 | } |
3329 | | |
3330 | | /* -------------------------------------------------------------------- */ |
3331 | | /* Do we have TOWGS84 parameters? */ |
3332 | | /* -------------------------------------------------------------------- */ |
3333 | 0 | #if !defined(GEO_NORMALIZE_DISABLE_TOWGS84) |
3334 | 0 | double adfTOWGS84[7] = {0.0}; |
3335 | |
|
3336 | 0 | if ((nGCS == KvUserDefined || eVersion == GEOTIFF_VERSION_1_0) && |
3337 | 0 | poSRS->GetTOWGS84(adfTOWGS84) == OGRERR_NONE) |
3338 | 0 | { |
3339 | | // If we are writing a SRS with a EPSG code, and that the EPSG code |
3340 | | // of the current SRS object and the one coming from the EPSG code |
3341 | | // are the same, then by default, do not write them. |
3342 | 0 | bool bUseReferenceTOWGS84 = false; |
3343 | 0 | const char *pszAuthName = poSRS->GetAuthorityName(nullptr); |
3344 | 0 | const char *pszAuthCode = poSRS->GetAuthorityCode(nullptr); |
3345 | 0 | if (pszAuthName && EQUAL(pszAuthName, "EPSG") && pszAuthCode) |
3346 | 0 | { |
3347 | 0 | CPLErrorHandlerPusher oErrorHandler(CPLQuietErrorHandler); |
3348 | 0 | OGRSpatialReference oRefSRS; |
3349 | 0 | double adfRefTOWGS84[7] = {0.0}; |
3350 | 0 | if (oRefSRS.importFromEPSG(atoi(pszAuthCode)) == OGRERR_NONE) |
3351 | 0 | { |
3352 | 0 | oRefSRS.AddGuessedTOWGS84(); |
3353 | 0 | if (oRefSRS.GetTOWGS84(adfRefTOWGS84) == OGRERR_NONE && |
3354 | 0 | memcmp(adfRefTOWGS84, adfTOWGS84, sizeof(adfTOWGS84)) == 0) |
3355 | 0 | { |
3356 | 0 | bUseReferenceTOWGS84 = true; |
3357 | 0 | } |
3358 | 0 | } |
3359 | 0 | } |
3360 | 0 | const char *pszWriteTOWGS84 = |
3361 | 0 | CPLGetConfigOption("GTIFF_WRITE_TOWGS84", "AUTO"); |
3362 | 0 | if ((EQUAL(pszWriteTOWGS84, "YES") || EQUAL(pszWriteTOWGS84, "TRUE") || |
3363 | 0 | EQUAL(pszWriteTOWGS84, "ON")) || |
3364 | 0 | (!bUseReferenceTOWGS84 && EQUAL(pszWriteTOWGS84, "AUTO"))) |
3365 | 0 | { |
3366 | 0 | if (adfTOWGS84[3] == 0.0 && adfTOWGS84[4] == 0.0 && |
3367 | 0 | adfTOWGS84[5] == 0.0 && adfTOWGS84[6] == 0.0) |
3368 | 0 | { |
3369 | 0 | if (nGCS == GCS_WGS_84 && adfTOWGS84[0] == 0.0 && |
3370 | 0 | adfTOWGS84[1] == 0.0 && adfTOWGS84[2] == 0.0) |
3371 | 0 | { |
3372 | 0 | ; // Do nothing. |
3373 | 0 | } |
3374 | 0 | else |
3375 | 0 | GTIFKeySet(psGTIF, GeogTOWGS84GeoKey, TYPE_DOUBLE, 3, |
3376 | 0 | adfTOWGS84); |
3377 | 0 | } |
3378 | 0 | else |
3379 | 0 | GTIFKeySet(psGTIF, GeogTOWGS84GeoKey, TYPE_DOUBLE, 7, |
3380 | 0 | adfTOWGS84); |
3381 | 0 | } |
3382 | 0 | } |
3383 | 0 | #endif |
3384 | | |
3385 | | /* -------------------------------------------------------------------- */ |
3386 | | /* Do we have vertical information to set? */ |
3387 | | /* -------------------------------------------------------------------- */ |
3388 | 0 | if (poSRS->GetAttrValue("COMPD_CS|VERT_CS") != nullptr) |
3389 | 0 | { |
3390 | 0 | bool bGotVertCSCode = false; |
3391 | 0 | const char *pszVertCSCode = poSRS->GetAuthorityCode("COMPD_CS|VERT_CS"); |
3392 | 0 | const char *pszVertCSAuthName = |
3393 | 0 | poSRS->GetAuthorityName("COMPD_CS|VERT_CS"); |
3394 | 0 | if (pszVertCSCode && pszVertCSAuthName && atoi(pszVertCSCode) && |
3395 | 0 | EQUAL(pszVertCSAuthName, "EPSG")) |
3396 | 0 | { |
3397 | 0 | bGotVertCSCode = true; |
3398 | 0 | GTIFKeySet(psGTIF, VerticalCSTypeGeoKey, TYPE_SHORT, 1, |
3399 | 0 | atoi(pszVertCSCode)); |
3400 | 0 | } |
3401 | 0 | else if (eVersion >= GEOTIFF_VERSION_1_1) |
3402 | 0 | { |
3403 | 0 | GTIFKeySet(psGTIF, VerticalCSTypeGeoKey, TYPE_SHORT, 1, |
3404 | 0 | KvUserDefined); |
3405 | 0 | } |
3406 | |
|
3407 | 0 | if (eVersion == GEOTIFF_VERSION_1_0 || !bGotVertCSCode) |
3408 | 0 | { |
3409 | 0 | oMapAsciiKeys[VerticalCitationGeoKey] = |
3410 | 0 | poSRS->GetAttrValue("COMPD_CS|VERT_CS"); |
3411 | |
|
3412 | 0 | const char *pszVertDatumCode = |
3413 | 0 | poSRS->GetAuthorityCode("COMPD_CS|VERT_CS|VERT_DATUM"); |
3414 | 0 | const char *pszVertDatumAuthName = |
3415 | 0 | poSRS->GetAuthorityName("COMPD_CS|VERT_CS|VERT_DATUM"); |
3416 | 0 | if (pszVertDatumCode && pszVertDatumAuthName && |
3417 | 0 | atoi(pszVertDatumCode) && EQUAL(pszVertDatumAuthName, "EPSG")) |
3418 | 0 | { |
3419 | 0 | GTIFKeySet(psGTIF, VerticalDatumGeoKey, TYPE_SHORT, 1, |
3420 | 0 | atoi(pszVertDatumCode)); |
3421 | 0 | } |
3422 | 0 | else if (eVersion >= GEOTIFF_VERSION_1_1) |
3423 | 0 | { |
3424 | 0 | GTIFKeySet(psGTIF, VerticalDatumGeoKey, TYPE_SHORT, 1, |
3425 | 0 | KvUserDefined); |
3426 | 0 | } |
3427 | |
|
3428 | 0 | const char *pszVertUnitCode = |
3429 | 0 | poSRS->GetAuthorityCode("COMPD_CS|VERT_CS|UNIT"); |
3430 | 0 | const char *pszVertUnitAuthName = |
3431 | 0 | poSRS->GetAuthorityName("COMPD_CS|VERT_CS|UNIT"); |
3432 | 0 | if (pszVertUnitCode && pszVertUnitAuthName && |
3433 | 0 | atoi(pszVertUnitCode) && EQUAL(pszVertUnitAuthName, "EPSG")) |
3434 | 0 | { |
3435 | 0 | GTIFKeySet(psGTIF, VerticalUnitsGeoKey, TYPE_SHORT, 1, |
3436 | 0 | atoi(pszVertUnitCode)); |
3437 | 0 | } |
3438 | 0 | else if (eVersion >= GEOTIFF_VERSION_1_1) |
3439 | 0 | { |
3440 | 0 | GTIFKeySet(psGTIF, VerticalUnitsGeoKey, TYPE_SHORT, 1, |
3441 | 0 | KvUserDefined); |
3442 | 0 | } |
3443 | 0 | } |
3444 | 0 | } |
3445 | 0 | else if (eVersion >= GEOTIFF_VERSION_1_1 && nVerticalCSKeyValue != 0) |
3446 | 0 | { |
3447 | 0 | GTIFKeySet(psGTIF, VerticalCSTypeGeoKey, TYPE_SHORT, 1, |
3448 | 0 | nVerticalCSKeyValue); |
3449 | 0 | } |
3450 | |
|
3451 | 0 | const double dfCoordinateEpoch = poSRS->GetCoordinateEpoch(); |
3452 | 0 | if (dfCoordinateEpoch > 0) |
3453 | 0 | { |
3454 | 0 | GTIFKeySet(psGTIF, CoordinateEpochGeoKey, TYPE_DOUBLE, 1, |
3455 | 0 | dfCoordinateEpoch); |
3456 | 0 | } |
3457 | | |
3458 | | /* -------------------------------------------------------------------- */ |
3459 | | /* Write all ascii keys */ |
3460 | | /* -------------------------------------------------------------------- */ |
3461 | 0 | for (const auto &oIter : oMapAsciiKeys) |
3462 | 0 | { |
3463 | 0 | GTIFKeySet(psGTIF, oIter.first, TYPE_ASCII, 0, oIter.second.c_str()); |
3464 | 0 | } |
3465 | |
|
3466 | 0 | return TRUE; |
3467 | 0 | } |
3468 | | |
3469 | | /************************************************************************/ |
3470 | | /* GTIFWktFromMemBuf() */ |
3471 | | /************************************************************************/ |
3472 | | |
3473 | | CPLErr GTIFWktFromMemBuf(int nSize, unsigned char *pabyBuffer, char **ppszWKT, |
3474 | | double *padfGeoTransform, int *pnGCPCount, |
3475 | | GDAL_GCP **ppasGCPList) |
3476 | 0 | { |
3477 | 0 | OGRSpatialReferenceH hSRS = nullptr; |
3478 | 0 | if (ppszWKT) |
3479 | 0 | *ppszWKT = nullptr; |
3480 | 0 | CPLErr eErr = |
3481 | 0 | GTIFWktFromMemBufEx(nSize, pabyBuffer, &hSRS, padfGeoTransform, |
3482 | 0 | pnGCPCount, ppasGCPList, nullptr, nullptr); |
3483 | 0 | if (eErr == CE_None) |
3484 | 0 | { |
3485 | 0 | if (hSRS && ppszWKT) |
3486 | 0 | { |
3487 | 0 | OSRExportToWkt(hSRS, ppszWKT); |
3488 | 0 | } |
3489 | 0 | } |
3490 | 0 | OSRDestroySpatialReference(hSRS); |
3491 | 0 | return eErr; |
3492 | 0 | } |
3493 | | |
3494 | | CPLErr GTIFWktFromMemBufEx(int nSize, unsigned char *pabyBuffer, |
3495 | | OGRSpatialReferenceH *phSRS, |
3496 | | double *padfGeoTransform, int *pnGCPCount, |
3497 | | GDAL_GCP **ppasGCPList, int *pbPixelIsPoint, |
3498 | | char ***ppapszRPCMD) |
3499 | | |
3500 | 0 | { |
3501 | 0 | const std::string osFilename( |
3502 | 0 | VSIMemGenerateHiddenFilename("wkt_from_mem_buf.tif")); |
3503 | | |
3504 | | /* -------------------------------------------------------------------- */ |
3505 | | /* Initialization of libtiff and libgeotiff. */ |
3506 | | /* -------------------------------------------------------------------- */ |
3507 | 0 | GTiffOneTimeInit(); // For RPC tag. |
3508 | 0 | LibgeotiffOneTimeInit(); |
3509 | | |
3510 | | /* -------------------------------------------------------------------- */ |
3511 | | /* Create a memory file from the buffer. */ |
3512 | | /* -------------------------------------------------------------------- */ |
3513 | 0 | VSILFILE *fp = |
3514 | 0 | VSIFileFromMemBuffer(osFilename.c_str(), pabyBuffer, nSize, FALSE); |
3515 | 0 | if (fp == nullptr) |
3516 | 0 | return CE_Failure; |
3517 | | |
3518 | | /* -------------------------------------------------------------------- */ |
3519 | | /* Initialize access to the memory geotiff structure. */ |
3520 | | /* -------------------------------------------------------------------- */ |
3521 | 0 | TIFF *hTIFF = VSI_TIFFOpen(osFilename.c_str(), "rc", fp); |
3522 | |
|
3523 | 0 | if (hTIFF == nullptr) |
3524 | 0 | { |
3525 | 0 | CPLError(CE_Failure, CPLE_AppDefined, |
3526 | 0 | "TIFF/GeoTIFF structure is corrupt."); |
3527 | 0 | VSIUnlink(osFilename.c_str()); |
3528 | 0 | CPL_IGNORE_RET_VAL(VSIFCloseL(fp)); |
3529 | 0 | return CE_Failure; |
3530 | 0 | } |
3531 | | |
3532 | | /* -------------------------------------------------------------------- */ |
3533 | | /* Get the projection definition. */ |
3534 | | /* -------------------------------------------------------------------- */ |
3535 | 0 | bool bPixelIsPoint = false; |
3536 | 0 | bool bPointGeoIgnore = false; |
3537 | 0 | unsigned short nRasterType = 0; |
3538 | |
|
3539 | 0 | GTIF *hGTIF = GTIFNew(hTIFF); |
3540 | 0 | if (hGTIF) |
3541 | 0 | GTIFAttachPROJContext(hGTIF, OSRGetProjTLSContext()); |
3542 | |
|
3543 | 0 | if (hGTIF != nullptr && |
3544 | 0 | GDALGTIFKeyGetSHORT(hGTIF, GTRasterTypeGeoKey, &nRasterType, 0, 1) == |
3545 | 0 | 1 && |
3546 | 0 | nRasterType == static_cast<unsigned short>(RasterPixelIsPoint)) |
3547 | 0 | { |
3548 | 0 | bPixelIsPoint = true; |
3549 | 0 | bPointGeoIgnore = |
3550 | 0 | CPLTestBool(CPLGetConfigOption("GTIFF_POINT_GEO_IGNORE", "FALSE")); |
3551 | 0 | } |
3552 | 0 | if (pbPixelIsPoint) |
3553 | 0 | *pbPixelIsPoint = bPixelIsPoint; |
3554 | 0 | if (ppapszRPCMD) |
3555 | 0 | *ppapszRPCMD = nullptr; |
3556 | |
|
3557 | 0 | if (phSRS) |
3558 | 0 | { |
3559 | 0 | *phSRS = nullptr; |
3560 | 0 | if (hGTIF != nullptr) |
3561 | 0 | { |
3562 | 0 | GTIFDefn *psGTIFDefn = GTIFAllocDefn(); |
3563 | 0 | if (GTIFGetDefn(hGTIF, psGTIFDefn)) |
3564 | 0 | { |
3565 | 0 | *phSRS = GTIFGetOGISDefnAsOSR(hGTIF, psGTIFDefn); |
3566 | 0 | } |
3567 | 0 | GTIFFreeDefn(psGTIFDefn); |
3568 | 0 | } |
3569 | 0 | } |
3570 | 0 | if (hGTIF) |
3571 | 0 | GTIFFree(hGTIF); |
3572 | | |
3573 | | /* -------------------------------------------------------------------- */ |
3574 | | /* Get geotransform or tiepoints. */ |
3575 | | /* -------------------------------------------------------------------- */ |
3576 | 0 | double *padfTiePoints = nullptr; |
3577 | 0 | double *padfScale = nullptr; |
3578 | 0 | double *padfMatrix = nullptr; |
3579 | 0 | int16_t nCount = 0; |
3580 | |
|
3581 | 0 | padfGeoTransform[0] = 0.0; |
3582 | 0 | padfGeoTransform[1] = 1.0; |
3583 | 0 | padfGeoTransform[2] = 0.0; |
3584 | 0 | padfGeoTransform[3] = 0.0; |
3585 | 0 | padfGeoTransform[4] = 0.0; |
3586 | 0 | padfGeoTransform[5] = 1.0; |
3587 | |
|
3588 | 0 | *pnGCPCount = 0; |
3589 | 0 | *ppasGCPList = nullptr; |
3590 | |
|
3591 | 0 | if (TIFFGetField(hTIFF, TIFFTAG_GEOPIXELSCALE, &nCount, &padfScale) && |
3592 | 0 | nCount >= 2) |
3593 | 0 | { |
3594 | 0 | padfGeoTransform[1] = padfScale[0]; |
3595 | 0 | padfGeoTransform[5] = -std::abs(padfScale[1]); |
3596 | |
|
3597 | 0 | if (TIFFGetField(hTIFF, TIFFTAG_GEOTIEPOINTS, &nCount, |
3598 | 0 | &padfTiePoints) && |
3599 | 0 | nCount >= 6) |
3600 | 0 | { |
3601 | 0 | padfGeoTransform[0] = |
3602 | 0 | padfTiePoints[3] - padfTiePoints[0] * padfGeoTransform[1]; |
3603 | 0 | padfGeoTransform[3] = |
3604 | 0 | padfTiePoints[4] - padfTiePoints[1] * padfGeoTransform[5]; |
3605 | | |
3606 | | // Adjust for pixel is point in transform. |
3607 | 0 | if (bPixelIsPoint && !bPointGeoIgnore) |
3608 | 0 | { |
3609 | 0 | padfGeoTransform[0] -= |
3610 | 0 | padfGeoTransform[1] * 0.5 + padfGeoTransform[2] * 0.5; |
3611 | 0 | padfGeoTransform[3] -= |
3612 | 0 | padfGeoTransform[4] * 0.5 + padfGeoTransform[5] * 0.5; |
3613 | 0 | } |
3614 | 0 | } |
3615 | 0 | } |
3616 | 0 | else if (TIFFGetField(hTIFF, TIFFTAG_GEOTIEPOINTS, &nCount, |
3617 | 0 | &padfTiePoints) && |
3618 | 0 | nCount >= 6) |
3619 | 0 | { |
3620 | 0 | *pnGCPCount = nCount / 6; |
3621 | 0 | *ppasGCPList = |
3622 | 0 | static_cast<GDAL_GCP *>(CPLCalloc(sizeof(GDAL_GCP), *pnGCPCount)); |
3623 | |
|
3624 | 0 | for (int iGCP = 0; iGCP < *pnGCPCount; iGCP++) |
3625 | 0 | { |
3626 | 0 | char szID[32] = {}; |
3627 | 0 | GDAL_GCP *psGCP = *ppasGCPList + iGCP; |
3628 | |
|
3629 | 0 | snprintf(szID, sizeof(szID), "%d", iGCP + 1); |
3630 | 0 | psGCP->pszId = CPLStrdup(szID); |
3631 | 0 | psGCP->pszInfo = CPLStrdup(""); |
3632 | 0 | psGCP->dfGCPPixel = padfTiePoints[iGCP * 6 + 0]; |
3633 | 0 | psGCP->dfGCPLine = padfTiePoints[iGCP * 6 + 1]; |
3634 | 0 | psGCP->dfGCPX = padfTiePoints[iGCP * 6 + 3]; |
3635 | 0 | psGCP->dfGCPY = padfTiePoints[iGCP * 6 + 4]; |
3636 | 0 | psGCP->dfGCPZ = padfTiePoints[iGCP * 6 + 5]; |
3637 | 0 | } |
3638 | 0 | } |
3639 | 0 | else if (TIFFGetField(hTIFF, TIFFTAG_GEOTRANSMATRIX, &nCount, |
3640 | 0 | &padfMatrix) && |
3641 | 0 | nCount == 16) |
3642 | 0 | { |
3643 | 0 | padfGeoTransform[0] = padfMatrix[3]; |
3644 | 0 | padfGeoTransform[1] = padfMatrix[0]; |
3645 | 0 | padfGeoTransform[2] = padfMatrix[1]; |
3646 | 0 | padfGeoTransform[3] = padfMatrix[7]; |
3647 | 0 | padfGeoTransform[4] = padfMatrix[4]; |
3648 | 0 | padfGeoTransform[5] = padfMatrix[5]; |
3649 | 0 | } |
3650 | | |
3651 | | /* -------------------------------------------------------------------- */ |
3652 | | /* Read RPC */ |
3653 | | /* -------------------------------------------------------------------- */ |
3654 | 0 | if (ppapszRPCMD != nullptr) |
3655 | 0 | { |
3656 | 0 | *ppapszRPCMD = GTiffDatasetReadRPCTag(hTIFF); |
3657 | 0 | } |
3658 | | |
3659 | | /* -------------------------------------------------------------------- */ |
3660 | | /* Cleanup. */ |
3661 | | /* -------------------------------------------------------------------- */ |
3662 | 0 | XTIFFClose(hTIFF); |
3663 | 0 | CPL_IGNORE_RET_VAL(VSIFCloseL(fp)); |
3664 | |
|
3665 | 0 | VSIUnlink(osFilename.c_str()); |
3666 | |
|
3667 | 0 | if (phSRS && *phSRS == nullptr) |
3668 | 0 | return CE_Failure; |
3669 | | |
3670 | 0 | return CE_None; |
3671 | 0 | } |
3672 | | |
3673 | | /************************************************************************/ |
3674 | | /* GTIFMemBufFromWkt() */ |
3675 | | /************************************************************************/ |
3676 | | |
3677 | | CPLErr GTIFMemBufFromWkt(const char *pszWKT, const double *padfGeoTransform, |
3678 | | int nGCPCount, const GDAL_GCP *pasGCPList, int *pnSize, |
3679 | | unsigned char **ppabyBuffer) |
3680 | 0 | { |
3681 | 0 | OGRSpatialReference oSRS; |
3682 | 0 | if (pszWKT != nullptr) |
3683 | 0 | oSRS.importFromWkt(pszWKT); |
3684 | 0 | return GTIFMemBufFromSRS(OGRSpatialReference::ToHandle(&oSRS), |
3685 | 0 | padfGeoTransform, nGCPCount, pasGCPList, pnSize, |
3686 | 0 | ppabyBuffer, FALSE, nullptr); |
3687 | 0 | } |
3688 | | |
3689 | | CPLErr GTIFMemBufFromSRS(OGRSpatialReferenceH hSRS, |
3690 | | const double *padfGeoTransform, int nGCPCount, |
3691 | | const GDAL_GCP *pasGCPList, int *pnSize, |
3692 | | unsigned char **ppabyBuffer, int bPixelIsPoint, |
3693 | | char **papszRPCMD) |
3694 | | |
3695 | 0 | { |
3696 | 0 | const std::string osFilename( |
3697 | 0 | VSIMemGenerateHiddenFilename("wkt_from_mem_buf.tif")); |
3698 | | |
3699 | | /* -------------------------------------------------------------------- */ |
3700 | | /* Initialization of libtiff and libgeotiff. */ |
3701 | | /* -------------------------------------------------------------------- */ |
3702 | 0 | GTiffOneTimeInit(); // For RPC tag. |
3703 | 0 | LibgeotiffOneTimeInit(); |
3704 | | |
3705 | | /* -------------------------------------------------------------------- */ |
3706 | | /* Initialize access to the memory geotiff structure. */ |
3707 | | /* -------------------------------------------------------------------- */ |
3708 | 0 | VSILFILE *fpL = VSIFOpenL(osFilename.c_str(), "w"); |
3709 | 0 | if (fpL == nullptr) |
3710 | 0 | return CE_Failure; |
3711 | | |
3712 | 0 | TIFF *hTIFF = VSI_TIFFOpen(osFilename.c_str(), "w", fpL); |
3713 | |
|
3714 | 0 | if (hTIFF == nullptr) |
3715 | 0 | { |
3716 | 0 | CPLError(CE_Failure, CPLE_AppDefined, |
3717 | 0 | "TIFF/GeoTIFF structure is corrupt."); |
3718 | 0 | CPL_IGNORE_RET_VAL(VSIFCloseL(fpL)); |
3719 | 0 | VSIUnlink(osFilename.c_str()); |
3720 | 0 | return CE_Failure; |
3721 | 0 | } |
3722 | | |
3723 | | /* -------------------------------------------------------------------- */ |
3724 | | /* Write some minimal set of image parameters. */ |
3725 | | /* -------------------------------------------------------------------- */ |
3726 | 0 | TIFFSetField(hTIFF, TIFFTAG_IMAGEWIDTH, 1); |
3727 | 0 | TIFFSetField(hTIFF, TIFFTAG_IMAGELENGTH, 1); |
3728 | 0 | TIFFSetField(hTIFF, TIFFTAG_BITSPERSAMPLE, 8); |
3729 | 0 | TIFFSetField(hTIFF, TIFFTAG_SAMPLESPERPIXEL, 1); |
3730 | 0 | TIFFSetField(hTIFF, TIFFTAG_ROWSPERSTRIP, 1); |
3731 | 0 | TIFFSetField(hTIFF, TIFFTAG_PLANARCONFIG, PLANARCONFIG_CONTIG); |
3732 | 0 | TIFFSetField(hTIFF, TIFFTAG_PHOTOMETRIC, PHOTOMETRIC_MINISBLACK); |
3733 | | |
3734 | | /* -------------------------------------------------------------------- */ |
3735 | | /* Get the projection definition. */ |
3736 | | /* -------------------------------------------------------------------- */ |
3737 | |
|
3738 | 0 | bool bPointGeoIgnore = false; |
3739 | 0 | if (bPixelIsPoint) |
3740 | 0 | { |
3741 | 0 | bPointGeoIgnore = |
3742 | 0 | CPLTestBool(CPLGetConfigOption("GTIFF_POINT_GEO_IGNORE", "FALSE")); |
3743 | 0 | } |
3744 | |
|
3745 | 0 | GTIF *hGTIF = nullptr; |
3746 | 0 | if (hSRS != nullptr || bPixelIsPoint) |
3747 | 0 | { |
3748 | 0 | hGTIF = GTIFNew(hTIFF); |
3749 | 0 | if (hGTIF) |
3750 | 0 | { |
3751 | 0 | GTIFAttachPROJContext(hGTIF, OSRGetProjTLSContext()); |
3752 | |
|
3753 | 0 | if (hSRS != nullptr) |
3754 | 0 | GTIFSetFromOGISDefnEx(hGTIF, hSRS, GEOTIFF_KEYS_STANDARD, |
3755 | 0 | GEOTIFF_VERSION_1_0); |
3756 | |
|
3757 | 0 | if (bPixelIsPoint) |
3758 | 0 | { |
3759 | 0 | GTIFKeySet(hGTIF, GTRasterTypeGeoKey, TYPE_SHORT, 1, |
3760 | 0 | RasterPixelIsPoint); |
3761 | 0 | } |
3762 | |
|
3763 | 0 | GTIFWriteKeys(hGTIF); |
3764 | 0 | GTIFFree(hGTIF); |
3765 | 0 | } |
3766 | 0 | } |
3767 | | |
3768 | | /* -------------------------------------------------------------------- */ |
3769 | | /* Set the geotransform, or GCPs. */ |
3770 | | /* -------------------------------------------------------------------- */ |
3771 | |
|
3772 | 0 | if (padfGeoTransform[0] != 0.0 || padfGeoTransform[1] != 1.0 || |
3773 | 0 | padfGeoTransform[2] != 0.0 || padfGeoTransform[3] != 0.0 || |
3774 | 0 | padfGeoTransform[4] != 0.0 || std::abs(padfGeoTransform[5]) != 1.0) |
3775 | 0 | { |
3776 | |
|
3777 | 0 | if (padfGeoTransform[2] == 0.0 && padfGeoTransform[4] == 0.0) |
3778 | 0 | { |
3779 | 0 | double adfPixelScale[3] = {padfGeoTransform[1], |
3780 | 0 | fabs(padfGeoTransform[5]), 0.0}; |
3781 | |
|
3782 | 0 | TIFFSetField(hTIFF, TIFFTAG_GEOPIXELSCALE, 3, adfPixelScale); |
3783 | |
|
3784 | 0 | double adfTiePoints[6] = { |
3785 | 0 | 0.0, 0.0, 0.0, padfGeoTransform[0], padfGeoTransform[3], 0.0}; |
3786 | |
|
3787 | 0 | if (bPixelIsPoint && !bPointGeoIgnore) |
3788 | 0 | { |
3789 | 0 | adfTiePoints[3] += |
3790 | 0 | padfGeoTransform[1] * 0.5 + padfGeoTransform[2] * 0.5; |
3791 | 0 | adfTiePoints[4] += |
3792 | 0 | padfGeoTransform[4] * 0.5 + padfGeoTransform[5] * 0.5; |
3793 | 0 | } |
3794 | |
|
3795 | 0 | TIFFSetField(hTIFF, TIFFTAG_GEOTIEPOINTS, 6, adfTiePoints); |
3796 | 0 | } |
3797 | 0 | else |
3798 | 0 | { |
3799 | 0 | double adfMatrix[16] = {0.0}; |
3800 | |
|
3801 | 0 | adfMatrix[0] = padfGeoTransform[1]; |
3802 | 0 | adfMatrix[1] = padfGeoTransform[2]; |
3803 | 0 | adfMatrix[3] = padfGeoTransform[0]; |
3804 | 0 | adfMatrix[4] = padfGeoTransform[4]; |
3805 | 0 | adfMatrix[5] = padfGeoTransform[5]; |
3806 | 0 | adfMatrix[7] = padfGeoTransform[3]; |
3807 | 0 | adfMatrix[15] = 1.0; |
3808 | |
|
3809 | 0 | if (bPixelIsPoint && !bPointGeoIgnore) |
3810 | 0 | { |
3811 | 0 | adfMatrix[3] += |
3812 | 0 | padfGeoTransform[1] * 0.5 + padfGeoTransform[2] * 0.5; |
3813 | 0 | adfMatrix[7] += |
3814 | 0 | padfGeoTransform[4] * 0.5 + padfGeoTransform[5] * 0.5; |
3815 | 0 | } |
3816 | |
|
3817 | 0 | TIFFSetField(hTIFF, TIFFTAG_GEOTRANSMATRIX, 16, adfMatrix); |
3818 | 0 | } |
3819 | 0 | } |
3820 | | |
3821 | | /* -------------------------------------------------------------------- */ |
3822 | | /* Otherwise write tiepoints if they are available. */ |
3823 | | /* -------------------------------------------------------------------- */ |
3824 | 0 | else if (nGCPCount > 0) |
3825 | 0 | { |
3826 | 0 | double *padfTiePoints = |
3827 | 0 | static_cast<double *>(CPLMalloc(6 * sizeof(double) * nGCPCount)); |
3828 | |
|
3829 | 0 | for (int iGCP = 0; iGCP < nGCPCount; iGCP++) |
3830 | 0 | { |
3831 | |
|
3832 | 0 | padfTiePoints[iGCP * 6 + 0] = pasGCPList[iGCP].dfGCPPixel; |
3833 | 0 | padfTiePoints[iGCP * 6 + 1] = pasGCPList[iGCP].dfGCPLine; |
3834 | 0 | padfTiePoints[iGCP * 6 + 2] = 0; |
3835 | 0 | padfTiePoints[iGCP * 6 + 3] = pasGCPList[iGCP].dfGCPX; |
3836 | 0 | padfTiePoints[iGCP * 6 + 4] = pasGCPList[iGCP].dfGCPY; |
3837 | 0 | padfTiePoints[iGCP * 6 + 5] = pasGCPList[iGCP].dfGCPZ; |
3838 | 0 | } |
3839 | |
|
3840 | 0 | TIFFSetField(hTIFF, TIFFTAG_GEOTIEPOINTS, 6 * nGCPCount, padfTiePoints); |
3841 | 0 | CPLFree(padfTiePoints); |
3842 | 0 | } |
3843 | | |
3844 | | /* -------------------------------------------------------------------- */ |
3845 | | /* Write RPC */ |
3846 | | /* -------------------------------------------------------------------- */ |
3847 | 0 | if (papszRPCMD != nullptr) |
3848 | 0 | { |
3849 | 0 | GTiffDatasetWriteRPCTag(hTIFF, papszRPCMD); |
3850 | 0 | } |
3851 | | |
3852 | | /* -------------------------------------------------------------------- */ |
3853 | | /* Cleanup and return the created memory buffer. */ |
3854 | | /* -------------------------------------------------------------------- */ |
3855 | 0 | GByte bySmallImage = 0; |
3856 | |
|
3857 | 0 | TIFFWriteEncodedStrip(hTIFF, 0, reinterpret_cast<char *>(&bySmallImage), 1); |
3858 | 0 | TIFFWriteCheck(hTIFF, TIFFIsTiled(hTIFF), "GTIFMemBufFromWkt"); |
3859 | 0 | TIFFWriteDirectory(hTIFF); |
3860 | |
|
3861 | 0 | XTIFFClose(hTIFF); |
3862 | 0 | CPL_IGNORE_RET_VAL(VSIFCloseL(fpL)); |
3863 | | |
3864 | | /* -------------------------------------------------------------------- */ |
3865 | | /* Read back from the memory buffer. It would be preferable */ |
3866 | | /* to be able to "steal" the memory buffer, but there isn't */ |
3867 | | /* currently any support for this. */ |
3868 | | /* -------------------------------------------------------------------- */ |
3869 | 0 | GUIntBig nBigLength = 0; |
3870 | |
|
3871 | 0 | *ppabyBuffer = VSIGetMemFileBuffer(osFilename.c_str(), &nBigLength, TRUE); |
3872 | 0 | *pnSize = static_cast<int>(nBigLength); |
3873 | |
|
3874 | 0 | return CE_None; |
3875 | 0 | } |