1 /******************************************************************************
2 * $Id: gt_wkt_srs.cpp 27182 2014-04-14 20:03:08Z rouault $
3 *
4 * Project: GeoTIFF Driver
5 * Purpose: Implements translation between GeoTIFF normalized projection
6 * definitions and OpenGIS WKT SRS format. This code is intended to
7 * be moved into libgeotiff someday if possible.
8 * Author: Frank Warmerdam, warmerdam@pobox.com
9 *
10 ******************************************************************************
11 * Copyright (c) 1999, Frank Warmerdam
12 * Copyright (c) 2008-2014, Even Rouault <even dot rouault at mines-paris dot org>
13 *
14 * Permission is hereby granted, free of charge, to any person obtaining a
15 * copy of this software and associated documentation files (the "Software"),
16 * to deal in the Software without restriction, including without limitation
17 * the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 * and/or sell copies of the Software, and to permit persons to whom the
19 * Software is furnished to do so, subject to the following conditions:
20 *
21 * The above copyright notice and this permission notice shall be included
22 * in all copies or substantial portions of the Software.
23 *
24 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 * DEALINGS IN THE SOFTWARE.
31 ****************************************************************************/
32
33 #include "cpl_error.h"
34 #include "cpl_conv.h"
35 #include "cpl_csv.h"
36 #include "gdal_csv.h"
37
38 #include "geovalues.h"
39 #include "ogr_spatialref.h"
40 #include "gdal.h"
41 #include "xtiffio.h"
42 #include "cpl_multiproc.h"
43 #include "tifvsi.h"
44 #include "gt_wkt_srs.h"
45 #include "gt_wkt_srs_for_gdal.h"
46 #include "gt_citation.h"
47
48 CPL_CVSID("$Id: gt_wkt_srs.cpp 27182 2014-04-14 20:03:08Z rouault $")
49
50 #define ProjLinearUnitsInterpCorrectGeoKey 3059
51
52 #ifndef CT_HotineObliqueMercatorAzimuthCenter
53 # define CT_HotineObliqueMercatorAzimuthCenter 9815
54 #endif
55
56 #if !defined(GTIFAtof)
57 # define GTIFAtof CPLAtof
58 #endif
59
60 CPL_C_START
61 void CPL_DLL LibgeotiffOneTimeInit();
62 void LibgeotiffOneTimeCleanupMutex();
63
64 CPL_C_END
65
66 // To remind myself not to use CPLString in this file!
67 #define CPLString Please_do_not_use_CPLString_in_this_file
68
69 static const char *papszDatumEquiv[] =
70 {
71 "Militar_Geographische_Institut",
72 "Militar_Geographische_Institute",
73 "World_Geodetic_System_1984",
74 "WGS_1984",
75 "WGS_72_Transit_Broadcast_Ephemeris",
76 "WGS_1972_Transit_Broadcast_Ephemeris",
77 "World_Geodetic_System_1972",
78 "WGS_1972",
79 "European_Terrestrial_Reference_System_89",
80 "European_Reference_System_1989",
81 NULL
82 };
83
84 // older libgeotiff's won't list this.
85 #ifndef CT_CylindricalEqualArea
86 # define CT_CylindricalEqualArea 28
87 #endif
88
89 /************************************************************************/
90 /* LibgeotiffOneTimeInit() */
91 /************************************************************************/
92
93 static void* hMutex = NULL;
94
LibgeotiffOneTimeInit()95 void LibgeotiffOneTimeInit()
96 {
97 static int bOneTimeInitDone = FALSE;
98 CPLMutexHolder oHolder( &hMutex);
99
100 if (bOneTimeInitDone)
101 return;
102
103 bOneTimeInitDone = TRUE;
104
105 // If linking with an external libgeotiff we hope this will call the
106 // SetCSVFilenameHook() in libgeotiff, not the one in gdal/port!
107 // SetCSVFilenameHook( GDALDefaultCSVFilename );
108 }
109
110 /************************************************************************/
111 /* LibgeotiffOneTimeCleanupMutex() */
112 /************************************************************************/
113
LibgeotiffOneTimeCleanupMutex()114 void LibgeotiffOneTimeCleanupMutex()
115 {
116 if( hMutex != NULL )
117 {
118 CPLDestroyMutex(hMutex);
119 hMutex = NULL;
120 }
121 }
122
123 /************************************************************************/
124 /* GTIFToCPLRecyleString() */
125 /* */
126 /* This changes a string from the libgeotiff heap to the GDAL */
127 /* heap. */
128 /************************************************************************/
129
GTIFToCPLRecycleString(char ** ppszTarget)130 static void GTIFToCPLRecycleString( char **ppszTarget )
131
132 {
133 if( *ppszTarget == NULL )
134 return;
135
136 char *pszTempString = CPLStrdup(*ppszTarget);
137 GTIFFreeMemory( *ppszTarget );
138 *ppszTarget = pszTempString;
139 }
140
141 /************************************************************************/
142 /* WKTMassageDatum() */
143 /* */
144 /* Massage an EPSG datum name into WMT format. Also transform */
145 /* specific exception cases into WKT versions. */
146 /************************************************************************/
147
WKTMassageDatum(char ** ppszDatum)148 static void WKTMassageDatum( char ** ppszDatum )
149
150 {
151 int i, j;
152 char *pszDatum;
153
154 pszDatum = *ppszDatum;
155 if (pszDatum[0] == '\0')
156 return;
157
158 /* -------------------------------------------------------------------- */
159 /* Translate non-alphanumeric values to underscores. */
160 /* -------------------------------------------------------------------- */
161 for( i = 0; pszDatum[i] != '\0'; i++ )
162 {
163 if( pszDatum[i] != '+'
164 && !(pszDatum[i] >= 'A' && pszDatum[i] <= 'Z')
165 && !(pszDatum[i] >= 'a' && pszDatum[i] <= 'z')
166 && !(pszDatum[i] >= '0' && pszDatum[i] <= '9') )
167 {
168 pszDatum[i] = '_';
169 }
170 }
171
172 /* -------------------------------------------------------------------- */
173 /* Remove repeated and trailing underscores. */
174 /* -------------------------------------------------------------------- */
175 for( i = 1, j = 0; pszDatum[i] != '\0'; i++ )
176 {
177 if( pszDatum[j] == '_' && pszDatum[i] == '_' )
178 continue;
179
180 pszDatum[++j] = pszDatum[i];
181 }
182 if( pszDatum[j] == '_' )
183 pszDatum[j] = '\0';
184 else
185 pszDatum[j+1] = '\0';
186
187 /* -------------------------------------------------------------------- */
188 /* Search for datum equivelences. Specific massaged names get */
189 /* mapped to OpenGIS specified names. */
190 /* -------------------------------------------------------------------- */
191 for( i = 0; papszDatumEquiv[i] != NULL; i += 2 )
192 {
193 if( EQUAL(*ppszDatum,papszDatumEquiv[i]) )
194 {
195 CPLFree( *ppszDatum );
196 *ppszDatum = CPLStrdup( papszDatumEquiv[i+1] );
197 return;
198 }
199 }
200 }
201
202 /************************************************************************/
203 /* GTIFCleanupImageineNames() */
204 /* */
205 /* Erdas Imagine sometimes emits big copyright messages, and */
206 /* other stuff into citations. These can be pretty messy when */
207 /* turned into WKT, so we try to trim and clean the strings */
208 /* somewhat. */
209 /************************************************************************/
210
211 /* For example:
212 GTCitationGeoKey (Ascii,215): "IMAGINE GeoTIFF Support\nCopyright 1991 - 2001 by ERDAS, Inc. All Rights Reserved\n@(#)$RCSfile$ $Revision: 27182 $ $Date: 2014-04-14 15:03:08 -0500 (Mon, 14 Apr 2014) $\nProjection Name = UTM\nUnits = meters\nGeoTIFF Units = meters"
213
214 GeogCitationGeoKey (Ascii,267): "IMAGINE GeoTIFF Support\nCopyright 1991 - 2001 by ERDAS, Inc. All Rights Reserved\n@(#)$RCSfile$ $Revision: 27182 $ $Date: 2014-04-14 15:03:08 -0500 (Mon, 14 Apr 2014) $\nUnable to match Ellipsoid (Datum) to a GeographicTypeGeoKey value\nEllipsoid = Clarke 1866\nDatum = NAD27 (CONUS)"
215
216 PCSCitationGeoKey (Ascii,214): "IMAGINE GeoTIFF Support\nCopyright 1991 - 2001 by ERDAS, Inc. All Rights Reserved\n@(#)$RCSfile$ $Revision: 27182 $ $Date: 2014-04-14 15:03:08 -0500 (Mon, 14 Apr 2014) $\nUTM Zone 10N\nEllipsoid = Clarke 1866\nDatum = NAD27 (CONUS)"
217
218 */
219
GTIFCleanupImagineNames(char * pszCitation)220 static void GTIFCleanupImagineNames( char *pszCitation )
221
222 {
223 if( strstr(pszCitation,"IMAGINE GeoTIFF") == NULL )
224 return;
225
226 /* -------------------------------------------------------------------- */
227 /* First, we skip past all the copyright, and RCS stuff. We */
228 /* assume that this will have a "$" at the end of it all. */
229 /* -------------------------------------------------------------------- */
230 char *pszSkip;
231
232 for( pszSkip = pszCitation + strlen(pszCitation) - 1;
233 pszSkip != pszCitation && *pszSkip != '$';
234 pszSkip-- ) {}
235
236 if( *pszSkip == '$' )
237 pszSkip++;
238
239 memmove( pszCitation, pszSkip, strlen(pszSkip)+1 );
240
241 /* -------------------------------------------------------------------- */
242 /* Convert any newlines into spaces, they really gum up the */
243 /* WKT. */
244 /* -------------------------------------------------------------------- */
245 int i;
246
247 for( i = 0; pszCitation[i] != '\0'; i++ )
248 {
249 if( pszCitation[i] == '\n' )
250 pszCitation[i] = ' ';
251 }
252 }
253
254 /************************************************************************/
255 /* GTIFGetOGISDefn() */
256 /************************************************************************/
257
GTIFGetOGISDefn(GTIF * hGTIF,GTIFDefn * psDefn)258 char *GTIFGetOGISDefn( GTIF *hGTIF, GTIFDefn * psDefn )
259
260 {
261 OGRSpatialReference oSRS;
262
263 /* -------------------------------------------------------------------- */
264 /* Make sure we have hooked CSVFilename(). */
265 /* -------------------------------------------------------------------- */
266 LibgeotiffOneTimeInit();
267
268 /* -------------------------------------------------------------------- */
269 /* Handle non-standard coordinate systems where GTModelTypeGeoKey */
270 /* is not defined, but ProjectedCSTypeGeoKey is defined (ticket #3019) */
271 /* -------------------------------------------------------------------- */
272 if( psDefn->Model == KvUserDefined && psDefn->PCS != KvUserDefined)
273 {
274 psDefn->Model = ModelTypeProjected;
275 }
276
277 /* -------------------------------------------------------------------- */
278 /* Handle non-standard coordinate systems as LOCAL_CS. */
279 /* -------------------------------------------------------------------- */
280 if( psDefn->Model != ModelTypeProjected
281 && psDefn->Model != ModelTypeGeographic
282 && psDefn->Model != ModelTypeGeocentric )
283 {
284 char *pszWKT;
285 char szPeStr[2400];
286
287 /** check if there is a pe string citation key **/
288 if( GTIFKeyGet( hGTIF, PCSCitationGeoKey, szPeStr, 0, sizeof(szPeStr) ) &&
289 strstr(szPeStr, "ESRI PE String = " ) )
290 {
291 pszWKT = CPLStrdup( szPeStr + strlen("ESRI PE String = ") );
292 return pszWKT;
293 }
294 else
295 {
296 char *pszUnitsName = NULL;
297 char szPCSName[300];
298 int nKeyCount = 0;
299 int anVersion[3];
300
301 if( hGTIF != NULL )
302 GTIFDirectoryInfo( hGTIF, anVersion, &nKeyCount );
303
304 if( nKeyCount > 0 ) // Use LOCAL_CS if we have any geokeys at all.
305 {
306 // Handle citation.
307 strcpy( szPCSName, "unnamed" );
308 if( !GTIFKeyGet( hGTIF, GTCitationGeoKey, szPCSName,
309 0, sizeof(szPCSName) ) )
310 GTIFKeyGet( hGTIF, GeogCitationGeoKey, szPCSName,
311 0, sizeof(szPCSName) );
312
313 GTIFCleanupImagineNames( szPCSName );
314 oSRS.SetLocalCS( szPCSName );
315
316 // Handle units
317 GTIFGetUOMLengthInfo( psDefn->UOMLength, &pszUnitsName, NULL );
318
319 if( pszUnitsName != NULL && psDefn->UOMLength != KvUserDefined )
320 {
321 oSRS.SetLinearUnits( pszUnitsName, psDefn->UOMLengthInMeters );
322 oSRS.SetAuthority( "LOCAL_CS|UNIT", "EPSG", psDefn->UOMLength);
323 }
324 else
325 oSRS.SetLinearUnits( "unknown", psDefn->UOMLengthInMeters );
326
327 GTIFFreeMemory( pszUnitsName );
328 }
329 oSRS.exportToWkt( &pszWKT );
330
331 return pszWKT;
332 }
333 }
334
335 /* -------------------------------------------------------------------- */
336 /* Handle Geocentric coordinate systems. */
337 /* -------------------------------------------------------------------- */
338 if( psDefn->Model == ModelTypeGeocentric )
339 {
340 char szName[300];
341
342 strcpy( szName, "unnamed" );
343 if( !GTIFKeyGet( hGTIF, GTCitationGeoKey, szName,
344 0, sizeof(szName) ) )
345 GTIFKeyGet( hGTIF, GeogCitationGeoKey, szName,
346 0, sizeof(szName) );
347
348 oSRS.SetGeocCS( szName );
349
350 char *pszUnitsName = NULL;
351
352 GTIFGetUOMLengthInfo( psDefn->UOMLength, &pszUnitsName, NULL );
353
354 if( pszUnitsName != NULL && psDefn->UOMLength != KvUserDefined )
355 {
356 oSRS.SetLinearUnits( pszUnitsName, psDefn->UOMLengthInMeters );
357 oSRS.SetAuthority( "GEOCCS|UNIT", "EPSG", psDefn->UOMLength );
358 }
359 else
360 oSRS.SetLinearUnits( "unknown", psDefn->UOMLengthInMeters );
361
362 GTIFFreeMemory( pszUnitsName );
363 }
364
365 /* -------------------------------------------------------------------- */
366 /* #3901: In libgeotiff 1.3.0 and earlier we incorrectly */
367 /* interpreted linear projection parameter geokeys (false */
368 /* easting/northing) as being in meters instead of the */
369 /* coordinate system of the file. The following code attempts */
370 /* to provide mechanisms for fixing the issue if we are linked */
371 /* with an older version of libgeotiff. */
372 /* -------------------------------------------------------------------- */
373 int iParm;
374 const char *pszLinearUnits =
375 CPLGetConfigOption( "GTIFF_LINEAR_UNITS", "DEFAULT" );
376
377 #if LIBGEOTIFF_VERSION <= 1300
378 if( EQUAL(pszLinearUnits,"DEFAULT") && psDefn->Projection == KvUserDefined )
379 {
380 for( iParm = 0; iParm < psDefn->nParms; iParm++ )
381 {
382 switch( psDefn->ProjParmId[iParm] )
383 {
384 case ProjFalseEastingGeoKey:
385 case ProjFalseNorthingGeoKey:
386 case ProjFalseOriginEastingGeoKey:
387 case ProjFalseOriginNorthingGeoKey:
388 case ProjCenterEastingGeoKey:
389 case ProjCenterNorthingGeoKey:
390 if( psDefn->UOMLengthInMeters != 0
391 && psDefn->UOMLengthInMeters != 1.0 )
392 {
393 psDefn->ProjParm[iParm] *= psDefn->UOMLengthInMeters;
394 CPLDebug( "GTIFF", "converting geokey to meters to fix bug in old libgeotiff" );
395 }
396 break;
397
398 default:
399 break;
400 }
401 }
402 }
403 #endif /* LIBGEOTIFF_VERSION <= 1300 */
404
405 /* -------------------------------------------------------------------- */
406 /* #3901: If folks have broken GeoTIFF files generated with */
407 /* older versions of GDAL+libgeotiff, then they may need a */
408 /* hack to allow them to be read properly. This is that */
409 /* hack. We basically try to undue the conversion applied by */
410 /* libgeotiff to meters (or above) to simulate the old */
411 /* behavior. */
412 /* -------------------------------------------------------------------- */
413 short bLinearUnitsMarkedCorrect = FALSE;
414
415 GTIFKeyGet(hGTIF, (geokey_t) ProjLinearUnitsInterpCorrectGeoKey,
416 &bLinearUnitsMarkedCorrect, 0, 1);
417
418 if( EQUAL(pszLinearUnits,"BROKEN")
419 && psDefn->Projection == KvUserDefined
420 && !bLinearUnitsMarkedCorrect )
421 {
422 for( iParm = 0; iParm < psDefn->nParms; iParm++ )
423 {
424 switch( psDefn->ProjParmId[iParm] )
425 {
426 case ProjFalseEastingGeoKey:
427 case ProjFalseNorthingGeoKey:
428 case ProjFalseOriginEastingGeoKey:
429 case ProjFalseOriginNorthingGeoKey:
430 case ProjCenterEastingGeoKey:
431 case ProjCenterNorthingGeoKey:
432 if( psDefn->UOMLengthInMeters != 0
433 && psDefn->UOMLengthInMeters != 1.0 )
434 {
435 psDefn->ProjParm[iParm] /= psDefn->UOMLengthInMeters;
436 CPLDebug( "GTIFF", "converting geokey to accommodate old broken file due to GTIFF_LINEAR_UNITS=BROKEN setting." );
437 }
438 break;
439
440 default:
441 break;
442 }
443 }
444 }
445
446 /* -------------------------------------------------------------------- */
447 /* If this is a projected SRS we set the PROJCS keyword first */
448 /* to ensure that the GEOGCS will be a child. */
449 /* -------------------------------------------------------------------- */
450 OGRBoolean linearUnitIsSet = FALSE;
451 if( psDefn->Model == ModelTypeProjected )
452 {
453 char szCTString[512];
454 strcpy( szCTString, "unnamed" );
455 if( psDefn->PCS != KvUserDefined )
456 {
457 char *pszPCSName = NULL;
458
459 GTIFGetPCSInfo( psDefn->PCS, &pszPCSName, NULL, NULL, NULL );
460
461 oSRS.SetNode( "PROJCS", pszPCSName ? pszPCSName : "unnamed" );
462 if ( pszPCSName )
463 GTIFFreeMemory( pszPCSName );
464
465 oSRS.SetAuthority( "PROJCS", "EPSG", psDefn->PCS );
466 }
467 else if(hGTIF && GTIFKeyGet( hGTIF, PCSCitationGeoKey, szCTString, 0,
468 sizeof(szCTString)) )
469 {
470 if (!SetCitationToSRS(hGTIF, szCTString, sizeof(szCTString),
471 PCSCitationGeoKey, &oSRS, &linearUnitIsSet))
472 oSRS.SetNode("PROJCS",szCTString);
473 }
474 else
475 {
476 if( hGTIF )
477 {
478 GTIFKeyGet( hGTIF, GTCitationGeoKey, szCTString, 0, sizeof(szCTString) );
479 if(!SetCitationToSRS(hGTIF, szCTString, sizeof(szCTString),
480 GTCitationGeoKey, &oSRS, &linearUnitIsSet))
481 oSRS.SetNode( "PROJCS", szCTString );
482 }
483 else
484 oSRS.SetNode( "PROJCS", szCTString );
485 }
486
487 /* Handle ESRI/Erdas style state plane and UTM in citation key */
488 if( CheckCitationKeyForStatePlaneUTM(hGTIF, psDefn, &oSRS, &linearUnitIsSet) )
489 {
490 char *pszWKT;
491 oSRS.morphFromESRI();
492 oSRS.FixupOrdering();
493 if( oSRS.exportToWkt( &pszWKT ) == OGRERR_NONE )
494 return pszWKT;
495 }
496
497 /* Handle ESRI PE string in citation */
498 szCTString[0] = '\0';
499 if( hGTIF && GTIFKeyGet( hGTIF, GTCitationGeoKey, szCTString, 0, sizeof(szCTString) ) )
500 SetCitationToSRS(hGTIF, szCTString, sizeof(szCTString), GTCitationGeoKey, &oSRS, &linearUnitIsSet);
501 }
502
503 /* ==================================================================== */
504 /* Setup the GeogCS */
505 /* ==================================================================== */
506 char *pszGeogName = NULL;
507 char *pszDatumName = NULL;
508 char *pszPMName = NULL;
509 char *pszSpheroidName = NULL;
510 char *pszAngularUnits = NULL;
511 double dfInvFlattening=0.0, dfSemiMajor=0.0;
512 char szGCSName[512];
513 OGRBoolean aUnitGot = FALSE;
514
515 if( !GTIFGetGCSInfo( psDefn->GCS, &pszGeogName, NULL, NULL, NULL )
516 && hGTIF != NULL
517 && GTIFKeyGet( hGTIF, GeogCitationGeoKey, szGCSName, 0,
518 sizeof(szGCSName)) )
519 {
520 GetGeogCSFromCitation(szGCSName, sizeof(szGCSName),
521 GeogCitationGeoKey,
522 &pszGeogName, &pszDatumName,
523 &pszPMName, &pszSpheroidName,
524 &pszAngularUnits);
525 }
526 else
527 GTIFToCPLRecycleString( &pszGeogName );
528
529 if( !pszDatumName )
530 {
531 GTIFGetDatumInfo( psDefn->Datum, &pszDatumName, NULL );
532 GTIFToCPLRecycleString( &pszDatumName );
533 }
534 if( !pszSpheroidName )
535 {
536 GTIFGetEllipsoidInfo( psDefn->Ellipsoid, &pszSpheroidName, NULL, NULL );
537 GTIFToCPLRecycleString( &pszSpheroidName );
538 }
539 else
540 {
541 GTIFKeyGet(hGTIF, GeogSemiMajorAxisGeoKey, &(psDefn->SemiMajor), 0, 1 );
542 GTIFKeyGet(hGTIF, GeogInvFlatteningGeoKey, &dfInvFlattening, 0, 1 );
543 }
544 if( !pszPMName )
545 {
546 GTIFGetPMInfo( psDefn->PM, &pszPMName, NULL );
547 GTIFToCPLRecycleString( &pszPMName );
548 }
549 else
550 GTIFKeyGet(hGTIF, GeogPrimeMeridianLongGeoKey, &(psDefn->PMLongToGreenwich), 0, 1 );
551
552 if( !pszAngularUnits )
553 {
554 GTIFGetUOMAngleInfo( psDefn->UOMAngle, &pszAngularUnits, NULL );
555 if( pszAngularUnits == NULL )
556 pszAngularUnits = CPLStrdup("unknown");
557 else
558 GTIFToCPLRecycleString( &pszAngularUnits );
559 }
560 else
561 {
562 GTIFKeyGet(hGTIF, GeogAngularUnitSizeGeoKey, &(psDefn->UOMAngleInDegrees), 0, 1 );
563 aUnitGot = TRUE;
564 }
565
566 if( pszDatumName != NULL )
567 WKTMassageDatum( &pszDatumName );
568
569 dfSemiMajor = psDefn->SemiMajor;
570 if( dfSemiMajor == 0.0 )
571 {
572 CPLFree(pszSpheroidName);
573 pszSpheroidName = CPLStrdup("unretrievable - using WGS84");
574 dfSemiMajor = SRS_WGS84_SEMIMAJOR;
575 dfInvFlattening = SRS_WGS84_INVFLATTENING;
576 }
577 else if( dfInvFlattening == 0.0 && ((psDefn->SemiMinor / psDefn->SemiMajor) < 0.99999999999999999
578 || (psDefn->SemiMinor / psDefn->SemiMajor) > 1.00000000000000001 ) )
579 {
580 dfInvFlattening = -1.0 / (psDefn->SemiMinor/psDefn->SemiMajor - 1.0);
581
582 /* Take official inverse flattening definition in the WGS84 case */
583 if (fabs(dfSemiMajor-SRS_WGS84_SEMIMAJOR) < 1e-10 &&
584 fabs(dfInvFlattening - SRS_WGS84_INVFLATTENING) < 1e-10)
585 dfInvFlattening = SRS_WGS84_INVFLATTENING;
586 }
587 if(!pszGeogName || strlen(pszGeogName) == 0)
588 {
589 CPLFree(pszGeogName);
590 pszGeogName = CPLStrdup( pszDatumName ? pszDatumName : "unknown" );
591 }
592
593 if(aUnitGot)
594 oSRS.SetGeogCS( pszGeogName, pszDatumName,
595 pszSpheroidName, dfSemiMajor, dfInvFlattening,
596 pszPMName,
597 psDefn->PMLongToGreenwich / psDefn->UOMAngleInDegrees,
598 pszAngularUnits,
599 psDefn->UOMAngleInDegrees );
600 else
601 oSRS.SetGeogCS( pszGeogName, pszDatumName,
602 pszSpheroidName, dfSemiMajor, dfInvFlattening,
603 pszPMName,
604 psDefn->PMLongToGreenwich / psDefn->UOMAngleInDegrees,
605 pszAngularUnits,
606 psDefn->UOMAngleInDegrees * 0.0174532925199433 );
607
608 if( psDefn->GCS != KvUserDefined && psDefn->GCS > 0 )
609 oSRS.SetAuthority( "GEOGCS", "EPSG", psDefn->GCS );
610
611 if( psDefn->Datum != KvUserDefined )
612 oSRS.SetAuthority( "DATUM", "EPSG", psDefn->Datum );
613
614 if( psDefn->Ellipsoid != KvUserDefined )
615 oSRS.SetAuthority( "SPHEROID", "EPSG", psDefn->Ellipsoid );
616
617 CPLFree( pszGeogName );
618 CPLFree( pszDatumName );
619 CPLFree( pszSpheroidName );
620 CPLFree( pszPMName );
621 CPLFree( pszAngularUnits );
622
623 #if LIBGEOTIFF_VERSION >= 1310 && !defined(GEO_NORMALIZE_DISABLE_TOWGS84)
624 if( psDefn->TOWGS84Count > 0 )
625 oSRS.SetTOWGS84( psDefn->TOWGS84[0],
626 psDefn->TOWGS84[1],
627 psDefn->TOWGS84[2],
628 psDefn->TOWGS84[3],
629 psDefn->TOWGS84[4],
630 psDefn->TOWGS84[5],
631 psDefn->TOWGS84[6] );
632 #endif
633
634 /* ==================================================================== */
635 /* Handle projection parameters. */
636 /* ==================================================================== */
637 if( psDefn->Model == ModelTypeProjected )
638 {
639 /* -------------------------------------------------------------------- */
640 /* Make a local copy of parms, and convert back into the */
641 /* angular units of the GEOGCS and the linear units of the */
642 /* projection. */
643 /* -------------------------------------------------------------------- */
644 double adfParm[10];
645 int i;
646
647 for( i = 0; i < MIN(10,psDefn->nParms); i++ )
648 adfParm[i] = psDefn->ProjParm[i];
649
650 for( ; i < 10; i++ )
651 adfParm[i] = 0.0;
652
653 if(!aUnitGot)
654 {
655 adfParm[0] *= psDefn->UOMAngleInDegrees;
656 adfParm[1] *= psDefn->UOMAngleInDegrees;
657 adfParm[2] *= psDefn->UOMAngleInDegrees;
658 adfParm[3] *= psDefn->UOMAngleInDegrees;
659 }
660 int unitCode = 0;
661 GTIFKeyGet(hGTIF, ProjLinearUnitsGeoKey, &unitCode, 0, 1 );
662 if(unitCode != KvUserDefined)
663 {
664 adfParm[5] /= psDefn->UOMLengthInMeters;
665 adfParm[6] /= psDefn->UOMLengthInMeters;
666 }
667
668 /* -------------------------------------------------------------------- */
669 /* Translation the fundamental projection. */
670 /* -------------------------------------------------------------------- */
671 switch( psDefn->CTProjection )
672 {
673 case CT_TransverseMercator:
674 oSRS.SetTM( adfParm[0], adfParm[1],
675 adfParm[4],
676 adfParm[5], adfParm[6] );
677 break;
678
679 case CT_TransvMercator_SouthOriented:
680 oSRS.SetTMSO( adfParm[0], adfParm[1],
681 adfParm[4],
682 adfParm[5], adfParm[6] );
683 break;
684
685 case CT_Mercator:
686 oSRS.SetMercator( adfParm[0], adfParm[1],
687 adfParm[4],
688 adfParm[5], adfParm[6] );
689
690 if (psDefn->Projection == 1024 || psDefn->Projection == 9841) // override hack for google mercator.
691 {
692 oSRS.SetExtension( "PROJCS", "PROJ4",
693 "+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext +no_defs" );
694 }
695 break;
696
697 case CT_ObliqueStereographic:
698 oSRS.SetOS( adfParm[0], adfParm[1],
699 adfParm[4],
700 adfParm[5], adfParm[6] );
701 break;
702
703 case CT_Stereographic:
704 oSRS.SetStereographic( adfParm[0], adfParm[1],
705 adfParm[4],
706 adfParm[5], adfParm[6] );
707 break;
708
709 case CT_ObliqueMercator: /* hotine */
710 oSRS.SetHOM( adfParm[0], adfParm[1],
711 adfParm[2], adfParm[3],
712 adfParm[4],
713 adfParm[5], adfParm[6] );
714 break;
715
716 case CT_HotineObliqueMercatorAzimuthCenter:
717 #if GDAL_VERSION_MINOR >= 11
718 oSRS.SetHOMAC( adfParm[0], adfParm[1],
719 adfParm[2], adfParm[3],
720 adfParm[4],
721 adfParm[5], adfParm[6] );
722 break;
723 #endif
724 case CT_EquidistantConic:
725 oSRS.SetEC( adfParm[0], adfParm[1],
726 adfParm[2], adfParm[3],
727 adfParm[5], adfParm[6] );
728 break;
729
730 case CT_CassiniSoldner:
731 oSRS.SetCS( adfParm[0], adfParm[1],
732 adfParm[5], adfParm[6] );
733 break;
734
735 case CT_Polyconic:
736 oSRS.SetPolyconic( adfParm[0], adfParm[1],
737 adfParm[5], adfParm[6] );
738 break;
739
740 case CT_AzimuthalEquidistant:
741 oSRS.SetAE( adfParm[0], adfParm[1],
742 adfParm[5], adfParm[6] );
743 break;
744
745 case CT_MillerCylindrical:
746 oSRS.SetMC( adfParm[0], adfParm[1],
747 adfParm[5], adfParm[6] );
748 break;
749
750 case CT_Equirectangular:
751 oSRS.SetEquirectangular2( adfParm[0], adfParm[1],
752 adfParm[2],
753 adfParm[5], adfParm[6] );
754 break;
755
756 case CT_Gnomonic:
757 oSRS.SetGnomonic( adfParm[0], adfParm[1],
758 adfParm[5], adfParm[6] );
759 break;
760
761 case CT_LambertAzimEqualArea:
762 oSRS.SetLAEA( adfParm[0], adfParm[1],
763 adfParm[5], adfParm[6] );
764 break;
765
766 case CT_Orthographic:
767 oSRS.SetOrthographic( adfParm[0], adfParm[1],
768 adfParm[5], adfParm[6] );
769 break;
770
771 case CT_Robinson:
772 oSRS.SetRobinson( adfParm[1],
773 adfParm[5], adfParm[6] );
774 break;
775
776 case CT_Sinusoidal:
777 oSRS.SetSinusoidal( adfParm[1],
778 adfParm[5], adfParm[6] );
779 break;
780
781 case CT_VanDerGrinten:
782 oSRS.SetVDG( adfParm[1],
783 adfParm[5], adfParm[6] );
784 break;
785
786 case CT_PolarStereographic:
787 oSRS.SetPS( adfParm[0], adfParm[1],
788 adfParm[4],
789 adfParm[5], adfParm[6] );
790 break;
791
792 case CT_LambertConfConic_2SP:
793 oSRS.SetLCC( adfParm[2], adfParm[3],
794 adfParm[0], adfParm[1],
795 adfParm[5], adfParm[6] );
796 break;
797
798 case CT_LambertConfConic_1SP:
799 oSRS.SetLCC1SP( adfParm[0], adfParm[1],
800 adfParm[4],
801 adfParm[5], adfParm[6] );
802 break;
803
804 case CT_AlbersEqualArea:
805 oSRS.SetACEA( adfParm[0], adfParm[1],
806 adfParm[2], adfParm[3],
807 adfParm[5], adfParm[6] );
808 break;
809
810 case CT_NewZealandMapGrid:
811 oSRS.SetNZMG( adfParm[0], adfParm[1],
812 adfParm[5], adfParm[6] );
813 break;
814
815 case CT_CylindricalEqualArea:
816 oSRS.SetCEA( adfParm[0], adfParm[1],
817 adfParm[5], adfParm[6] );
818 break;
819 default:
820 if( oSRS.IsProjected() )
821 oSRS.GetRoot()->SetValue( "LOCAL_CS" );
822 break;
823 }
824
825 /* -------------------------------------------------------------------- */
826 /* Set projection units. */
827 /* -------------------------------------------------------------------- */
828 if(!linearUnitIsSet)
829 {
830 char *pszUnitsName = NULL;
831
832 GTIFGetUOMLengthInfo( psDefn->UOMLength, &pszUnitsName, NULL );
833
834 if( pszUnitsName != NULL && psDefn->UOMLength != KvUserDefined )
835 {
836 oSRS.SetLinearUnits( pszUnitsName, psDefn->UOMLengthInMeters );
837 oSRS.SetAuthority( "PROJCS|UNIT", "EPSG", psDefn->UOMLength );
838 }
839 else
840 oSRS.SetLinearUnits( "unknown", psDefn->UOMLengthInMeters );
841
842 GTIFFreeMemory( pszUnitsName );
843 }
844 }
845
846 /* ==================================================================== */
847 /* Handle vertical coordinate system information if we have it. */
848 /* ==================================================================== */
849 short verticalCSType = -1;
850 short verticalDatum = -1;
851 short verticalUnits = -1;
852 const char *pszFilename = NULL;
853 const char *pszValue;
854 char szSearchKey[128];
855 bool bNeedManualVertCS = false;
856 char citation[2048];
857
858 // Don't do anything if there is no apparent vertical information.
859 GTIFKeyGet( hGTIF, VerticalCSTypeGeoKey, &verticalCSType, 0, 1 );
860 GTIFKeyGet( hGTIF, VerticalDatumGeoKey, &verticalDatum, 0, 1 );
861 GTIFKeyGet( hGTIF, VerticalUnitsGeoKey, &verticalUnits, 0, 1 );
862
863 if( (verticalCSType != -1 || verticalDatum != -1 || verticalUnits != -1)
864 && (oSRS.IsGeographic() || oSRS.IsProjected() || oSRS.IsLocal()) )
865 {
866 if( !GTIFKeyGet( hGTIF, VerticalCitationGeoKey, &citation,
867 0, sizeof(citation) ) )
868 strcpy( citation, "unknown" );
869
870 /* -------------------------------------------------------------------- */
871 /* The original geotiff specification appears to have */
872 /* misconstrued the EPSG codes 5101 to 5106 to be vertical */
873 /* coordinate system codes, when in fact they are vertical */
874 /* datum codes. So if these are found in the */
875 /* VerticalCSTypeGeoKey move them to the VerticalDatumGeoKey */
876 /* and insert the "normal" corresponding VerticalCSTypeGeoKey */
877 /* value. */
878 /* -------------------------------------------------------------------- */
879 if( (verticalCSType >= 5101 && verticalCSType <= 5112)
880 && verticalDatum == -1 )
881 {
882 verticalDatum = verticalCSType;
883 verticalCSType = verticalDatum + 600;
884 }
885
886 /* -------------------------------------------------------------------- */
887 /* This addresses another case where the EGM96 Vertical Datum code */
888 /* is mis-used as a Vertical CS code (#4922) */
889 /* -------------------------------------------------------------------- */
890 if( verticalCSType == 5171 )
891 {
892 verticalDatum = 5171;
893 verticalCSType = 5773;
894 }
895
896 /* -------------------------------------------------------------------- */
897 /* Somewhat similarly, codes 5001 to 5033 were treated as */
898 /* vertical coordinate systems based on ellipsoidal heights. */
899 /* We use the corresponding 2d geodetic datum as the vertical */
900 /* datum and clear the vertical coordinate system code since */
901 /* there isn't one in epsg. */
902 /* -------------------------------------------------------------------- */
903 if( (verticalCSType >= 5001 && verticalCSType <= 5033)
904 && verticalDatum == -1 )
905 {
906 verticalDatum = verticalCSType+1000;
907 verticalCSType = -1;
908 }
909
910 /* -------------------------------------------------------------------- */
911 /* Promote to being a compound coordinate system. */
912 /* -------------------------------------------------------------------- */
913 OGR_SRSNode *poOldRoot = oSRS.GetRoot()->Clone();
914
915 oSRS.Clear();
916 oSRS.SetNode( "COMPD_CS", "unknown" );
917 oSRS.GetRoot()->AddChild( poOldRoot );
918
919 /* -------------------------------------------------------------------- */
920 /* If we have the vertical cs, try to look it up using the */
921 /* vertcs.csv file, and use the definition provided by that. */
922 /* -------------------------------------------------------------------- */
923 bNeedManualVertCS = true;
924
925 if( verticalCSType != KvUserDefined && verticalCSType > 0 )
926 {
927 OGRSpatialReference oVertSRS;
928 if( oVertSRS.importFromEPSG( verticalCSType ) == OGRERR_NONE )
929 {
930 oSRS.GetRoot()->AddChild( oVertSRS.GetRoot()->Clone() );
931 bNeedManualVertCS = false;
932 }
933 }
934 }
935
936 /* -------------------------------------------------------------------- */
937 /* Collect some information from the VerticalCS if not provided */
938 /* via geokeys. */
939 /* -------------------------------------------------------------------- */
940 if( bNeedManualVertCS )
941 {
942 if( verticalCSType > 0 && verticalCSType != KvUserDefined )
943 {
944 pszFilename = CSVFilename( "coordinate_reference_system.csv" );
945 sprintf( szSearchKey, "%d", verticalCSType );
946
947 if( verticalDatum < 1 || verticalDatum == KvUserDefined )
948 {
949 pszValue = CSVGetField( pszFilename,
950 "coord_ref_sys_code",
951 szSearchKey, CC_Integer,
952 "datum_code" );
953 if( pszValue != NULL )
954 verticalDatum = (short) atoi(pszValue);
955 }
956
957 if( EQUAL(citation,"unknown") )
958 {
959 pszValue = CSVGetField( pszFilename,
960 "coord_ref_sys_code",
961 szSearchKey, CC_Integer,
962 "coord_ref_sys_name" );
963 if( pszValue != NULL && *pszValue != '\0' )
964 strncpy( citation, pszValue, sizeof(citation) );
965 }
966
967 if( verticalUnits < 1 || verticalUnits == KvUserDefined )
968 {
969 pszValue = CSVGetField( pszFilename,
970 "coord_ref_sys_code",
971 szSearchKey, CC_Integer,
972 "coord_sys_code" );
973 if( pszValue != NULL )
974 {
975 pszFilename = CSVFilename( "coordinate_axis.csv" );
976 pszValue = CSVGetField( pszFilename,
977 "coord_sys_code",
978 pszValue, CC_Integer,
979 "uom_code" );
980 if( pszValue != NULL )
981 verticalUnits = (short) atoi(pszValue);
982 }
983 }
984 }
985
986 /* -------------------------------------------------------------------- */
987 /* Setup VERT_CS with citation if present. */
988 /* -------------------------------------------------------------------- */
989 oSRS.SetNode( "COMPD_CS|VERT_CS", citation );
990
991 /* -------------------------------------------------------------------- */
992 /* Setup the vertical datum. */
993 /* -------------------------------------------------------------------- */
994 const char *pszVDatumName = "unknown";
995 const char *pszVDatumType = "2005"; // CS_VD_GeoidModelDerived
996
997 if( verticalDatum > 0 && verticalDatum != KvUserDefined )
998 {
999 pszFilename = CSVFilename( "datum.csv" );
1000 if( EQUAL(pszFilename,"datum.csv") )
1001 pszFilename = CSVFilename( "gdal_datum.csv" );
1002
1003 sprintf( szSearchKey, "%d", verticalDatum );
1004
1005 pszValue = CSVGetField( pszFilename,
1006 "DATUM_CODE", szSearchKey, CC_Integer,
1007 "DATUM_NAME" );
1008 if( pszValue != NULL && *pszValue != '\0' )
1009 pszVDatumName = pszValue;
1010
1011 pszValue = CSVGetField( pszFilename,
1012 "DATUM_CODE", szSearchKey, CC_Integer,
1013 "DATUM_TYPE" );
1014 if( pszValue != NULL && EQUALN(pszValue,"geodetic",8) )
1015 pszVDatumType = "2002"; // CS_VD_Ellipsoidal
1016
1017 // We unfortunately don't know how to identify other
1018 // vertical datum types, particularly orthometric (2001).
1019 }
1020
1021 oSRS.SetNode( "COMPD_CS|VERT_CS|VERT_DATUM", pszVDatumName );
1022 oSRS.GetAttrNode( "COMPD_CS|VERT_CS|VERT_DATUM" )
1023 ->AddChild( new OGR_SRSNode( pszVDatumType ) );
1024 if( verticalDatum > 0 && verticalDatum != KvUserDefined )
1025 oSRS.SetAuthority( "COMPD_CS|VERT_CS|VERT_DATUM", "EPSG",
1026 verticalDatum );
1027
1028 /* -------------------------------------------------------------------- */
1029 /* Set the vertical units. */
1030 /* -------------------------------------------------------------------- */
1031 if( verticalUnits > 0 && verticalUnits != KvUserDefined
1032 && verticalUnits != 9001 )
1033 {
1034 char szInMeters[128];
1035
1036 pszFilename = CSVFilename("unit_of_measure.csv");
1037
1038 // Name
1039 sprintf( szSearchKey, "%d", verticalUnits );
1040 pszValue = CSVGetField( pszFilename,
1041 "uom_code", szSearchKey, CC_Integer,
1042 "unit_of_meas_name" );
1043 if( pszValue == NULL )
1044 pszValue = "unknown";
1045
1046 oSRS.SetNode( "COMPD_CS|VERT_CS|UNIT", pszValue );
1047
1048 // Value
1049 double dfFactorB, dfFactorC;
1050 dfFactorB = GTIFAtof(
1051 CSVGetField( pszFilename,
1052 "uom_code", szSearchKey, CC_Integer,
1053 "factor_b" ));
1054 dfFactorC = GTIFAtof(
1055 CSVGetField( pszFilename,
1056 "uom_code", szSearchKey, CC_Integer,
1057 "factor_c" ));
1058 if( dfFactorB != 0.0 && dfFactorC != 0.0 )
1059 sprintf( szInMeters, "%.16g", dfFactorB / dfFactorC );
1060 else
1061 strcpy( szInMeters, "1" );
1062
1063
1064 oSRS.GetAttrNode( "COMPD_CS|VERT_CS|UNIT" )
1065 ->AddChild( new OGR_SRSNode( szInMeters ) );
1066
1067 oSRS.SetAuthority( "COMPD_CS|VERT_CS|UNIT", "EPSG", verticalUnits);
1068 }
1069 else
1070 {
1071 oSRS.SetNode( "COMPD_CS|VERT_CS|UNIT", "metre" );
1072 oSRS.GetAttrNode( "COMPD_CS|VERT_CS|UNIT" )
1073 ->AddChild( new OGR_SRSNode( "1.0" ) );
1074 oSRS.SetAuthority( "COMPD_CS|VERT_CS|UNIT", "EPSG", 9001 );
1075 }
1076
1077 /* -------------------------------------------------------------------- */
1078 /* Set the axis and VERT_CS authority. */
1079 /* -------------------------------------------------------------------- */
1080 oSRS.SetNode( "COMPD_CS|VERT_CS|AXIS", "Up" );
1081 oSRS.GetAttrNode( "COMPD_CS|VERT_CS|AXIS" )
1082 ->AddChild( new OGR_SRSNode( "UP" ) );
1083
1084 if( verticalCSType > 0 && verticalCSType != KvUserDefined )
1085 oSRS.SetAuthority( "COMPD_CS|VERT_CS", "EPSG", verticalCSType );
1086 }
1087
1088 /* ==================================================================== */
1089 /* Return the WKT serialization of the object. */
1090 /* ==================================================================== */
1091 char *pszWKT;
1092
1093 oSRS.FixupOrdering();
1094
1095 if( oSRS.exportToWkt( &pszWKT ) == OGRERR_NONE )
1096 return pszWKT;
1097 else
1098 return NULL;
1099 }
1100
1101 /************************************************************************/
1102 /* OGCDatumName2EPSGDatumCode() */
1103 /************************************************************************/
1104
OGCDatumName2EPSGDatumCode(const char * pszOGCName)1105 static int OGCDatumName2EPSGDatumCode( const char * pszOGCName )
1106
1107 {
1108 FILE *fp;
1109 char **papszTokens;
1110 int nReturn = KvUserDefined;
1111
1112
1113 /* -------------------------------------------------------------------- */
1114 /* Do we know it as a built in? */
1115 /* -------------------------------------------------------------------- */
1116 if( EQUAL(pszOGCName,"NAD27")
1117 || EQUAL(pszOGCName,"North_American_Datum_1927") )
1118 return Datum_North_American_Datum_1927;
1119 else if( EQUAL(pszOGCName,"NAD83")
1120 || EQUAL(pszOGCName,"North_American_Datum_1983") )
1121 return Datum_North_American_Datum_1983;
1122 else if( EQUAL(pszOGCName,"WGS84") || EQUAL(pszOGCName,"WGS_1984")
1123 || EQUAL(pszOGCName,"WGS 84"))
1124 return Datum_WGS84;
1125 else if( EQUAL(pszOGCName,"WGS72") || EQUAL(pszOGCName,"WGS_1972") )
1126 return Datum_WGS72;
1127
1128 /* -------------------------------------------------------------------- */
1129 /* Open the table if possible. */
1130 /* -------------------------------------------------------------------- */
1131 fp = VSIFOpen( CSVFilename("gdal_datum.csv"), "r" );
1132 if( fp == NULL )
1133 fp = VSIFOpen( CSVFilename("datum.csv"), "r" );
1134
1135 if( fp == NULL )
1136 return nReturn;
1137
1138 /* -------------------------------------------------------------------- */
1139 /* Discard the first line with field names. */
1140 /* -------------------------------------------------------------------- */
1141 CSLDestroy( CSVReadParseLine( fp ) );
1142
1143 /* -------------------------------------------------------------------- */
1144 /* Read lines looking for our datum. */
1145 /* -------------------------------------------------------------------- */
1146 for( papszTokens = CSVReadParseLine( fp );
1147 CSLCount(papszTokens) > 2 && nReturn == KvUserDefined;
1148 papszTokens = CSVReadParseLine( fp ) )
1149 {
1150 WKTMassageDatum( papszTokens + 1 );
1151
1152 if( EQUAL(papszTokens[1], pszOGCName) )
1153 nReturn = atoi(papszTokens[0]);
1154
1155 CSLDestroy( papszTokens );
1156 }
1157
1158 CSLDestroy( papszTokens );
1159 VSIFClose( fp );
1160
1161 return nReturn;
1162 }
1163
1164 /************************************************************************/
1165 /* GTIFSetFromOGISDefn() */
1166 /* */
1167 /* Write GeoTIFF projection tags from an OGC WKT definition. */
1168 /************************************************************************/
1169
GTIFSetFromOGISDefn(GTIF * psGTIF,const char * pszOGCWKT)1170 int GTIFSetFromOGISDefn( GTIF * psGTIF, const char *pszOGCWKT )
1171
1172 {
1173 OGRSpatialReference *poSRS;
1174 int nPCS = KvUserDefined;
1175 OGRErr eErr;
1176 OGRBoolean peStrStored = FALSE;
1177
1178 GTIFKeySet(psGTIF, GTRasterTypeGeoKey, TYPE_SHORT, 1,
1179 RasterPixelIsArea);
1180
1181 /* -------------------------------------------------------------------- */
1182 /* Create an OGRSpatialReference object corresponding to the */
1183 /* string. */
1184 /* -------------------------------------------------------------------- */
1185 // poSRS = new OGRSpatialReference();
1186 poSRS = (OGRSpatialReference*)OSRNewSpatialReference(NULL);
1187 if( poSRS->importFromWkt((char **) &pszOGCWKT) != OGRERR_NONE )
1188 {
1189 OGRSpatialReference::DestroySpatialReference(poSRS);
1190 return FALSE;
1191 }
1192
1193 /* -------------------------------------------------------------------- */
1194 /* Get the ellipsoid definition. */
1195 /* -------------------------------------------------------------------- */
1196 short nSpheroid = KvUserDefined;
1197 double dfSemiMajor, dfInvFlattening;
1198
1199 if( poSRS->GetAuthorityName("PROJCS|GEOGCS|DATUM|SPHEROID") != NULL
1200 && EQUAL(poSRS->GetAuthorityName("PROJCS|GEOGCS|DATUM|SPHEROID"),
1201 "EPSG"))
1202 {
1203 nSpheroid = (short)
1204 atoi(poSRS->GetAuthorityCode("PROJCS|GEOGCS|DATUM|SPHEROID"));
1205 }
1206 else if( poSRS->GetAuthorityName("GEOGCS|DATUM|SPHEROID") != NULL
1207 && EQUAL(poSRS->GetAuthorityName("GEOGCS|DATUM|SPHEROID"),"EPSG"))
1208 {
1209 nSpheroid = (short)
1210 atoi(poSRS->GetAuthorityCode("GEOGCS|DATUM|SPHEROID"));
1211 }
1212
1213 dfSemiMajor = poSRS->GetSemiMajor( &eErr );
1214 dfInvFlattening = poSRS->GetInvFlattening( &eErr );
1215 if( eErr != OGRERR_NONE )
1216 {
1217 dfSemiMajor = 0.0;
1218 dfInvFlattening = 0.0;
1219 }
1220
1221 /* -------------------------------------------------------------------- */
1222 /* Get the Datum so we can special case a few PCS codes. */
1223 /* -------------------------------------------------------------------- */
1224 int nDatum = KvUserDefined;
1225
1226 if( poSRS->GetAuthorityName("PROJCS|GEOGCS|DATUM") != NULL
1227 && EQUAL(poSRS->GetAuthorityName("PROJCS|GEOGCS|DATUM"),"EPSG") )
1228 nDatum = atoi(poSRS->GetAuthorityCode("PROJCS|GEOGCS|DATUM"));
1229 else if( poSRS->GetAuthorityName("GEOGCS|DATUM") != NULL
1230 && EQUAL(poSRS->GetAuthorityName("GEOGCS|DATUM"),"EPSG") )
1231 nDatum = atoi(poSRS->GetAuthorityCode("GEOGCS|DATUM"));
1232 else if( poSRS->GetAttrValue("DATUM") != NULL )
1233 nDatum = OGCDatumName2EPSGDatumCode( poSRS->GetAttrValue("DATUM") );
1234
1235 /* -------------------------------------------------------------------- */
1236 /* Get the GCS if possible. */
1237 /* -------------------------------------------------------------------- */
1238 int nGCS = KvUserDefined;
1239
1240 if( poSRS->GetAuthorityName("PROJCS|GEOGCS") != NULL
1241 && EQUAL(poSRS->GetAuthorityName("PROJCS|GEOGCS"),"EPSG") )
1242 nGCS = atoi(poSRS->GetAuthorityCode("PROJCS|GEOGCS"));
1243 else if( poSRS->GetAuthorityName("GEOGCS") != NULL
1244 && EQUAL(poSRS->GetAuthorityName("GEOGCS"),"EPSG") )
1245 nGCS = atoi(poSRS->GetAuthorityCode("GEOGCS"));
1246
1247 if( nGCS > 32767 )
1248 nGCS = KvUserDefined;
1249
1250 /* -------------------------------------------------------------------- */
1251 /* Get the linear units. */
1252 /* -------------------------------------------------------------------- */
1253 char *pszLinearUOMName = NULL;
1254 double dfLinearUOM = poSRS->GetLinearUnits( &pszLinearUOMName );
1255 int nUOMLengthCode = 9001; /* meters */
1256
1257 if( poSRS->GetAuthorityName("PROJCS|UNIT") != NULL
1258 && EQUAL(poSRS->GetAuthorityName("PROJCS|UNIT"),"EPSG")
1259 && poSRS->GetAttrNode( "PROJCS|UNIT" ) != poSRS->GetAttrNode("GEOGCS|UNIT") )
1260 nUOMLengthCode = atoi(poSRS->GetAuthorityCode("PROJCS|UNIT"));
1261 else if( (pszLinearUOMName != NULL
1262 && EQUAL(pszLinearUOMName,SRS_UL_FOOT))
1263 || fabs(dfLinearUOM-GTIFAtof(SRS_UL_FOOT_CONV)) < 0.0000001 )
1264 nUOMLengthCode = 9002; /* international foot */
1265 else if( (pszLinearUOMName != NULL
1266 && EQUAL(pszLinearUOMName,SRS_UL_US_FOOT))
1267 || ABS(dfLinearUOM-GTIFAtof(SRS_UL_US_FOOT_CONV)) < 0.0000001 )
1268 nUOMLengthCode = 9003; /* us survey foot */
1269 else if( fabs(dfLinearUOM-1.0) > 0.00000001 )
1270 nUOMLengthCode = KvUserDefined;
1271
1272 /* -------------------------------------------------------------------- */
1273 /* Get some authority values. */
1274 /* -------------------------------------------------------------------- */
1275 if( poSRS->GetAuthorityName("PROJCS") != NULL
1276 && EQUAL(poSRS->GetAuthorityName("PROJCS"),"EPSG") )
1277 {
1278 nPCS = atoi(poSRS->GetAuthorityCode("PROJCS"));
1279 if( nPCS > 32767 )
1280 nPCS = KvUserDefined;
1281 }
1282
1283 /* -------------------------------------------------------------------- */
1284 /* Handle the projection transformation. */
1285 /* -------------------------------------------------------------------- */
1286 const char *pszProjection = poSRS->GetAttrValue( "PROJECTION" );
1287 int bWritePEString = FALSE;
1288
1289 if( nPCS != KvUserDefined )
1290 {
1291 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1292 ModelTypeProjected);
1293 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, nPCS );
1294 }
1295 else if( poSRS->IsGeocentric() )
1296 {
1297 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1298 ModelTypeGeocentric );
1299 }
1300 else if( pszProjection == NULL )
1301 {
1302 if( poSRS->IsGeographic() )
1303 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1304 ModelTypeGeographic);
1305 // otherwise, presumably something like LOCAL_CS.
1306 }
1307 else if( EQUAL(pszProjection,SRS_PT_ALBERS_CONIC_EQUAL_AREA) )
1308 {
1309 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1310 ModelTypeProjected);
1311 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1312 KvUserDefined );
1313 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1314 KvUserDefined );
1315
1316 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1317 CT_AlbersEqualArea );
1318
1319 GTIFKeySet(psGTIF, ProjStdParallelGeoKey, TYPE_DOUBLE, 1,
1320 poSRS->GetNormProjParm( SRS_PP_STANDARD_PARALLEL_1, 0.0 ) );
1321
1322 GTIFKeySet(psGTIF, ProjStdParallel2GeoKey, TYPE_DOUBLE, 1,
1323 poSRS->GetNormProjParm( SRS_PP_STANDARD_PARALLEL_2, 0.0 ) );
1324
1325 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
1326 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_CENTER, 0.0 ) );
1327
1328 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
1329 poSRS->GetNormProjParm( SRS_PP_LONGITUDE_OF_CENTER, 0.0 ) );
1330
1331 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1332 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1333
1334 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1335 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1336 }
1337
1338 else if( poSRS->GetUTMZone() != 0 )
1339 {
1340 int bNorth, nZone, nProjection;
1341
1342 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1343 ModelTypeProjected);
1344
1345 nZone = poSRS->GetUTMZone( &bNorth );
1346
1347 if( nDatum == Datum_North_American_Datum_1983 && nZone >= 3
1348 && nZone <= 22 && bNorth && nUOMLengthCode == 9001 )
1349 {
1350 nPCS = 26900 + nZone;
1351
1352 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, nPCS );
1353 }
1354 else if( nDatum == Datum_North_American_Datum_1927 && nZone >= 3
1355 && nZone <= 22 && bNorth && nUOMLengthCode == 9001 )
1356 {
1357 nPCS = 26700 + nZone;
1358
1359 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, nPCS );
1360 }
1361 else if( nDatum == Datum_WGS84 && nUOMLengthCode == 9001 )
1362 {
1363 if( bNorth )
1364 nPCS = 32600 + nZone;
1365 else
1366 nPCS = 32700 + nZone;
1367
1368 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1, nPCS );
1369 }
1370 else
1371 {
1372 if( bNorth )
1373 nProjection = 16000 + nZone;
1374 else
1375 nProjection = 16100 + nZone;
1376
1377
1378 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1379 KvUserDefined );
1380
1381 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1, nProjection );
1382 }
1383 }
1384
1385 else if( EQUAL(pszProjection,SRS_PT_TRANSVERSE_MERCATOR) )
1386 {
1387 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1388 ModelTypeProjected);
1389 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1390 KvUserDefined );
1391 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1392 KvUserDefined );
1393
1394 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1395 CT_TransverseMercator );
1396
1397 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
1398 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1399
1400 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
1401 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1402
1403 GTIFKeySet(psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1,
1404 poSRS->GetNormProjParm( SRS_PP_SCALE_FACTOR, 1.0 ) );
1405
1406 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1407 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1408
1409 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1410 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1411 }
1412
1413 else if( EQUAL(pszProjection,SRS_PT_TRANSVERSE_MERCATOR_SOUTH_ORIENTED) )
1414 {
1415 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1416 ModelTypeProjected);
1417 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1418 KvUserDefined );
1419 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1420 KvUserDefined );
1421
1422 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1423 CT_TransvMercator_SouthOriented );
1424
1425 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
1426 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1427
1428 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
1429 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1430
1431 GTIFKeySet(psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1,
1432 poSRS->GetNormProjParm( SRS_PP_SCALE_FACTOR, 1.0 ) );
1433
1434 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1435 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1436
1437 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1438 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1439 }
1440
1441 else if( EQUAL(pszProjection,SRS_PT_MERCATOR_2SP)
1442 || EQUAL(pszProjection,SRS_PT_MERCATOR_1SP) )
1443
1444 {
1445 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1446 ModelTypeProjected);
1447 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1448 KvUserDefined );
1449 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1450 KvUserDefined );
1451
1452 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1453 CT_Mercator );
1454
1455 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
1456 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1457
1458 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
1459 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1460
1461 GTIFKeySet(psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1,
1462 poSRS->GetNormProjParm( SRS_PP_SCALE_FACTOR, 1.0 ) );
1463
1464 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1465 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1466
1467 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1468 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1469 }
1470
1471 else if( EQUAL(pszProjection,SRS_PT_OBLIQUE_STEREOGRAPHIC) )
1472 {
1473 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1474 ModelTypeProjected);
1475 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1476 KvUserDefined );
1477 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1478 KvUserDefined );
1479
1480 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1481 CT_ObliqueStereographic );
1482
1483 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
1484 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1485
1486 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
1487 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1488
1489 GTIFKeySet(psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1,
1490 poSRS->GetNormProjParm( SRS_PP_SCALE_FACTOR, 1.0 ) );
1491
1492 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1493 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1494
1495 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1496 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1497 }
1498
1499 else if( EQUAL(pszProjection,SRS_PT_STEREOGRAPHIC) )
1500 {
1501 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1502 ModelTypeProjected);
1503 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1504 KvUserDefined );
1505 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1506 KvUserDefined );
1507
1508 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1509 CT_Stereographic );
1510
1511 GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1,
1512 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1513
1514 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1515 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1516
1517 GTIFKeySet(psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1,
1518 poSRS->GetNormProjParm( SRS_PP_SCALE_FACTOR, 1.0 ) );
1519
1520 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1521 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1522
1523 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1524 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1525 }
1526
1527 else if( EQUAL(pszProjection,SRS_PT_POLAR_STEREOGRAPHIC) )
1528 {
1529 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1530 ModelTypeProjected);
1531 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1532 KvUserDefined );
1533 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1534 KvUserDefined );
1535
1536 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1537 CT_PolarStereographic );
1538
1539 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
1540 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1541
1542 GTIFKeySet(psGTIF, ProjStraightVertPoleLongGeoKey, TYPE_DOUBLE, 1,
1543 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1544
1545 GTIFKeySet(psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1,
1546 poSRS->GetNormProjParm( SRS_PP_SCALE_FACTOR, 1.0 ) );
1547
1548 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1549 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1550
1551 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1552 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1553 }
1554
1555 else if( EQUAL(pszProjection,SRS_PT_HOTINE_OBLIQUE_MERCATOR) )
1556 {
1557 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1558 ModelTypeProjected);
1559 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1560 KvUserDefined );
1561 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1562 KvUserDefined );
1563
1564 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1565 CT_ObliqueMercator );
1566
1567 GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1,
1568 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_CENTER, 0.0 ) );
1569
1570 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1571 poSRS->GetNormProjParm( SRS_PP_LONGITUDE_OF_CENTER, 0.0 ) );
1572
1573 GTIFKeySet(psGTIF, ProjAzimuthAngleGeoKey, TYPE_DOUBLE, 1,
1574 poSRS->GetNormProjParm( SRS_PP_AZIMUTH, 0.0 ) );
1575
1576 GTIFKeySet(psGTIF, ProjRectifiedGridAngleGeoKey, TYPE_DOUBLE, 1,
1577 poSRS->GetNormProjParm( SRS_PP_RECTIFIED_GRID_ANGLE, 0.0 ) );
1578
1579 GTIFKeySet(psGTIF, ProjScaleAtCenterGeoKey, TYPE_DOUBLE, 1,
1580 poSRS->GetNormProjParm( SRS_PP_SCALE_FACTOR, 1.0 ) );
1581
1582 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1583 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1584
1585 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1586 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1587 }
1588
1589 #ifdef SRS_PT_HOTINE_OBLIQUE_MERCATOR_AZIMUTH_CENTER
1590 else if( EQUAL(pszProjection,SRS_PT_HOTINE_OBLIQUE_MERCATOR_AZIMUTH_CENTER) )
1591 {
1592 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1593 ModelTypeProjected);
1594 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1595 KvUserDefined );
1596 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1597 KvUserDefined );
1598
1599 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1600 CT_HotineObliqueMercatorAzimuthCenter );
1601
1602 GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1,
1603 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_CENTER, 0.0 ) );
1604
1605 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1606 poSRS->GetNormProjParm( SRS_PP_LONGITUDE_OF_CENTER, 0.0 ) );
1607
1608 GTIFKeySet(psGTIF, ProjAzimuthAngleGeoKey, TYPE_DOUBLE, 1,
1609 poSRS->GetNormProjParm( SRS_PP_AZIMUTH, 0.0 ) );
1610
1611 GTIFKeySet(psGTIF, ProjRectifiedGridAngleGeoKey, TYPE_DOUBLE, 1,
1612 poSRS->GetNormProjParm( SRS_PP_RECTIFIED_GRID_ANGLE, 0.0 ) );
1613
1614 GTIFKeySet(psGTIF, ProjScaleAtCenterGeoKey, TYPE_DOUBLE, 1,
1615 poSRS->GetNormProjParm( SRS_PP_SCALE_FACTOR, 1.0 ) );
1616
1617 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1618 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1619
1620 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1621 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1622 }
1623 #endif
1624
1625 else if( EQUAL(pszProjection,SRS_PT_CASSINI_SOLDNER) )
1626 {
1627 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1628 ModelTypeProjected);
1629 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1630 KvUserDefined );
1631 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1632 KvUserDefined );
1633
1634 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1635 CT_CassiniSoldner );
1636
1637 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
1638 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1639
1640 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
1641 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1642
1643 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1644 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1645
1646 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1647 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1648 }
1649
1650 else if( EQUAL(pszProjection,SRS_PT_EQUIDISTANT_CONIC) )
1651 {
1652 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1653 ModelTypeProjected);
1654 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1655 KvUserDefined );
1656 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1657 KvUserDefined );
1658
1659 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1660 CT_EquidistantConic );
1661
1662 GTIFKeySet(psGTIF, ProjStdParallel1GeoKey, TYPE_DOUBLE, 1,
1663 poSRS->GetNormProjParm( SRS_PP_STANDARD_PARALLEL_1, 0.0 ) );
1664
1665 GTIFKeySet(psGTIF, ProjStdParallel2GeoKey, TYPE_DOUBLE, 1,
1666 poSRS->GetNormProjParm( SRS_PP_STANDARD_PARALLEL_2, 0.0 ) );
1667
1668 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
1669 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_CENTER, 0.0 ) );
1670
1671 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
1672 poSRS->GetNormProjParm( SRS_PP_LONGITUDE_OF_CENTER, 0.0 ) );
1673
1674 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1675 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1676
1677 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1678 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1679 }
1680
1681 else if( EQUAL(pszProjection,SRS_PT_POLYCONIC) )
1682 {
1683 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1684 ModelTypeProjected);
1685 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1686 KvUserDefined );
1687 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1688 KvUserDefined );
1689
1690 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1691 CT_Polyconic );
1692
1693 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
1694 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1695
1696 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
1697 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1698
1699 GTIFKeySet(psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1,
1700 poSRS->GetNormProjParm( SRS_PP_SCALE_FACTOR, 1.0 ) );
1701
1702 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1703 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1704
1705 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1706 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1707 }
1708
1709 else if( EQUAL(pszProjection,SRS_PT_AZIMUTHAL_EQUIDISTANT) )
1710 {
1711 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1712 ModelTypeProjected);
1713 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1714 KvUserDefined );
1715 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1716 KvUserDefined );
1717
1718 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1719 CT_AzimuthalEquidistant );
1720
1721 GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1,
1722 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_CENTER, 0.0 ) );
1723
1724 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1725 poSRS->GetNormProjParm( SRS_PP_LONGITUDE_OF_CENTER, 0.0 ) );
1726
1727 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1728 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1729
1730 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1731 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1732 }
1733
1734 else if( EQUAL(pszProjection,SRS_PT_MILLER_CYLINDRICAL) )
1735 {
1736 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1737 ModelTypeProjected);
1738 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1739 KvUserDefined );
1740 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1741 KvUserDefined );
1742
1743 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1744 CT_MillerCylindrical );
1745
1746 GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1,
1747 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_CENTER, 0.0 ) );
1748
1749 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1750 poSRS->GetNormProjParm( SRS_PP_LONGITUDE_OF_CENTER, 0.0 ) );
1751
1752 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1753 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1754
1755 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1756 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1757 }
1758
1759 else if( EQUAL(pszProjection,SRS_PT_EQUIRECTANGULAR) )
1760 {
1761 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1762 ModelTypeProjected);
1763 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1764 KvUserDefined );
1765 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1766 KvUserDefined );
1767
1768 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1769 CT_Equirectangular );
1770
1771 GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1,
1772 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1773
1774 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1775 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1776
1777 GTIFKeySet(psGTIF, ProjStdParallel1GeoKey, TYPE_DOUBLE, 1,
1778 poSRS->GetNormProjParm( SRS_PP_STANDARD_PARALLEL_1, 0.0 ) );
1779
1780 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1781 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1782
1783 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1784 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1785 }
1786
1787 else if( EQUAL(pszProjection,SRS_PT_GNOMONIC) )
1788 {
1789 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1790 ModelTypeProjected);
1791 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1792 KvUserDefined );
1793 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1794 KvUserDefined );
1795
1796 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1797 CT_Gnomonic );
1798
1799 GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1,
1800 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1801
1802 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1803 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1804
1805 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1806 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1807
1808 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1809 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1810 }
1811
1812 else if( EQUAL(pszProjection,SRS_PT_LAMBERT_AZIMUTHAL_EQUAL_AREA) )
1813 {
1814 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1815 ModelTypeProjected);
1816 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1817 KvUserDefined );
1818 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1819 KvUserDefined );
1820
1821 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1822 CT_LambertAzimEqualArea );
1823
1824 GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1,
1825 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_CENTER, 0.0 ) );
1826
1827 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1828 poSRS->GetNormProjParm( SRS_PP_LONGITUDE_OF_CENTER, 0.0 ) );
1829
1830 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1831 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1832
1833 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1834 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1835 }
1836
1837 else if( EQUAL(pszProjection,SRS_PT_ORTHOGRAPHIC) )
1838 {
1839 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1840 ModelTypeProjected);
1841 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1842 KvUserDefined );
1843 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1844 KvUserDefined );
1845
1846 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1847 CT_Orthographic );
1848
1849 GTIFKeySet(psGTIF, ProjCenterLatGeoKey, TYPE_DOUBLE, 1,
1850 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1851
1852 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1853 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1854
1855 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1856 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1857
1858 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1859 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1860 }
1861
1862 else if( EQUAL(pszProjection,SRS_PT_NEW_ZEALAND_MAP_GRID) )
1863 {
1864 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1865 ModelTypeProjected);
1866 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1867 KvUserDefined );
1868 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1869 KvUserDefined );
1870
1871 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1872 CT_NewZealandMapGrid );
1873
1874 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
1875 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1876
1877 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
1878 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1879
1880 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1881 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1882
1883 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1884 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1885 }
1886
1887 else if( EQUAL(pszProjection,SRS_PT_ROBINSON) )
1888 {
1889 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1890 ModelTypeProjected);
1891 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1892 KvUserDefined );
1893 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1894 KvUserDefined );
1895
1896 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1897 CT_Robinson );
1898
1899 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1900 poSRS->GetNormProjParm( SRS_PP_LONGITUDE_OF_CENTER, 0.0 ) );
1901
1902 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1903 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1904
1905 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1906 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1907 }
1908
1909 else if( EQUAL(pszProjection,SRS_PT_SINUSOIDAL) )
1910 {
1911 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1912 ModelTypeProjected);
1913 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1914 KvUserDefined );
1915 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1916 KvUserDefined );
1917
1918 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1919 CT_Sinusoidal );
1920
1921 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1922 poSRS->GetNormProjParm( SRS_PP_LONGITUDE_OF_CENTER, 0.0 ) );
1923
1924 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1925 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1926
1927 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1928 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1929 }
1930
1931 else if( EQUAL(pszProjection,SRS_PT_VANDERGRINTEN) )
1932 {
1933 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1934 ModelTypeProjected);
1935 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1936 KvUserDefined );
1937 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1938 KvUserDefined );
1939
1940 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1941 CT_VanDerGrinten );
1942
1943 GTIFKeySet(psGTIF, ProjCenterLongGeoKey, TYPE_DOUBLE, 1,
1944 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
1945
1946 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1947 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1948
1949 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1950 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1951 }
1952
1953 else if( EQUAL(pszProjection,SRS_PT_ALBERS_CONIC_EQUAL_AREA) )
1954 {
1955 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1956 ModelTypeProjected);
1957 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1958 KvUserDefined );
1959 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1960 KvUserDefined );
1961
1962 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1963 CT_AlbersEqualArea );
1964
1965 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
1966 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_CENTER, 0.0 ) );
1967
1968 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
1969 poSRS->GetNormProjParm( SRS_PP_LONGITUDE_OF_CENTER, 0.0 ) );
1970
1971 GTIFKeySet(psGTIF, ProjStdParallel1GeoKey, TYPE_DOUBLE, 1,
1972 poSRS->GetNormProjParm( SRS_PP_STANDARD_PARALLEL_1, 0.0 ) );
1973
1974 GTIFKeySet(psGTIF, ProjStdParallel2GeoKey, TYPE_DOUBLE, 1,
1975 poSRS->GetNormProjParm( SRS_PP_STANDARD_PARALLEL_2, 0.0 ) );
1976
1977 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
1978 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
1979
1980 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
1981 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
1982 }
1983
1984 else if( EQUAL(pszProjection,SRS_PT_LAMBERT_CONFORMAL_CONIC_2SP) )
1985 {
1986 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
1987 ModelTypeProjected);
1988 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
1989 KvUserDefined );
1990 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
1991 KvUserDefined );
1992
1993 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
1994 CT_LambertConfConic_2SP );
1995
1996 GTIFKeySet(psGTIF, ProjFalseOriginLatGeoKey, TYPE_DOUBLE, 1,
1997 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
1998
1999 GTIFKeySet(psGTIF, ProjFalseOriginLongGeoKey, TYPE_DOUBLE, 1,
2000 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
2001
2002 GTIFKeySet(psGTIF, ProjStdParallel1GeoKey, TYPE_DOUBLE, 1,
2003 poSRS->GetNormProjParm( SRS_PP_STANDARD_PARALLEL_1, 0.0 ) );
2004
2005 GTIFKeySet(psGTIF, ProjStdParallel2GeoKey, TYPE_DOUBLE, 1,
2006 poSRS->GetNormProjParm( SRS_PP_STANDARD_PARALLEL_2, 0.0 ) );
2007
2008 GTIFKeySet(psGTIF, ProjFalseOriginEastingGeoKey, TYPE_DOUBLE, 1,
2009 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
2010
2011 GTIFKeySet(psGTIF, ProjFalseOriginNorthingGeoKey, TYPE_DOUBLE, 1,
2012 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
2013 }
2014
2015 else if( EQUAL(pszProjection,SRS_PT_LAMBERT_CONFORMAL_CONIC_1SP) )
2016 {
2017 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
2018 ModelTypeProjected);
2019 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
2020 KvUserDefined );
2021 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
2022 KvUserDefined );
2023
2024 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
2025 CT_LambertConfConic_1SP );
2026
2027 GTIFKeySet(psGTIF, ProjNatOriginLatGeoKey, TYPE_DOUBLE, 1,
2028 poSRS->GetNormProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0 ) );
2029
2030 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
2031 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
2032
2033 GTIFKeySet(psGTIF, ProjScaleAtNatOriginGeoKey, TYPE_DOUBLE, 1,
2034 poSRS->GetNormProjParm( SRS_PP_SCALE_FACTOR, 1.0 ) );
2035
2036 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
2037 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
2038
2039 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
2040 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
2041 }
2042
2043 else if( EQUAL(pszProjection,SRS_PT_CYLINDRICAL_EQUAL_AREA) )
2044 {
2045 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
2046 ModelTypeProjected);
2047 GTIFKeySet(psGTIF, ProjectedCSTypeGeoKey, TYPE_SHORT, 1,
2048 KvUserDefined );
2049 GTIFKeySet(psGTIF, ProjectionGeoKey, TYPE_SHORT, 1,
2050 KvUserDefined );
2051
2052 GTIFKeySet(psGTIF, ProjCoordTransGeoKey, TYPE_SHORT, 1,
2053 CT_CylindricalEqualArea );
2054
2055 GTIFKeySet(psGTIF, ProjNatOriginLongGeoKey, TYPE_DOUBLE, 1,
2056 poSRS->GetNormProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0 ) );
2057
2058 GTIFKeySet(psGTIF, ProjStdParallel1GeoKey, TYPE_DOUBLE, 1,
2059 poSRS->GetNormProjParm( SRS_PP_STANDARD_PARALLEL_1, 0.0 ) );
2060
2061 GTIFKeySet(psGTIF, ProjFalseEastingGeoKey, TYPE_DOUBLE, 1,
2062 poSRS->GetProjParm( SRS_PP_FALSE_EASTING, 0.0 ) );
2063
2064 GTIFKeySet(psGTIF, ProjFalseNorthingGeoKey, TYPE_DOUBLE, 1,
2065 poSRS->GetProjParm( SRS_PP_FALSE_NORTHING, 0.0 ) );
2066 }
2067
2068 else
2069 {
2070 bWritePEString = TRUE;
2071 }
2072
2073 // Note that VERTCS is an ESRI "spelling" of VERT_CS so we assume if
2074 // we find it that we should try to treat this as a PE string.
2075 bWritePEString |= (poSRS->GetAttrValue("VERTCS") != NULL);
2076
2077 if( bWritePEString
2078 && CSLTestBoolean( CPLGetConfigOption("GTIFF_ESRI_CITATION",
2079 "YES") ) )
2080 {
2081 /* Anyhing we can't map, we store as an ESRI PE string with a citation key */
2082 char *pszPEString = NULL;
2083 poSRS->morphToESRI();
2084 poSRS->exportToWkt( &pszPEString );
2085 int peStrLen = strlen(pszPEString);
2086 if(peStrLen > 0)
2087 {
2088 char *outPeStr = (char *) CPLMalloc( peStrLen + strlen("ESRI PE String = ")+1 );
2089 strcpy(outPeStr, "ESRI PE String = ");
2090 strcat(outPeStr, pszPEString);
2091 GTIFKeySet( psGTIF, PCSCitationGeoKey, TYPE_ASCII, 0, outPeStr );
2092 peStrStored = TRUE;
2093 CPLFree( outPeStr );
2094 }
2095 if(pszPEString)
2096 CPLFree( pszPEString );
2097 GTIFKeySet(psGTIF, GTModelTypeGeoKey, TYPE_SHORT, 1,
2098 KvUserDefined );
2099 }
2100
2101 /* -------------------------------------------------------------------- */
2102 /* Is there a false easting/northing set? If so, write out a */
2103 /* special geokey tag to indicate that GDAL has written these */
2104 /* with the proper interpretation of the linear units. */
2105 /* -------------------------------------------------------------------- */
2106 double dfFE = 0.0, dfFN = 0.0;
2107
2108 if( (GTIFKeyGet(psGTIF, ProjFalseEastingGeoKey, &dfFE, 0, 1)
2109 || GTIFKeyGet(psGTIF, ProjFalseNorthingGeoKey, &dfFN, 0, 1)
2110 || GTIFKeyGet(psGTIF, ProjFalseOriginEastingGeoKey, &dfFE, 0, 1)
2111 || GTIFKeyGet(psGTIF, ProjFalseOriginNorthingGeoKey, &dfFN, 0, 1))
2112 && (dfFE != 0.0 || dfFN != 0.0)
2113 && nUOMLengthCode != 9001 )
2114 {
2115 GTIFKeySet(psGTIF, (geokey_t) ProjLinearUnitsInterpCorrectGeoKey,
2116 TYPE_SHORT, 1, (short) 1 );
2117 }
2118
2119 /* -------------------------------------------------------------------- */
2120 /* Write linear units information. */
2121 /* -------------------------------------------------------------------- */
2122 if( poSRS->IsGeocentric() )
2123 {
2124 GTIFKeySet(psGTIF, GeogLinearUnitsGeoKey, TYPE_SHORT, 1,
2125 nUOMLengthCode );
2126 if( nUOMLengthCode == KvUserDefined )
2127 GTIFKeySet( psGTIF, GeogLinearUnitSizeGeoKey, TYPE_DOUBLE, 1,
2128 dfLinearUOM);
2129 }
2130 else if( !poSRS->IsGeographic() )
2131 {
2132 GTIFKeySet(psGTIF, ProjLinearUnitsGeoKey, TYPE_SHORT, 1,
2133 nUOMLengthCode );
2134 if( nUOMLengthCode == KvUserDefined )
2135 GTIFKeySet( psGTIF, ProjLinearUnitSizeGeoKey, TYPE_DOUBLE, 1,
2136 dfLinearUOM);
2137
2138 /* if linear units name is available and user defined, store it as citation */
2139 if(!peStrStored
2140 && nUOMLengthCode == KvUserDefined
2141 && pszLinearUOMName
2142 && strlen(pszLinearUOMName)>0
2143 && CSLTestBoolean( CPLGetConfigOption("GTIFF_ESRI_CITATION",
2144 "YES") ) )
2145 {
2146 SetLinearUnitCitation(psGTIF, pszLinearUOMName);
2147 }
2148 }
2149
2150 /* -------------------------------------------------------------------- */
2151 /* Write angular units. Always Degrees for now. */
2152 /* Changed to support different angular units */
2153 /* -------------------------------------------------------------------- */
2154
2155 char* angUnitName = NULL;
2156 double angUnitValue = poSRS->GetAngularUnits(&angUnitName);
2157 if(EQUAL(angUnitName, "Degree"))
2158 GTIFKeySet(psGTIF, GeogAngularUnitsGeoKey, TYPE_SHORT, 1,
2159 Angular_Degree );
2160 else if(angUnitName)
2161 {
2162 GTIFKeySet(psGTIF, GeogCitationGeoKey, TYPE_ASCII, 0,
2163 angUnitName ); // it may be rewritten if the gcs is userdefined
2164 GTIFKeySet(psGTIF, GeogAngularUnitSizeGeoKey, TYPE_DOUBLE, 1,
2165 angUnitValue );
2166 }
2167
2168 /* -------------------------------------------------------------------- */
2169 /* Try to write a citation from the main coordinate system */
2170 /* name. */
2171 /* -------------------------------------------------------------------- */
2172 if( poSRS->GetRoot() != NULL
2173 && poSRS->GetRoot()->GetChild(0) != NULL
2174 && (poSRS->IsProjected() || poSRS->IsLocal() || poSRS->IsGeocentric()) )
2175 {
2176 GTIFKeySet( psGTIF, GTCitationGeoKey, TYPE_ASCII, 0,
2177 poSRS->GetRoot()->GetChild(0)->GetValue() );
2178 }
2179
2180 /* -------------------------------------------------------------------- */
2181 /* Try to write a GCS citation. */
2182 /* -------------------------------------------------------------------- */
2183 OGR_SRSNode *poGCS = poSRS->GetAttrNode( "GEOGCS" );
2184
2185 if( poGCS != NULL && poGCS->GetChild(0) != NULL )
2186 {
2187 GTIFKeySet( psGTIF, GeogCitationGeoKey, TYPE_ASCII, 0,
2188 poGCS->GetChild(0)->GetValue() );
2189 }
2190
2191 /* -------------------------------------------------------------------- */
2192 /* Try to identify the GCS/datum, scanning the EPSG datum file for */
2193 /* a match. */
2194 /* -------------------------------------------------------------------- */
2195 if( nPCS == KvUserDefined )
2196 {
2197 if( nGCS == KvUserDefined )
2198 {
2199 if( nDatum == Datum_North_American_Datum_1927 )
2200 nGCS = GCS_NAD27;
2201 else if( nDatum == Datum_North_American_Datum_1983 )
2202 nGCS = GCS_NAD83;
2203 else if( nDatum == Datum_WGS84 || nDatum == DatumE_WGS84 )
2204 nGCS = GCS_WGS_84;
2205 }
2206
2207 if( nGCS != KvUserDefined )
2208 {
2209 GTIFKeySet( psGTIF, GeographicTypeGeoKey, TYPE_SHORT,
2210 1, nGCS );
2211 }
2212 else if( nDatum != KvUserDefined )
2213 {
2214 GTIFKeySet( psGTIF, GeographicTypeGeoKey, TYPE_SHORT, 1,
2215 KvUserDefined );
2216 GTIFKeySet( psGTIF, GeogGeodeticDatumGeoKey, TYPE_SHORT,
2217 1, nDatum );
2218 }
2219 else if( nSpheroid != KvUserDefined )
2220 {
2221 GTIFKeySet( psGTIF, GeographicTypeGeoKey, TYPE_SHORT, 1,
2222 KvUserDefined );
2223 GTIFKeySet( psGTIF, GeogGeodeticDatumGeoKey, TYPE_SHORT,
2224 1, KvUserDefined );
2225 GTIFKeySet( psGTIF, GeogEllipsoidGeoKey, TYPE_SHORT, 1,
2226 nSpheroid );
2227 }
2228 else if( dfSemiMajor != 0.0 )
2229 {
2230 GTIFKeySet( psGTIF, GeographicTypeGeoKey, TYPE_SHORT, 1,
2231 KvUserDefined );
2232 GTIFKeySet( psGTIF, GeogGeodeticDatumGeoKey, TYPE_SHORT,
2233 1, KvUserDefined );
2234 GTIFKeySet( psGTIF, GeogEllipsoidGeoKey, TYPE_SHORT, 1,
2235 KvUserDefined );
2236 GTIFKeySet( psGTIF, GeogSemiMajorAxisGeoKey, TYPE_DOUBLE, 1,
2237 dfSemiMajor );
2238 if( dfInvFlattening == 0.0 )
2239 GTIFKeySet( psGTIF, GeogSemiMinorAxisGeoKey, TYPE_DOUBLE, 1,
2240 dfSemiMajor );
2241 else
2242 GTIFKeySet( psGTIF, GeogInvFlatteningGeoKey, TYPE_DOUBLE, 1,
2243 dfInvFlattening );
2244 }
2245 else if( poSRS->GetAttrValue("DATUM") != NULL
2246 && strstr(poSRS->GetAttrValue("DATUM"),"unknown") == NULL
2247 && strstr(poSRS->GetAttrValue("DATUM"),"unnamed") == NULL )
2248
2249 {
2250 CPLError( CE_Warning, CPLE_AppDefined,
2251 "Couldn't translate `%s' to a GeoTIFF datum.\n",
2252 poSRS->GetAttrValue("DATUM") );
2253 }
2254
2255 /* Always set InvFlattening if it is avaliable. */
2256 /* So that it doesn'tneed to calculate from SemiMinor */
2257 if( dfInvFlattening != 0.0 )
2258 GTIFKeySet( psGTIF, GeogInvFlatteningGeoKey, TYPE_DOUBLE, 1,
2259 dfInvFlattening );
2260 /* Always set SemiMajor to keep the precision and in case of editing */
2261 if( dfSemiMajor != 0.0 )
2262 GTIFKeySet( psGTIF, GeogSemiMajorAxisGeoKey, TYPE_DOUBLE, 1,
2263 dfSemiMajor );
2264
2265 if( nGCS == KvUserDefined
2266 && CSLTestBoolean( CPLGetConfigOption("GTIFF_ESRI_CITATION",
2267 "YES") ) )
2268 SetGeogCSCitation(psGTIF, poSRS, angUnitName, nDatum, nSpheroid);
2269 }
2270
2271 /* -------------------------------------------------------------------- */
2272 /* Do we have TOWGS84 parameters? */
2273 /* -------------------------------------------------------------------- */
2274
2275 #if LIBGEOTIFF_VERSION >= 1310 && !defined(GEO_NORMALIZE_DISABLE_TOWGS84)
2276 double adfTOWGS84[7];
2277
2278 if( poSRS->GetTOWGS84( adfTOWGS84 ) == OGRERR_NONE )
2279 {
2280 if( adfTOWGS84[3] == 0.0 && adfTOWGS84[4] == 0.0
2281 && adfTOWGS84[5] == 0.0 && adfTOWGS84[6] == 0.0 )
2282 {
2283 if( nGCS == GCS_WGS_84 && adfTOWGS84[0] == 0.0
2284 && adfTOWGS84[1] == 0.0 && adfTOWGS84[2] == 0.0 )
2285 {
2286 ; /* do nothing */
2287 }
2288 else
2289 GTIFKeySet( psGTIF, GeogTOWGS84GeoKey, TYPE_DOUBLE, 3,
2290 adfTOWGS84 );
2291 }
2292 else
2293 GTIFKeySet( psGTIF, GeogTOWGS84GeoKey, TYPE_DOUBLE, 7,
2294 adfTOWGS84 );
2295 }
2296 #endif
2297
2298 /* -------------------------------------------------------------------- */
2299 /* Do we have vertical datum information to set? */
2300 /* -------------------------------------------------------------------- */
2301 if( poSRS->GetAttrValue( "COMPD_CS|VERT_CS" ) != NULL )
2302 {
2303 const char *pszValue;
2304
2305 GTIFKeySet( psGTIF, VerticalCitationGeoKey, TYPE_ASCII, 0,
2306 poSRS->GetAttrValue( "COMPD_CS|VERT_CS" ) );
2307
2308 pszValue = poSRS->GetAuthorityCode( "COMPD_CS|VERT_CS" );
2309 if( pszValue && atoi(pszValue) )
2310 GTIFKeySet( psGTIF, VerticalCSTypeGeoKey, TYPE_SHORT, 1,
2311 atoi(pszValue) );
2312
2313 pszValue = poSRS->GetAuthorityCode( "COMPD_CS|VERT_CS|VERT_DATUM" );
2314 if( pszValue && atoi(pszValue) )
2315 GTIFKeySet( psGTIF, VerticalDatumGeoKey, TYPE_SHORT, 1,
2316 atoi(pszValue) );
2317
2318 pszValue = poSRS->GetAuthorityCode( "COMPD_CS|VERT_CS|UNIT" );
2319 if( pszValue && atoi(pszValue) )
2320 GTIFKeySet( psGTIF, VerticalUnitsGeoKey, TYPE_SHORT, 1,
2321 atoi(pszValue) );
2322 }
2323
2324 /* -------------------------------------------------------------------- */
2325 /* Cleanup */
2326 /* -------------------------------------------------------------------- */
2327 OGRSpatialReference::DestroySpatialReference(poSRS);
2328 return TRUE;
2329 }
2330
2331 /************************************************************************/
2332 /* GTIFWktFromMemBuf() */
2333 /************************************************************************/
2334
GTIFWktFromMemBuf(int nSize,unsigned char * pabyBuffer,char ** ppszWKT,double * padfGeoTransform,int * pnGCPCount,GDAL_GCP ** ppasGCPList)2335 CPLErr GTIFWktFromMemBuf( int nSize, unsigned char *pabyBuffer,
2336 char **ppszWKT, double *padfGeoTransform,
2337 int *pnGCPCount, GDAL_GCP **ppasGCPList )
2338 {
2339 return GTIFWktFromMemBufEx(nSize, pabyBuffer, ppszWKT, padfGeoTransform,
2340 pnGCPCount, ppasGCPList, NULL);
2341 }
2342
GTIFWktFromMemBufEx(int nSize,unsigned char * pabyBuffer,char ** ppszWKT,double * padfGeoTransform,int * pnGCPCount,GDAL_GCP ** ppasGCPList,int * pbPixelIsPoint)2343 CPLErr GTIFWktFromMemBufEx( int nSize, unsigned char *pabyBuffer,
2344 char **ppszWKT, double *padfGeoTransform,
2345 int *pnGCPCount, GDAL_GCP **ppasGCPList,
2346 int *pbPixelIsPoint )
2347
2348 {
2349 bool bPixelIsPoint = false;
2350 int bPointGeoIgnore = FALSE;
2351 short nRasterType;
2352 char szFilename[100];
2353
2354 sprintf( szFilename, "/vsimem/wkt_from_mem_buf_%ld.tif",
2355 (long) CPLGetPID() );
2356
2357 /* -------------------------------------------------------------------- */
2358 /* Make sure we have hooked CSVFilename(). */
2359 /* -------------------------------------------------------------------- */
2360 LibgeotiffOneTimeInit();
2361
2362 /* -------------------------------------------------------------------- */
2363 /* Create a memory file from the buffer. */
2364 /* -------------------------------------------------------------------- */
2365 VSILFILE *fp = VSIFileFromMemBuffer( szFilename, pabyBuffer, nSize, FALSE );
2366 if( fp == NULL )
2367 return CE_Failure;
2368 VSIFCloseL( fp );
2369
2370 /* -------------------------------------------------------------------- */
2371 /* Initialize access to the memory geotiff structure. */
2372 /* -------------------------------------------------------------------- */
2373 TIFF *hTIFF;
2374 hTIFF = VSI_TIFFOpen( szFilename, "rc" );
2375
2376 if( hTIFF == NULL )
2377 {
2378 CPLError( CE_Failure, CPLE_AppDefined,
2379 "TIFF/GeoTIFF structure is corrupt." );
2380 VSIUnlink( szFilename );
2381 return CE_Failure;
2382 }
2383
2384 /* -------------------------------------------------------------------- */
2385 /* Get the projection definition. */
2386 /* -------------------------------------------------------------------- */
2387 GTIF *hGTIF;
2388 GTIFDefn *psGTIFDefn;
2389
2390 hGTIF = GTIFNew(hTIFF);
2391
2392 if( hGTIF != NULL && GTIFKeyGet(hGTIF, GTRasterTypeGeoKey, &nRasterType,
2393 0, 1 ) == 1
2394 && nRasterType == (short) RasterPixelIsPoint )
2395 {
2396 bPixelIsPoint = true;
2397 bPointGeoIgnore =
2398 CSLTestBoolean( CPLGetConfigOption("GTIFF_POINT_GEO_IGNORE",
2399 "FALSE") );
2400 }
2401 if( pbPixelIsPoint )
2402 *pbPixelIsPoint = bPixelIsPoint;
2403
2404 #if LIBGEOTIFF_VERSION >= 1410
2405 psGTIFDefn = GTIFAllocDefn();
2406 #else
2407 psGTIFDefn = (GTIFDefn *) CPLCalloc(1,sizeof(GTIFDefn));
2408 #endif
2409
2410
2411 if( hGTIF != NULL && GTIFGetDefn( hGTIF, psGTIFDefn ) )
2412 *ppszWKT = GTIFGetOGISDefn( hGTIF, psGTIFDefn );
2413 else
2414 *ppszWKT = NULL;
2415
2416 if( hGTIF )
2417 GTIFFree( hGTIF );
2418
2419 #if LIBGEOTIFF_VERSION >= 1410
2420 GTIFFreeDefn(psGTIFDefn);
2421 #else
2422 CPLFree(psGTIFDefn);
2423 #endif
2424
2425 /* -------------------------------------------------------------------- */
2426 /* Get geotransform or tiepoints. */
2427 /* -------------------------------------------------------------------- */
2428 double *padfTiePoints, *padfScale, *padfMatrix;
2429 int16 nCount;
2430
2431 padfGeoTransform[0] = 0.0;
2432 padfGeoTransform[1] = 1.0;
2433 padfGeoTransform[2] = 0.0;
2434 padfGeoTransform[3] = 0.0;
2435 padfGeoTransform[4] = 0.0;
2436 padfGeoTransform[5] = 1.0;
2437
2438 *pnGCPCount = 0;
2439 *ppasGCPList = NULL;
2440
2441 if( TIFFGetField(hTIFF,TIFFTAG_GEOPIXELSCALE,&nCount,&padfScale )
2442 && nCount >= 2 )
2443 {
2444 padfGeoTransform[1] = padfScale[0];
2445 padfGeoTransform[5] = - ABS(padfScale[1]);
2446
2447 if( TIFFGetField(hTIFF,TIFFTAG_GEOTIEPOINTS,&nCount,&padfTiePoints )
2448 && nCount >= 6 )
2449 {
2450 padfGeoTransform[0] =
2451 padfTiePoints[3] - padfTiePoints[0] * padfGeoTransform[1];
2452 padfGeoTransform[3] =
2453 padfTiePoints[4] - padfTiePoints[1] * padfGeoTransform[5];
2454
2455 // adjust for pixel is point in transform
2456 if( bPixelIsPoint && !bPointGeoIgnore )
2457 {
2458 padfGeoTransform[0] -= (padfGeoTransform[1] * 0.5 + padfGeoTransform[2] * 0.5);
2459 padfGeoTransform[3] -= (padfGeoTransform[4] * 0.5 + padfGeoTransform[5] * 0.5);
2460 }
2461 }
2462 }
2463
2464 else if( TIFFGetField(hTIFF,TIFFTAG_GEOTIEPOINTS,&nCount,&padfTiePoints )
2465 && nCount >= 6 )
2466 {
2467 *pnGCPCount = nCount / 6;
2468 *ppasGCPList = (GDAL_GCP *) CPLCalloc(sizeof(GDAL_GCP),*pnGCPCount);
2469
2470 for( int iGCP = 0; iGCP < *pnGCPCount; iGCP++ )
2471 {
2472 char szID[32];
2473 GDAL_GCP *psGCP = *ppasGCPList + iGCP;
2474
2475 sprintf( szID, "%d", iGCP+1 );
2476 psGCP->pszId = CPLStrdup( szID );
2477 psGCP->pszInfo = CPLStrdup("");
2478 psGCP->dfGCPPixel = padfTiePoints[iGCP*6+0];
2479 psGCP->dfGCPLine = padfTiePoints[iGCP*6+1];
2480 psGCP->dfGCPX = padfTiePoints[iGCP*6+3];
2481 psGCP->dfGCPY = padfTiePoints[iGCP*6+4];
2482 psGCP->dfGCPZ = padfTiePoints[iGCP*6+5];
2483 }
2484 }
2485
2486 else if( TIFFGetField(hTIFF,TIFFTAG_GEOTRANSMATRIX,&nCount,&padfMatrix )
2487 && nCount == 16 )
2488 {
2489 padfGeoTransform[0] = padfMatrix[3];
2490 padfGeoTransform[1] = padfMatrix[0];
2491 padfGeoTransform[2] = padfMatrix[1];
2492 padfGeoTransform[3] = padfMatrix[7];
2493 padfGeoTransform[4] = padfMatrix[4];
2494 padfGeoTransform[5] = padfMatrix[5];
2495 }
2496
2497 /* -------------------------------------------------------------------- */
2498 /* Cleanup. */
2499 /* -------------------------------------------------------------------- */
2500 XTIFFClose( hTIFF );
2501
2502 VSIUnlink( szFilename );
2503
2504 if( *ppszWKT == NULL )
2505 return CE_Failure;
2506 else
2507 return CE_None;
2508 }
2509
2510 /************************************************************************/
2511 /* GTIFMemBufFromWkt() */
2512 /************************************************************************/
2513
GTIFMemBufFromWkt(const char * pszWKT,const double * padfGeoTransform,int nGCPCount,const GDAL_GCP * pasGCPList,int * pnSize,unsigned char ** ppabyBuffer)2514 CPLErr GTIFMemBufFromWkt( const char *pszWKT, const double *padfGeoTransform,
2515 int nGCPCount, const GDAL_GCP *pasGCPList,
2516 int *pnSize, unsigned char **ppabyBuffer )
2517 {
2518 return GTIFMemBufFromWktEx(pszWKT, padfGeoTransform,
2519 nGCPCount,pasGCPList,
2520 pnSize, ppabyBuffer, FALSE);
2521 }
2522
GTIFMemBufFromWktEx(const char * pszWKT,const double * padfGeoTransform,int nGCPCount,const GDAL_GCP * pasGCPList,int * pnSize,unsigned char ** ppabyBuffer,int bPixelIsPoint)2523 CPLErr GTIFMemBufFromWktEx( const char *pszWKT, const double *padfGeoTransform,
2524 int nGCPCount, const GDAL_GCP *pasGCPList,
2525 int *pnSize, unsigned char **ppabyBuffer,
2526 int bPixelIsPoint )
2527
2528 {
2529 TIFF *hTIFF;
2530 GTIF *hGTIF;
2531 char szFilename[100];
2532
2533 sprintf( szFilename, "/vsimem/wkt_from_mem_buf_%ld.tif",
2534 (long) CPLGetPID() );
2535
2536 /* -------------------------------------------------------------------- */
2537 /* Make sure we have hooked CSVFilename(). */
2538 /* -------------------------------------------------------------------- */
2539 LibgeotiffOneTimeInit();
2540
2541 /* -------------------------------------------------------------------- */
2542 /* Initialize access to the memory geotiff structure. */
2543 /* -------------------------------------------------------------------- */
2544 hTIFF = VSI_TIFFOpen( szFilename, "w" );
2545
2546 if( hTIFF == NULL )
2547 {
2548 CPLError( CE_Failure, CPLE_AppDefined,
2549 "TIFF/GeoTIFF structure is corrupt." );
2550 return CE_Failure;
2551 }
2552
2553 /* -------------------------------------------------------------------- */
2554 /* Write some minimal set of image parameters. */
2555 /* -------------------------------------------------------------------- */
2556 TIFFSetField( hTIFF, TIFFTAG_IMAGEWIDTH, 1 );
2557 TIFFSetField( hTIFF, TIFFTAG_IMAGELENGTH, 1 );
2558 TIFFSetField( hTIFF, TIFFTAG_BITSPERSAMPLE, 8 );
2559 TIFFSetField( hTIFF, TIFFTAG_SAMPLESPERPIXEL, 1 );
2560 TIFFSetField( hTIFF, TIFFTAG_ROWSPERSTRIP, 1 );
2561 TIFFSetField( hTIFF, TIFFTAG_PLANARCONFIG, PLANARCONFIG_CONTIG );
2562 TIFFSetField( hTIFF, TIFFTAG_PHOTOMETRIC, PHOTOMETRIC_MINISBLACK );
2563
2564 /* -------------------------------------------------------------------- */
2565 /* Get the projection definition. */
2566 /* -------------------------------------------------------------------- */
2567
2568 int bPointGeoIgnore = FALSE;
2569 if( bPixelIsPoint )
2570 {
2571 bPointGeoIgnore =
2572 CSLTestBoolean( CPLGetConfigOption("GTIFF_POINT_GEO_IGNORE",
2573 "FALSE") );
2574 }
2575
2576 if( pszWKT != NULL || bPixelIsPoint )
2577 {
2578 hGTIF = GTIFNew(hTIFF);
2579 if( pszWKT != NULL )
2580 GTIFSetFromOGISDefn( hGTIF, pszWKT );
2581
2582 if( bPixelIsPoint )
2583 {
2584 GTIFKeySet(hGTIF, GTRasterTypeGeoKey, TYPE_SHORT, 1,
2585 RasterPixelIsPoint);
2586 }
2587
2588 GTIFWriteKeys( hGTIF );
2589 GTIFFree( hGTIF );
2590 }
2591
2592 /* -------------------------------------------------------------------- */
2593 /* Set the geotransform, or GCPs. */
2594 /* -------------------------------------------------------------------- */
2595
2596 if( padfGeoTransform[0] != 0.0 || padfGeoTransform[1] != 1.0
2597 || padfGeoTransform[2] != 0.0 || padfGeoTransform[3] != 0.0
2598 || padfGeoTransform[4] != 0.0 || ABS(padfGeoTransform[5]) != 1.0 )
2599 {
2600
2601 if( padfGeoTransform[2] == 0.0 && padfGeoTransform[4] == 0.0 )
2602 {
2603 double adfPixelScale[3], adfTiePoints[6];
2604
2605 adfPixelScale[0] = padfGeoTransform[1];
2606 adfPixelScale[1] = fabs(padfGeoTransform[5]);
2607 adfPixelScale[2] = 0.0;
2608
2609 TIFFSetField( hTIFF, TIFFTAG_GEOPIXELSCALE, 3, adfPixelScale );
2610
2611 adfTiePoints[0] = 0.0;
2612 adfTiePoints[1] = 0.0;
2613 adfTiePoints[2] = 0.0;
2614 adfTiePoints[3] = padfGeoTransform[0];
2615 adfTiePoints[4] = padfGeoTransform[3];
2616 adfTiePoints[5] = 0.0;
2617
2618 if( bPixelIsPoint && !bPointGeoIgnore )
2619 {
2620 adfTiePoints[3] += padfGeoTransform[1] * 0.5 + padfGeoTransform[2] * 0.5;
2621 adfTiePoints[4] += padfGeoTransform[4] * 0.5 + padfGeoTransform[5] * 0.5;
2622 }
2623
2624 TIFFSetField( hTIFF, TIFFTAG_GEOTIEPOINTS, 6, adfTiePoints );
2625 }
2626 else
2627 {
2628 double adfMatrix[16];
2629
2630 memset(adfMatrix,0,sizeof(double) * 16);
2631
2632 adfMatrix[0] = padfGeoTransform[1];
2633 adfMatrix[1] = padfGeoTransform[2];
2634 adfMatrix[3] = padfGeoTransform[0];
2635 adfMatrix[4] = padfGeoTransform[4];
2636 adfMatrix[5] = padfGeoTransform[5];
2637 adfMatrix[7] = padfGeoTransform[3];
2638 adfMatrix[15] = 1.0;
2639
2640 if( bPixelIsPoint && !bPointGeoIgnore )
2641 {
2642 adfMatrix[3] += padfGeoTransform[1] * 0.5 + padfGeoTransform[2] * 0.5;
2643 adfMatrix[7] += padfGeoTransform[4] * 0.5 + padfGeoTransform[5] * 0.5;
2644 }
2645
2646 TIFFSetField( hTIFF, TIFFTAG_GEOTRANSMATRIX, 16, adfMatrix );
2647 }
2648 }
2649
2650 /* -------------------------------------------------------------------- */
2651 /* Otherwise write tiepoints if they are available. */
2652 /* -------------------------------------------------------------------- */
2653 else if( nGCPCount > 0 )
2654 {
2655 double *padfTiePoints;
2656
2657 padfTiePoints = (double *) CPLMalloc(6*sizeof(double)*nGCPCount);
2658
2659 for( int iGCP = 0; iGCP < nGCPCount; iGCP++ )
2660 {
2661
2662 padfTiePoints[iGCP*6+0] = pasGCPList[iGCP].dfGCPPixel;
2663 padfTiePoints[iGCP*6+1] = pasGCPList[iGCP].dfGCPLine;
2664 padfTiePoints[iGCP*6+2] = 0;
2665 padfTiePoints[iGCP*6+3] = pasGCPList[iGCP].dfGCPX;
2666 padfTiePoints[iGCP*6+4] = pasGCPList[iGCP].dfGCPY;
2667 padfTiePoints[iGCP*6+5] = pasGCPList[iGCP].dfGCPZ;
2668 }
2669
2670 TIFFSetField( hTIFF, TIFFTAG_GEOTIEPOINTS, 6*nGCPCount, padfTiePoints);
2671 CPLFree( padfTiePoints );
2672 }
2673
2674 /* -------------------------------------------------------------------- */
2675 /* Cleanup and return the created memory buffer. */
2676 /* -------------------------------------------------------------------- */
2677 GByte bySmallImage = 0;
2678
2679 TIFFWriteEncodedStrip( hTIFF, 0, (char *) &bySmallImage, 1 );
2680 TIFFWriteCheck( hTIFF, TIFFIsTiled(hTIFF), "GTIFMemBufFromWkt");
2681 TIFFWriteDirectory( hTIFF );
2682
2683 XTIFFClose( hTIFF );
2684
2685 /* -------------------------------------------------------------------- */
2686 /* Read back from the memory buffer. It would be preferrable */
2687 /* to be able to "steal" the memory buffer, but there isn't */
2688 /* currently any support for this. */
2689 /* -------------------------------------------------------------------- */
2690 GUIntBig nBigLength;
2691
2692 *ppabyBuffer = VSIGetMemFileBuffer( szFilename, &nBigLength, TRUE );
2693 *pnSize = (int) nBigLength;
2694
2695 return CE_None;
2696 }
2697