1 : /*****************************************************************************
2 : * $Id: IdrisiDataset.cpp 19780 2010-05-28 14:53:15Z ilucena $
3 : *
4 : * Project: Idrisi Raster Image File Driver
5 : * Purpose: Read/write Idrisi Raster Image Format RST
6 : * Author: Ivan Lucena, ivan.lucena@pmldnet.com
7 : *
8 : ******************************************************************************
9 : * Copyright( c ) 2006, Ivan Lucena
10 : *
11 : * Permission is hereby granted, free of charge, to any person obtaining a
12 : * copy of this software and associated documentation files( the "Software" ),
13 : * to deal in the Software without restriction, including without limitation
14 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
15 : * and/or sell copies of the Software, and to permit persons to whom the
16 : * Software is furnished to do so, subject to the following conditions:
17 : *
18 : * The above copyright notice and this permission notice shall be included
19 : * in all copies or substantial portions of the Software.
20 : *
21 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
22 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
23 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
24 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
25 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
26 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
27 : * DEALINGS IN THE SOFTWARE.
28 : ****************************************************************************/
29 :
30 : #include "gdal_priv.h"
31 : #include "cpl_conv.h"
32 : #include "cpl_string.h"
33 : #include "cpl_csv.h"
34 : #include "ogr_spatialref.h"
35 : #include "gdal_pam.h"
36 : #include "gdal_alg.h"
37 : #include "gdal_rat.h"
38 :
39 : CPL_CVSID( "$Id: IdrisiDataset.cpp 19780 2010-05-28 14:53:15Z ilucena $" );
40 :
41 : CPL_C_START
42 : void GDALRegister_IDRISI( void);
43 : CPL_C_END
44 :
45 : #ifdef WIN32
46 : #define PATHDELIM '\\'
47 : #else
48 : #define PATHDELIM '/'
49 : #endif
50 :
51 : //----- Safe numeric conversion, NULL as zero
52 : #define atoi_nz(s) (s == NULL ? (int) 0 : atoi(s))
53 : #define atof_nz(s) (s == NULL ? (double) 0.0 : atof(s))
54 :
55 : //----- file extensions:
56 : #define extRST "rst"
57 : #define extRDC "rdc"
58 : #define extSMP "smp"
59 : #define extREF "ref"
60 : #define extRSTu "RST"
61 : #define extRDCu "RDC"
62 : #define extSMPu "SMP"
63 : #define extREFu "REF"
64 :
65 : //----- field names on rdc file:
66 : #define rdcFILE_FORMAT "file format "
67 : #define rdcFILE_TITLE "file title "
68 : #define rdcDATA_TYPE "data type "
69 : #define rdcFILE_TYPE "file type "
70 : #define rdcCOLUMNS "columns "
71 : #define rdcROWS "rows "
72 : #define rdcREF_SYSTEM "ref. system "
73 : #define rdcREF_UNITS "ref. units "
74 : #define rdcUNIT_DIST "unit dist. "
75 : #define rdcMIN_X "min. X "
76 : #define rdcMAX_X "max. X "
77 : #define rdcMIN_Y "min. Y "
78 : #define rdcMAX_Y "max. Y "
79 : #define rdcPOSN_ERROR "pos'n error "
80 : #define rdcRESOLUTION "resolution "
81 : #define rdcMIN_VALUE "min. value "
82 : #define rdcMAX_VALUE "max. value "
83 : #define rdcDISPLAY_MIN "display min "
84 : #define rdcDISPLAY_MAX "display max "
85 : #define rdcVALUE_UNITS "value units "
86 : #define rdcVALUE_ERROR "value error "
87 : #define rdcFLAG_VALUE "flag value "
88 : #define rdcFLAG_DEFN "flag def'n "
89 : #define rdcFLAG_DEFN2 "flag def`n "
90 : #define rdcLEGEND_CATS "legend cats "
91 : #define rdcLINEAGES "lineage "
92 : #define rdcCOMMENTS "comment "
93 : #define rdcCODE_N "code %6d "
94 :
95 : //----- ".ref" file field names:
96 : #define refREF_SYSTEM "ref. system "
97 : #define refREF_SYSTEM2 "ref.system "
98 : #define refPROJECTION "projection "
99 : #define refDATUM "datum "
100 : #define refDELTA_WGS84 "delta WGS84 "
101 : #define refELLIPSOID "ellipsoid "
102 : #define refMAJOR_SAX "major s-ax "
103 : #define refMINOR_SAX "minor s-ax "
104 : #define refORIGIN_LONG "origin long "
105 : #define refORIGIN_LAT "origin lat "
106 : #define refORIGIN_X "origin X "
107 : #define refORIGIN_Y "origin Y "
108 : #define refSCALE_FAC "scale fac "
109 : #define refUNITS "units "
110 : #define refPARAMETERS "parameters "
111 : #define refSTANDL_1 "stand ln 1 "
112 : #define refSTANDL_2 "stand ln 2 "
113 :
114 : //----- standard values:
115 : #define rstVERSION "Idrisi Raster A.1"
116 : #define rstBYTE "byte"
117 : #define rstINTEGER "integer"
118 : #define rstREAL "real"
119 : #define rstRGB24 "rgb24"
120 : #define rstDEGREE "deg"
121 : #define rstMETER "m"
122 : #define rstLATLONG "latlong"
123 : #define rstLATLONG2 "lat/long"
124 : #define rstPLANE "plane"
125 : #define rstUTM "utm-%d%c"
126 : #define rstSPC "spc%2d%2s%d"
127 :
128 : //----- palette file( .smp ) header size:
129 : #define smpHEADERSIZE 18
130 :
131 : //----- check if file exists:
132 : bool FileExists( const char *pszPath );
133 :
134 : //----- Reference Table
135 : struct ReferenceTab {
136 : int nCode;
137 : const char *pszName;
138 : };
139 :
140 : //----- USA State's reference table to USGS PCS Code
141 : static const ReferenceTab aoUSStateTable[] = {
142 : {101, "al"},
143 : {201, "az"},
144 : {301, "ar"},
145 : {401, "ca"},
146 : {501, "co"},
147 : {600, "ct"},
148 : {700, "de"},
149 : {901, "fl"},
150 : {1001, "ga"},
151 : {1101, "id"},
152 : {1201, "il"},
153 : {1301, "in"},
154 : {1401, "ia"},
155 : {1501, "ks"},
156 : {1601, "ky"},
157 : {1701, "la"},
158 : {1801, "me"},
159 : {1900, "md"},
160 : {2001, "ma"},
161 : {2111, "mi"},
162 : {2201, "mn"},
163 : {2301, "ms"},
164 : {2401, "mo"},
165 : {2500, "mt"},
166 : {2600, "ne"},
167 : {2701, "nv"},
168 : {2800, "nh"},
169 : {2900, "nj"},
170 : {3001, "nm"},
171 : {3101, "ny"},
172 : {3200, "nc"},
173 : {3301, "nd"},
174 : {3401, "oh"},
175 : {3501, "ok"},
176 : {3601, "or"},
177 : {3701, "pa"},
178 : {3800, "ri"},
179 : {3900, "sc"},
180 : {4001, "sd"},
181 : {4100, "tn"},
182 : {4201, "tx"},
183 : {4301, "ut"},
184 : {4400, "vt"},
185 : {4501, "va"},
186 : {4601, "wa"},
187 : {4701, "wv"},
188 : {4801, "wv"},
189 : {4901, "wy"},
190 : {5001, "ak"},
191 : {5101, "hi"},
192 : {5200, "pr"}
193 : };
194 : #define US_STATE_COUNT ( sizeof( aoUSStateTable ) / sizeof( ReferenceTab ) )
195 :
196 : //----- Get the Code of a US State
197 : int GetStateCode ( const char *pszState );
198 :
199 : //----- Get the state name of a Code
200 : const char *GetStateName( int nCode );
201 :
202 : //----- Conversion Table definition
203 : struct ConvertionTab {
204 : const char *pszName;
205 : int nDefaultI;
206 : int nDefaultG;
207 : double dfConv;
208 : };
209 :
210 : //----- Linear Unit Conversion Table
211 : static const ConvertionTab aoLinearUnitsConv[] = {
212 : {"m", /* 0 */ 0, 1, 1.0},
213 : {SRS_UL_METER, /* 1 */ 0, 1, 1.0},
214 : {"meters", /* 2 */ 0, 1, 1.0},
215 : {"metre", /* 3 */ 0, 1, 1.0},
216 :
217 : {"ft", /* 4 */ 4, 5, 0.3048},
218 : {SRS_UL_FOOT, /* 5 */ 4, 5, 0.3048},
219 : {"feet", /* 6 */ 4, 5, 0.3048},
220 : {"foot_us", /* 7 */ 4, 5, 0.3048006},
221 : {"u.s. foot", /* 8 */ 4, 5, 0.3048006},
222 :
223 : {"mi", /* 9 */ 9, 10, 1612.9},
224 : {"mile", /* 10 */ 9, 10, 1612.9},
225 : {"miles", /* 11 */ 9, 10, 1612.9},
226 :
227 : {"km", /* 12 */ 12, 13, 1000.0},
228 : {"kilometers", /* 13 */ 12, 13, 1000.0},
229 : {"kilometer", /* 14 */ 12, 13, 1000.0},
230 : {"kilometre", /* 15 */ 12, 13, 1000.0},
231 :
232 : {"deg", /* 16 */ 16, 17, 0.0},
233 : {SRS_UA_DEGREE, /* 17 */ 16, 17, 0.0},
234 : {"degrees", /* 18 */ 16, 17, 0.0},
235 :
236 : {"rad", /* 19 */ 19, 20, 0.0},
237 : {SRS_UA_RADIAN, /* 20 */ 19, 20, 0.0},
238 : {"radians", /* 21 */ 19, 20, 0.0}
239 : };
240 : #define LINEAR_UNITS_COUNT (sizeof(aoLinearUnitsConv) / sizeof(ConvertionTab))
241 :
242 : //----- Get the index of a given linear unit
243 : int GetUnitIndex( const char *pszUnitName );
244 :
245 : //----- Get the defaut name
246 : char *GetUnitDefault( const char *pszUnitName, const char *pszToMeter = NULL );
247 :
248 : //----- Get the "to meter"
249 : int GetToMeterIndex( const char *pszToMeter );
250 :
251 : //----- CSLSaveCRLF
252 : int SaveAsCRLF(char **papszStrList, const char *pszFname);
253 :
254 : //----- Classes pre-definition:
255 : class IdrisiDataset;
256 : class IdrisiRasterBand;
257 :
258 : // ----------------------------------------------------------------------------
259 : // Idrisi GDALDataset
260 : // ----------------------------------------------------------------------------
261 :
262 : class IdrisiDataset : public GDALPamDataset
263 : {
264 : friend class IdrisiRasterBand;
265 :
266 : private:
267 : FILE *fp;
268 :
269 : char *pszFilename;
270 : char *pszDocFilename;
271 : char **papszRDC;
272 : double adfGeoTransform[6];
273 :
274 : char *pszProjection;
275 : char **papszCategories;
276 : char *pszUnitType;
277 :
278 : CPLErr GeoReference2Wkt( const char *pszRefSystem,
279 : const char *pszRefUnits,
280 : char **pszProjString );
281 :
282 : CPLErr Wkt2GeoReference( const char *pszProjString,
283 : char **pszRefSystem,
284 : char **pszRefUnit );
285 :
286 : protected:
287 : GDALColorTable *poColorTable;
288 :
289 : public:
290 : IdrisiDataset();
291 : ~IdrisiDataset();
292 :
293 : static GDALDataset *Open( GDALOpenInfo *poOpenInfo );
294 : static GDALDataset *Create( const char *pszFilename,
295 : int nXSize,
296 : int nYSize,
297 : int nBands,
298 : GDALDataType eType,
299 : char **papszOptions );
300 : static GDALDataset *CreateCopy( const char *pszFilename,
301 : GDALDataset *poSrcDS,
302 : int bStrict,
303 : char **papszOptions,
304 : GDALProgressFunc pfnProgress,
305 : void * pProgressData );
306 : virtual char **GetFileList(void);
307 : virtual CPLErr GetGeoTransform( double *padfTransform );
308 : virtual CPLErr SetGeoTransform( double *padfTransform );
309 : virtual const char *GetProjectionRef( void );
310 : virtual CPLErr SetProjection( const char *pszProjString );
311 : };
312 :
313 : // ----------------------------------------------------------------------------
314 : // Idrisi GDALPamRasterBand
315 : // ----------------------------------------------------------------------------
316 :
317 : class IdrisiRasterBand : public GDALPamRasterBand
318 : {
319 : friend class IdrisiDataset;
320 :
321 : GDALRasterAttributeTable *poDefaultRAT;
322 :
323 : private:
324 : int nRecordSize;
325 : GByte *pabyScanLine;
326 :
327 : public:
328 : IdrisiRasterBand( IdrisiDataset *poDS,
329 : int nBand,
330 : GDALDataType eDataType );
331 : ~IdrisiRasterBand();
332 :
333 : virtual double GetNoDataValue( int *pbSuccess = NULL );
334 : virtual double GetMinimum( int *pbSuccess = NULL );
335 : virtual double GetMaximum( int *pbSuccess = NULL );
336 : virtual CPLErr IReadBlock( int nBlockXOff, int nBlockYOff, void *pImage );
337 : virtual CPLErr IWriteBlock( int nBlockXOff, int nBlockYOff, void *pImage );
338 : virtual GDALColorTable *GetColorTable();
339 : virtual GDALColorInterp GetColorInterpretation();
340 : virtual char **GetCategoryNames();
341 : virtual const char *GetUnitType();
342 :
343 : virtual CPLErr SetCategoryNames( char **papszCategoryNames );
344 : virtual CPLErr SetNoDataValue( double dfNoDataValue );
345 : virtual CPLErr SetColorTable( GDALColorTable *poColorTable );
346 : virtual CPLErr SetUnitType( const char *pszUnitType );
347 : virtual CPLErr SetStatistics( double dfMin, double dfMax,
348 : double dfMean, double dfStdDev );
349 : virtual const GDALRasterAttributeTable *GetDefaultRAT();
350 : virtual CPLErr SetDefaultRAT( const GDALRasterAttributeTable * );
351 :
352 : float dfMaximum;
353 : float dfMinimum;
354 : bool bFirstVal;
355 : };
356 :
357 : // ------------------------------------------------------------------------ //
358 : // Implementation of IdrisiDataset //
359 : // ------------------------------------------------------------------------ //
360 :
361 : /************************************************************************/
362 : /* IdrisiDataset() */
363 : /************************************************************************/
364 :
365 16 : IdrisiDataset::IdrisiDataset()
366 : {
367 16 : pszFilename = NULL;
368 16 : fp = NULL;
369 16 : papszRDC = NULL;
370 16 : pszDocFilename = NULL;
371 16 : pszProjection = NULL;
372 16 : poColorTable = new GDALColorTable();
373 16 : papszCategories = NULL;
374 16 : pszUnitType = NULL;
375 :
376 16 : adfGeoTransform[0] = 0.0;
377 16 : adfGeoTransform[1] = 1.0;
378 16 : adfGeoTransform[2] = 0.0;
379 16 : adfGeoTransform[3] = 0.0;
380 16 : adfGeoTransform[4] = 0.0;
381 16 : adfGeoTransform[5] = 1.0;
382 16 : }
383 :
384 : /************************************************************************/
385 : /* ~IdrisiDataset() */
386 : /************************************************************************/
387 :
388 16 : IdrisiDataset::~IdrisiDataset()
389 : {
390 16 : FlushCache();
391 :
392 16 : if( papszRDC != NULL )
393 : {
394 16 : if( eAccess == GA_Update )
395 : {
396 : int i;
397 :
398 16 : for( i = 0; i < nBands; i++ )
399 : {
400 10 : IdrisiRasterBand *poBand = (IdrisiRasterBand*) GetRasterBand( i + 1 );
401 10 : poBand->SetStatistics( poBand->dfMinimum, poBand->dfMaximum, 0.0, 0.0 );
402 : }
403 :
404 6 : CSLSetNameValueSeparator( papszRDC, ": " );
405 6 : SaveAsCRLF( papszRDC, pszDocFilename );
406 : }
407 16 : CSLDestroy( papszRDC );
408 : }
409 :
410 16 : if( poColorTable )
411 : {
412 16 : delete poColorTable;
413 : }
414 16 : CPLFree( pszFilename );
415 16 : CPLFree( pszDocFilename );
416 16 : CPLFree( pszProjection );
417 16 : CSLDestroy( papszCategories );
418 16 : CPLFree( pszUnitType );
419 :
420 16 : if( fp != NULL )
421 16 : VSIFCloseL( fp );
422 16 : }
423 :
424 : /************************************************************************/
425 : /* Open() */
426 : /************************************************************************/
427 :
428 10165 : GDALDataset *IdrisiDataset::Open( GDALOpenInfo *poOpenInfo )
429 : {
430 10165 : if( ( poOpenInfo->fp == NULL ) ||
431 : ( EQUAL( CPLGetExtension( poOpenInfo->pszFilename ), extRST ) == FALSE ) )
432 10149 : return NULL;
433 :
434 : // --------------------------------------------------------------------
435 : // Check the documentation file .rdc
436 : // --------------------------------------------------------------------
437 :
438 16 : const char *pszLDocFilename = CPLResetExtension( poOpenInfo->pszFilename, extRDC );
439 :
440 16 : if( ! FileExists( pszLDocFilename ) )
441 : {
442 0 : pszLDocFilename = CPLResetExtension( poOpenInfo->pszFilename, extRDCu );
443 :
444 0 : if( ! FileExists( pszLDocFilename ) )
445 : {
446 0 : return NULL;
447 : }
448 : }
449 :
450 16 : char **papszLRDC = CSLLoad( pszLDocFilename );
451 :
452 16 : CSLSetNameValueSeparator( papszLRDC, ":" );
453 :
454 16 : const char *pszVersion = CSLFetchNameValue( papszLRDC, rdcFILE_FORMAT );
455 :
456 16 : if( pszVersion == NULL || !EQUAL( pszVersion, rstVERSION ) )
457 : {
458 0 : CSLDestroy( papszLRDC );
459 0 : return NULL;
460 : }
461 :
462 : // --------------------------------------------------------------------
463 : // Create a corresponding GDALDataset
464 : // --------------------------------------------------------------------
465 :
466 : IdrisiDataset *poDS;
467 :
468 16 : poDS = new IdrisiDataset();
469 16 : poDS->eAccess = poOpenInfo->eAccess;
470 16 : poDS->pszFilename = CPLStrdup( poOpenInfo->pszFilename );
471 :
472 16 : if( poOpenInfo->eAccess == GA_ReadOnly )
473 : {
474 10 : poDS->fp = VSIFOpenL( poDS->pszFilename, "rb" );
475 : }
476 : else
477 : {
478 6 : poDS->fp = VSIFOpenL( poDS->pszFilename, "r+b" );
479 : }
480 :
481 16 : if( poDS->fp == NULL )
482 : {
483 0 : CSLDestroy( papszLRDC );
484 0 : return NULL;
485 : }
486 :
487 16 : poDS->pszDocFilename = CPLStrdup( pszLDocFilename );
488 16 : poDS->papszRDC = CSLDuplicate( papszLRDC );
489 16 : CSLDestroy( papszLRDC );
490 :
491 : // --------------------------------------------------------------------
492 : // Load information from rdc
493 : // --------------------------------------------------------------------
494 :
495 16 : poDS->nRasterXSize = atoi_nz( CSLFetchNameValue( poDS->papszRDC, rdcCOLUMNS ) );
496 16 : poDS->nRasterYSize = atoi_nz( CSLFetchNameValue( poDS->papszRDC, rdcROWS ) );
497 :
498 : // --------------------------------------------------------------------
499 : // Create band information
500 : // --------------------------------------------------------------------
501 :
502 16 : const char *pszDataType = CSLFetchNameValue( poDS->papszRDC, rdcDATA_TYPE );
503 :
504 16 : if( EQUAL( pszDataType, rstBYTE ) )
505 : {
506 3 : poDS->nBands = 1;
507 3 : poDS->SetBand( 1, new IdrisiRasterBand( poDS, 1, GDT_Byte ) );
508 : }
509 13 : else if( EQUAL( pszDataType, rstINTEGER ) )
510 : {
511 2 : poDS->nBands = 1;
512 2 : poDS->SetBand( 1, new IdrisiRasterBand( poDS, 1, GDT_Int16 ) );
513 : }
514 11 : else if( EQUAL( pszDataType, rstREAL ) )
515 : {
516 6 : poDS->nBands = 1;
517 6 : poDS->SetBand( 1, new IdrisiRasterBand( poDS, 1, GDT_Float32 ) );
518 : }
519 5 : else if( EQUAL( pszDataType, rstRGB24 ) )
520 : {
521 5 : poDS->nBands = 3;
522 5 : poDS->SetBand( 1, new IdrisiRasterBand( poDS, 1, GDT_Byte ) );
523 10 : poDS->SetBand( 2, new IdrisiRasterBand( poDS, 2, GDT_Byte ) );
524 10 : poDS->SetBand( 3, new IdrisiRasterBand( poDS, 3, GDT_Byte ) );
525 : }
526 : else
527 : {
528 : CPLError(CE_Failure, CPLE_AppDefined,
529 0 : "Unknown data type : %s", pszDataType);
530 0 : delete poDS;
531 0 : return NULL;
532 : }
533 :
534 : int i;
535 42 : for(i=0;i<poDS->nBands;i++)
536 : {
537 26 : IdrisiRasterBand* band = (IdrisiRasterBand*) poDS->GetRasterBand(i+1);
538 26 : if (band->pabyScanLine == NULL)
539 : {
540 0 : delete poDS;
541 0 : return NULL;
542 : }
543 : }
544 :
545 : // --------------------------------------------------------------------
546 : // Load the transformation matrix
547 : // --------------------------------------------------------------------
548 :
549 16 : const char *pszMinX = CSLFetchNameValue( poDS->papszRDC, rdcMIN_X );
550 :
551 16 : if( strlen( pszMinX ) > 0 )
552 : {
553 : double dfMinX, dfMaxX, dfMinY, dfMaxY, dfUnit, dfXPixSz, dfYPixSz;
554 :
555 16 : dfMinX = atof_nz( CSLFetchNameValue( poDS->papszRDC, rdcMIN_X ) );
556 16 : dfMaxX = atof_nz( CSLFetchNameValue( poDS->papszRDC, rdcMAX_X ) );
557 16 : dfMinY = atof_nz( CSLFetchNameValue( poDS->papszRDC, rdcMIN_Y ) );
558 16 : dfMaxY = atof_nz( CSLFetchNameValue( poDS->papszRDC, rdcMAX_Y ) );
559 16 : dfUnit = atof_nz( CSLFetchNameValue( poDS->papszRDC, rdcUNIT_DIST ) );
560 :
561 16 : dfMinX = dfMinX * dfUnit;
562 16 : dfMaxX = dfMaxX * dfUnit;
563 16 : dfMinY = dfMinY * dfUnit;
564 16 : dfMaxY = dfMaxY * dfUnit;
565 :
566 16 : dfYPixSz = ( dfMinY - dfMaxY ) / poDS->nRasterYSize;
567 16 : dfXPixSz = ( dfMaxX - dfMinX ) / poDS->nRasterXSize;
568 :
569 16 : poDS->adfGeoTransform[0] = dfMinX;
570 16 : poDS->adfGeoTransform[1] = dfXPixSz;
571 16 : poDS->adfGeoTransform[2] = 0.0;
572 16 : poDS->adfGeoTransform[3] = dfMaxY;
573 16 : poDS->adfGeoTransform[4] = 0.0;
574 16 : poDS->adfGeoTransform[5] = dfYPixSz;
575 : }
576 :
577 : // --------------------------------------------------------------------
578 : // Set Color Table in the presence of a smp file
579 : // --------------------------------------------------------------------
580 :
581 16 : if( poDS->nBands != 3 )
582 : {
583 11 : const char *pszSMPFilename = CPLResetExtension( poDS->pszFilename, extSMP );
584 : FILE *fpSMP;
585 11 : if( ( fpSMP = VSIFOpenL( pszSMPFilename, "rb" ) ) != NULL )
586 : {
587 0 : int dfMaxValue = atoi_nz( CSLFetchNameValue( poDS->papszRDC, rdcMAX_VALUE ) );
588 0 : int nCatCount = atoi_nz( CSLFetchNameValue( poDS->papszRDC, rdcLEGEND_CATS ) );
589 0 : if( nCatCount == 0 )
590 0 : dfMaxValue = 255;
591 0 : VSIFSeekL( fpSMP, smpHEADERSIZE, SEEK_SET );
592 : GDALColorEntry oEntry;
593 : unsigned char aucRGB[3];
594 0 : int i = 0;
595 0 : while( ( VSIFReadL( &aucRGB, sizeof( aucRGB ), 1, fpSMP ) ) &&( i <= dfMaxValue ) )
596 : {
597 0 : oEntry.c1 = (short) aucRGB[0];
598 0 : oEntry.c2 = (short) aucRGB[1];
599 0 : oEntry.c3 = (short) aucRGB[2];
600 0 : oEntry.c4 = (short) 255;
601 0 : poDS->poColorTable->SetColorEntry( i, &oEntry );
602 0 : i++;
603 : }
604 0 : VSIFCloseL( fpSMP );
605 : }
606 : }
607 :
608 : // --------------------------------------------------------------------
609 : // Check for Unit Type
610 : // --------------------------------------------------------------------
611 :
612 16 : const char *pszValueUnit = CSLFetchNameValue( poDS->papszRDC, rdcVALUE_UNITS );
613 :
614 16 : if( pszValueUnit == NULL )
615 0 : poDS->pszUnitType = CPLStrdup( "unspecified" );
616 : else
617 : {
618 16 : if( EQUALN( pszValueUnit, "meter", 4 ) )
619 : {
620 0 : poDS->pszUnitType = CPLStrdup( "m" );
621 : }
622 16 : else if( EQUALN( pszValueUnit, "feet", 4 ) )
623 : {
624 0 : poDS->pszUnitType = CPLStrdup( "ft" );
625 : }
626 : else
627 16 : poDS->pszUnitType = CPLStrdup( pszValueUnit );
628 : }
629 :
630 : // --------------------------------------------------------------------
631 : // Check for category names.
632 : // --------------------------------------------------------------------
633 :
634 16 : int nCatCount = atoi_nz( CSLFetchNameValue( poDS->papszRDC, rdcLEGEND_CATS ) );
635 :
636 16 : if( nCatCount > 0 )
637 : {
638 : // ----------------------------------------------------------------
639 : // Sequentialize categories names, from 0 to the last "code n"
640 : // ----------------------------------------------------------------
641 :
642 0 : int nCode = 0;
643 0 : int nCount = 0;
644 0 : int nLine = -1;
645 0 : for( int i = 0;( i < CSLCount( poDS->papszRDC ) ) &&( nLine == -1 ); i++ )
646 0 : if( EQUALN( poDS->papszRDC[i], rdcLEGEND_CATS, 11 ) )
647 0 : nLine = i;
648 0 : if( nLine > 0 )
649 : {
650 0 : sscanf( poDS->papszRDC[++nLine], rdcCODE_N, &nCode );
651 0 : for( int i = 0;( i < 255 ) &&( nCount < nCatCount ); i++ )
652 : {
653 0 : if( i == nCode )
654 : {
655 : poDS->papszCategories =
656 : CSLAddString( poDS->papszCategories,
657 0 : CPLParseNameValue( poDS->papszRDC[nLine], NULL ) );
658 0 : nCount++;
659 0 : if( nCount < nCatCount )
660 0 : sscanf( poDS->papszRDC[++nLine], rdcCODE_N, &nCode );
661 : }
662 : else
663 0 : poDS->papszCategories = CSLAddString( poDS->papszCategories, "" );
664 : }
665 : }
666 : }
667 :
668 : /* -------------------------------------------------------------------- */
669 : /* Automatic Generated Color Table */
670 : /* -------------------------------------------------------------------- */
671 :
672 16 : if( poDS->papszCategories != NULL &&
673 : ( poDS->poColorTable->GetColorEntryCount() == 0 ) )
674 : {
675 0 : int nEntryCount = CSLCount(poDS->papszCategories);
676 :
677 : GDALColorEntry sFromColor;
678 0 : sFromColor.c1 = (short) ( 255 );
679 0 : sFromColor.c2 = (short) ( 0 );
680 0 : sFromColor.c3 = (short) ( 0 );
681 0 : sFromColor.c4 = (short) ( 255 );
682 :
683 : GDALColorEntry sToColor;
684 0 : sToColor.c1 = (short) ( 0 );
685 0 : sToColor.c2 = (short) ( 0 );
686 0 : sToColor.c3 = (short) ( 255 );
687 0 : sToColor.c4 = (short) ( 255 );
688 :
689 : poDS->poColorTable->CreateColorRamp(
690 0 : 0, &sFromColor, ( nEntryCount - 1 ), &sToColor );
691 : }
692 :
693 : /* -------------------------------------------------------------------- */
694 : /* Initialize any PAM information. */
695 : /* -------------------------------------------------------------------- */
696 :
697 16 : poDS->SetDescription( poOpenInfo->pszFilename );
698 16 : poDS->TryLoadXML();
699 :
700 : /* -------------------------------------------------------------------- */
701 : /* Check for external overviews. */
702 : /* -------------------------------------------------------------------- */
703 :
704 16 : poDS->oOvManager.Initialize( poDS, poOpenInfo->pszFilename );
705 :
706 16 : return( poDS );
707 : }
708 :
709 : /************************************************************************/
710 : /* Create() */
711 : /************************************************************************/
712 :
713 : GDALDataset *IdrisiDataset::Create( const char *pszFilename,
714 : int nXSize,
715 : int nYSize,
716 : int nBands,
717 : GDALDataType eType,
718 32 : char **papszOptions )
719 : {
720 : (void) papszOptions;
721 :
722 : // --------------------------------------------------------------------
723 : // Check input options
724 : // --------------------------------------------------------------------
725 :
726 32 : if( nBands != 1 )
727 : {
728 17 : if( !( nBands == 3 && eType == GDT_Byte ) )
729 : {
730 : CPLError( CE_Failure, CPLE_AppDefined,
731 : "Attempt to create IDRISI dataset with an illegal "
732 : "number of bands(%d) or data type(%s).\n",
733 14 : nBands, GDALGetDataTypeName( eType ) );
734 14 : return NULL;
735 : }
736 : }
737 :
738 : // ----------------------------------------------------------------
739 : // Create the header file with minimun information
740 : // ----------------------------------------------------------------
741 :
742 : const char *pszLDataType;
743 :
744 18 : switch( eType )
745 : {
746 : case GDT_Byte:
747 5 : if( nBands == 1 )
748 2 : pszLDataType = rstBYTE;
749 : else
750 3 : pszLDataType = rstRGB24;
751 5 : break;
752 : case GDT_Int16:
753 2 : pszLDataType = rstINTEGER;
754 2 : break;
755 : case GDT_Float32:
756 3 : pszLDataType = rstREAL;
757 3 : break;
758 : default:
759 : CPLError( CE_Failure, CPLE_AppDefined,
760 : "Attempt to create IDRISI dataset with an illegal "
761 : "data type(%s).\n",
762 8 : GDALGetDataTypeName( eType ) );
763 8 : return NULL;
764 :
765 : };
766 :
767 : char **papszLRDC;
768 :
769 10 : papszLRDC = NULL;
770 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcFILE_FORMAT, rstVERSION );
771 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcFILE_TITLE, "" );
772 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcDATA_TYPE, pszLDataType );
773 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcFILE_TYPE, "binary" );
774 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcCOLUMNS, CPLSPrintf( "%d", nXSize ) );
775 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcROWS, CPLSPrintf( "%d", nYSize ) );
776 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcREF_SYSTEM, "plane" );
777 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcREF_UNITS, "m" );
778 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcUNIT_DIST, "1" );
779 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcMIN_X, "0" );
780 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcMAX_X, CPLSPrintf( "%d", nXSize ) );
781 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcMIN_Y, "0" );
782 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcMAX_Y, CPLSPrintf( "%d", nYSize ) );
783 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcPOSN_ERROR, "unspecified" );
784 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcRESOLUTION, "1.0" );
785 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcMIN_VALUE, "0" );
786 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcMAX_VALUE, "0" );
787 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcDISPLAY_MIN, "0" );
788 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcDISPLAY_MAX, "0" );
789 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcVALUE_UNITS, "unspecified" );
790 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcVALUE_ERROR, "unspecified" );
791 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcFLAG_VALUE, "none" );
792 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcFLAG_DEFN, "none" );
793 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcLEGEND_CATS, "0" );
794 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcLINEAGES, "" );
795 10 : papszLRDC = CSLAddNameValue( papszLRDC, rdcCOMMENTS, "" );
796 :
797 : const char *pszLDocFilename;
798 10 : pszLDocFilename = CPLResetExtension( pszFilename, extRDC );
799 :
800 10 : CSLSetNameValueSeparator( papszLRDC, ": " );
801 10 : SaveAsCRLF( papszLRDC, pszLDocFilename );
802 10 : CSLDestroy( papszLRDC );
803 :
804 : // ----------------------------------------------------------------
805 : // Create an empty data file
806 : // ----------------------------------------------------------------
807 :
808 : FILE *fp;
809 :
810 10 : fp = VSIFOpenL( pszFilename, "wb+" );
811 :
812 10 : if( fp == NULL )
813 : {
814 : CPLError( CE_Failure, CPLE_OpenFailed,
815 0 : "Attempt to create file %s' failed.\n", pszFilename );
816 0 : return NULL;
817 : }
818 10 : VSIFCloseL( fp );
819 :
820 10 : return (IdrisiDataset *) GDALOpen( pszFilename, GA_Update );
821 : }
822 :
823 : /************************************************************************/
824 : /* CreateCopy() */
825 : /************************************************************************/
826 :
827 : GDALDataset *IdrisiDataset::CreateCopy( const char *pszFilename,
828 : GDALDataset *poSrcDS,
829 : int bStrict,
830 : char **papszOptions,
831 : GDALProgressFunc pfnProgress,
832 17 : void *pProgressData )
833 : {
834 17 : if( !pfnProgress( 0.0, NULL, pProgressData ) )
835 0 : return NULL;
836 :
837 : // ------------------------------------------------------------------------
838 : // Check number of bands
839 : // ------------------------------------------------------------------------
840 :
841 17 : if ( ( ( poSrcDS->GetRasterCount() == 1 ) ||
842 : ( ( poSrcDS->GetRasterCount() == 3 ) &&
843 : ( ( poSrcDS->GetRasterBand( 1 )->GetRasterDataType() == GDT_Byte ) &&
844 : ( poSrcDS->GetRasterBand( 2 )->GetRasterDataType() == GDT_Byte ) &&
845 : ( poSrcDS->GetRasterBand( 3 )->GetRasterDataType() == GDT_Byte ) ) ) ) == FALSE )
846 : {
847 :
848 : CPLError( CE_Failure, CPLE_AppDefined,
849 : "Attempt to create IDRISI dataset with an illegal "
850 : "number of bands(%d).\n",
851 4 : poSrcDS->GetRasterCount() );
852 4 : return NULL;
853 : }
854 :
855 : // ------------------------------------------------------------------------
856 : // Check Data types
857 : // ------------------------------------------------------------------------
858 :
859 : int i;
860 :
861 22 : for( i = 1; i <= poSrcDS->GetRasterCount(); i++ )
862 : {
863 17 : GDALDataType eType = poSrcDS->GetRasterBand( i )->GetRasterDataType();
864 :
865 17 : if( bStrict )
866 : {
867 14 : if( eType != GDT_Byte &&
868 : eType != GDT_Int16 &&
869 : eType != GDT_Float32 )
870 : {
871 : CPLError( CE_Failure, CPLE_AppDefined,
872 : "Attempt to create IDRISI dataset in strict mode "
873 : "with an illegal data type(%s).\n",
874 8 : GDALGetDataTypeName( eType ) );
875 8 : return NULL;
876 : }
877 : }
878 : else
879 : {
880 3 : if( eType != GDT_Byte &&
881 : eType != GDT_Int16 &&
882 : eType != GDT_UInt16 &&
883 : eType != GDT_UInt32 &&
884 : eType != GDT_Int32 &&
885 : eType != GDT_Float32 &&
886 : eType != GDT_Float64 )
887 : {
888 : CPLError( CE_Failure, CPLE_AppDefined,
889 : "Attempt to create IDRISI dataset with an illegal "
890 : "data type(%s).\n",
891 0 : GDALGetDataTypeName( eType ) );
892 0 : return NULL;
893 : }
894 : }
895 : }
896 :
897 : // --------------------------------------------------------------------
898 : // Define data type
899 : // --------------------------------------------------------------------
900 :
901 5 : GDALRasterBand *poBand = poSrcDS->GetRasterBand( 1 );
902 5 : GDALDataType eType = poBand->GetRasterDataType();
903 :
904 5 : int bSuccessMin = FALSE;
905 5 : int bSuccessMax = FALSE;
906 :
907 : double dfMin;
908 : double dfMax;
909 :
910 5 : dfMin = poBand->GetMinimum( &bSuccessMin );
911 5 : dfMax = poBand->GetMaximum( &bSuccessMax );
912 :
913 5 : if( ! ( bSuccessMin && bSuccessMax ) )
914 : {
915 5 : poBand->GetStatistics( false, false, &dfMin, &dfMax, NULL, NULL );
916 : }
917 :
918 5 : if(!( ( eType == GDT_Byte ) ||
919 : ( eType == GDT_Int16 ) ||
920 : ( eType == GDT_Float32 ) ) )
921 : {
922 0 : if( eType == GDT_Float64 )
923 : {
924 0 : eType = GDT_Float32;
925 : }
926 : else
927 : {
928 0 : if( ( dfMin < (double) SHRT_MIN ) ||
929 : ( dfMax > (double) SHRT_MAX ) )
930 : {
931 0 : eType = GDT_Float32;
932 : }
933 : else
934 : {
935 0 : eType = GDT_Int16;
936 : }
937 : }
938 : }
939 :
940 : // --------------------------------------------------------------------
941 : // Create the dataset
942 : // --------------------------------------------------------------------
943 :
944 : IdrisiDataset *poDS;
945 :
946 : poDS = (IdrisiDataset *) IdrisiDataset::Create( pszFilename,
947 : poSrcDS->GetRasterXSize(),
948 : poSrcDS->GetRasterYSize(),
949 : poSrcDS->GetRasterCount(),
950 : eType,
951 5 : papszOptions );
952 :
953 5 : if( poDS == NULL )
954 4 : return NULL;
955 :
956 : // --------------------------------------------------------------------
957 : // Copy information to the dataset
958 : // --------------------------------------------------------------------
959 :
960 : double adfGeoTransform[6];
961 :
962 1 : poDS->SetProjection( poSrcDS->GetProjectionRef() );
963 :
964 1 : if( poSrcDS->GetGeoTransform( adfGeoTransform ) == CE_None )
965 : {
966 1 : poDS->SetGeoTransform( adfGeoTransform );
967 : }
968 : else
969 : {
970 0 : CSLSetNameValue( poDS->papszRDC, rdcMIN_X, CPLSPrintf( "%.1f", -0.5 ));
971 : CSLSetNameValue( poDS->papszRDC, rdcMAX_X, CPLSPrintf( "%.1f",
972 0 : poSrcDS->GetRasterXSize() - 0.5 ) );
973 : CSLSetNameValue( poDS->papszRDC, rdcMIN_Y, CPLSPrintf( "%.1f",
974 0 : poSrcDS->GetRasterYSize() + 0.5 ) );
975 0 : CSLSetNameValue( poDS->papszRDC, rdcMAX_Y, CPLSPrintf( "%.1f", 0.5 ));
976 0 : CSLSetNameValue( poDS->papszRDC, rdcRESOLUTION, "1" );
977 : }
978 :
979 : // --------------------------------------------------------------------
980 : // Copy information to the raster band
981 : // --------------------------------------------------------------------
982 :
983 : GDALRasterBand *poSrcBand;
984 : int bHasNoDataValue;
985 : double dfNoDataValue;
986 :
987 4 : for( i = 1; i <= poDS->nBands; i++ )
988 : {
989 3 : poSrcBand = poSrcDS->GetRasterBand( i );
990 3 : poBand = poDS->GetRasterBand( i );
991 :
992 3 : if( poDS->nBands == 1 )
993 : {
994 0 : poBand->SetUnitType( poSrcBand->GetUnitType() );
995 0 : poBand->SetColorTable( poSrcBand->GetColorTable() );
996 0 : poBand->SetCategoryNames( poSrcBand->GetCategoryNames() );
997 :
998 0 : const GDALRasterAttributeTable *poRAT = poSrcBand->GetDefaultRAT();
999 :
1000 0 : if( poRAT != NULL )
1001 : {
1002 0 : poBand->SetDefaultRAT( poRAT );
1003 : }
1004 : }
1005 :
1006 3 : dfMin = poSrcBand->GetMinimum( NULL );
1007 3 : dfMax = poSrcBand->GetMaximum( NULL );
1008 3 : poBand->SetStatistics( dfMin, dfMax, 0.0, 0.0 );
1009 3 : dfNoDataValue = poSrcBand->GetNoDataValue( &bHasNoDataValue );
1010 3 : if( bHasNoDataValue )
1011 0 : poBand->SetNoDataValue( dfNoDataValue );
1012 : }
1013 :
1014 : // --------------------------------------------------------------------
1015 : // Copy image data
1016 : // --------------------------------------------------------------------
1017 :
1018 : GDALDatasetCopyWholeRaster( (GDALDatasetH) poSrcDS,
1019 : (GDALDatasetH) poDS, NULL,
1020 1 : pfnProgress, pProgressData );
1021 :
1022 : // --------------------------------------------------------------------
1023 : // Finalize
1024 : // --------------------------------------------------------------------
1025 :
1026 1 : poDS->FlushCache();
1027 :
1028 1 : return poDS;
1029 : }
1030 :
1031 : /************************************************************************/
1032 : /* GetFileList() */
1033 : /************************************************************************/
1034 :
1035 2 : char **IdrisiDataset::GetFileList()
1036 : {
1037 2 : char **papszFileList = GDALPamDataset::GetFileList();
1038 : const char *pszAssociated;
1039 :
1040 : // --------------------------------------------------------------------
1041 : // Symbol table file
1042 : // --------------------------------------------------------------------
1043 :
1044 2 : pszAssociated = CPLResetExtension( pszFilename, extSMP );
1045 :
1046 2 : if( FileExists( pszAssociated ) )
1047 : {
1048 0 : papszFileList = CSLAddString( papszFileList, pszAssociated );
1049 : }
1050 : else
1051 : {
1052 2 : pszAssociated = CPLResetExtension( pszFilename, extSMPu );
1053 :
1054 2 : if( FileExists( pszAssociated ) )
1055 : {
1056 0 : papszFileList = CSLAddString( papszFileList, pszAssociated );
1057 : }
1058 : }
1059 :
1060 : // --------------------------------------------------------------------
1061 : // Documentation file
1062 : // --------------------------------------------------------------------
1063 :
1064 2 : pszAssociated = CPLResetExtension( pszFilename, extRDC );
1065 :
1066 2 : if( FileExists( pszAssociated ) )
1067 : {
1068 2 : papszFileList = CSLAddString( papszFileList, pszAssociated );
1069 : }
1070 : else
1071 : {
1072 0 : pszAssociated = CPLResetExtension( pszFilename, extRDCu );
1073 :
1074 0 : if( FileExists( pszAssociated ) )
1075 : {
1076 0 : papszFileList = CSLAddString( papszFileList, pszAssociated );
1077 : }
1078 : }
1079 :
1080 : // --------------------------------------------------------------------
1081 : // Reference file
1082 : // --------------------------------------------------------------------
1083 :
1084 2 : pszAssociated = CPLResetExtension( pszFilename, extREF );
1085 :
1086 2 : if( FileExists( pszAssociated ) )
1087 : {
1088 0 : papszFileList = CSLAddString( papszFileList, pszAssociated );
1089 : }
1090 : else
1091 : {
1092 2 : pszAssociated = CPLResetExtension( pszFilename, extREFu );
1093 :
1094 2 : if( FileExists( pszAssociated ) )
1095 : {
1096 0 : papszFileList = CSLAddString( papszFileList, pszAssociated );
1097 : }
1098 : }
1099 :
1100 2 : return papszFileList;
1101 : }
1102 :
1103 : /************************************************************************/
1104 : /* GetGeoTransform() */
1105 : /************************************************************************/
1106 :
1107 1 : CPLErr IdrisiDataset::GetGeoTransform( double * padfTransform )
1108 : {
1109 1 : if( GDALPamDataset::GetGeoTransform( padfTransform ) != CE_None )
1110 : {
1111 1 : memcpy( padfTransform, adfGeoTransform, sizeof( double ) * 6 );
1112 :
1113 1 : if( adfGeoTransform[0] == 0.0
1114 : && adfGeoTransform[1] == 1.0
1115 : && adfGeoTransform[2] == 0.0
1116 : && adfGeoTransform[3] == 0.0
1117 : && adfGeoTransform[4] == 0.0
1118 : && adfGeoTransform[5] == 1.0 )
1119 0 : return CE_Failure;
1120 : }
1121 :
1122 1 : return CE_None;
1123 : }
1124 :
1125 : /************************************************************************/
1126 : /* SetGeoTransform() */
1127 : /************************************************************************/
1128 :
1129 5 : CPLErr IdrisiDataset::SetGeoTransform( double * padfTransform )
1130 : {
1131 5 : if( padfTransform[2] != 0.0 || padfTransform[4] != 0.0 )
1132 : {
1133 : CPLError( CE_Failure, CPLE_AppDefined,
1134 : "Attempt to set rotated geotransform on Idrisi Raster file.\n"
1135 0 : "Idrisi Raster does not support rotation.\n" );
1136 0 : return CE_Failure;
1137 : }
1138 :
1139 : // --------------------------------------------------------------------
1140 : // Update the .rdc file
1141 : // --------------------------------------------------------------------
1142 :
1143 : double dfMinX, dfMaxX, dfMinY, dfMaxY, dfXPixSz, dfYPixSz;
1144 :
1145 5 : dfXPixSz = padfTransform[1];
1146 5 : dfYPixSz = padfTransform[5];
1147 5 : dfMinX = padfTransform[0];
1148 5 : dfMaxX = ( dfXPixSz * nRasterXSize ) + dfMinX;
1149 :
1150 5 : if( dfYPixSz < 0 )
1151 : {
1152 5 : dfMaxY = padfTransform[3];
1153 5 : dfMinY = ( dfYPixSz * nRasterYSize ) + padfTransform[3];
1154 : }
1155 : else
1156 : {
1157 0 : dfMaxY = ( dfYPixSz * nRasterYSize ) + padfTransform[3];
1158 0 : dfMinY = padfTransform[3];
1159 : }
1160 :
1161 5 : CSLSetNameValue( papszRDC, rdcMIN_X, CPLSPrintf( "%.7f", dfMinX ) );
1162 5 : CSLSetNameValue( papszRDC, rdcMAX_X, CPLSPrintf( "%.7f", dfMaxX ) );
1163 5 : CSLSetNameValue( papszRDC, rdcMIN_Y, CPLSPrintf( "%.7f", dfMinY ) );
1164 5 : CSLSetNameValue( papszRDC, rdcMAX_Y, CPLSPrintf( "%.7f", dfMaxY ) );
1165 5 : CSLSetNameValue( papszRDC, rdcRESOLUTION, CPLSPrintf( "%.7f", fabs( dfYPixSz ) ) );
1166 :
1167 : // --------------------------------------------------------------------
1168 : // Update the Dataset attribute
1169 : // --------------------------------------------------------------------
1170 :
1171 5 : memcpy( adfGeoTransform, padfTransform, sizeof( double ) * 6 );
1172 :
1173 5 : return CE_None;
1174 : }
1175 :
1176 : /************************************************************************/
1177 : /* GetProjectionRef() */
1178 : /************************************************************************/
1179 :
1180 1 : const char *IdrisiDataset::GetProjectionRef( void )
1181 : {
1182 1 : const char *pszPamSRS = GDALPamDataset::GetProjectionRef();
1183 :
1184 1 : if( pszPamSRS != NULL && strlen( pszPamSRS ) > 0 )
1185 0 : return pszPamSRS;
1186 :
1187 1 : if( pszProjection == NULL )
1188 : {
1189 1 : const char *pszRefSystem = CSLFetchNameValue( papszRDC, rdcREF_SYSTEM );
1190 1 : const char *pszRefUnit = CSLFetchNameValue( papszRDC, rdcREF_UNITS );
1191 :
1192 1 : GeoReference2Wkt( pszRefSystem, pszRefUnit, &pszProjection );
1193 : }
1194 1 : return pszProjection;
1195 : }
1196 :
1197 : /************************************************************************/
1198 : /* SetProjection() */
1199 : /************************************************************************/
1200 :
1201 5 : CPLErr IdrisiDataset::SetProjection( const char *pszProjString )
1202 : {
1203 5 : CPLFree( pszProjection );
1204 5 : pszProjection = CPLStrdup( pszProjString );
1205 5 : CPLErr eResult = CE_None;
1206 :
1207 5 : char *pszRefSystem = NULL;
1208 5 : char *pszRefUnit = NULL;
1209 :
1210 5 : eResult = Wkt2GeoReference( pszProjString, &pszRefSystem, &pszRefUnit );
1211 :
1212 5 : CSLSetNameValue( papszRDC, rdcREF_SYSTEM, pszRefSystem );
1213 5 : CSLSetNameValue( papszRDC, rdcREF_UNITS, pszRefUnit );
1214 :
1215 5 : CPLFree( pszRefSystem );
1216 5 : CPLFree( pszRefUnit );
1217 :
1218 5 : return eResult;
1219 : }
1220 :
1221 : /************************************************************************/
1222 : /* IdrisiRasterBand() */
1223 : /************************************************************************/
1224 :
1225 : IdrisiRasterBand::IdrisiRasterBand( IdrisiDataset *poDS,
1226 : int nBand,
1227 26 : GDALDataType eDataType )
1228 : {
1229 26 : this->poDS = poDS;
1230 26 : this->nBand = nBand;
1231 26 : this->eDataType = eDataType;
1232 26 : this->poDefaultRAT = NULL;
1233 26 : this->dfMinimum = 0.0;
1234 26 : this->dfMaximum = 0.0;
1235 26 : this->bFirstVal = true;
1236 :
1237 : // --------------------------------------------------------------------
1238 : // Set Dimension
1239 : // --------------------------------------------------------------------
1240 :
1241 26 : nBlockYSize = 1;
1242 26 : nBlockXSize = poDS->GetRasterXSize();
1243 :
1244 : // --------------------------------------------------------------------
1245 : // Get ready for reading and writing
1246 : // --------------------------------------------------------------------
1247 :
1248 26 : nRecordSize = poDS->GetRasterXSize() * GDALGetDataTypeSize( eDataType ) / 8;
1249 :
1250 26 : pabyScanLine = (GByte*) VSIMalloc2( nRecordSize, poDS->nBands );
1251 :
1252 26 : if( pabyScanLine == NULL )
1253 : {
1254 : CPLError(CE_Failure, CPLE_OutOfMemory,
1255 : "IdrisiRasterBand::IdrisiRasterBand : Out of memory (nRasterXSize = %d)",
1256 0 : poDS->GetRasterXSize());
1257 : }
1258 :
1259 26 : nRecordSize *= poDS->nBands;
1260 26 : }
1261 :
1262 : /************************************************************************/
1263 : /* ~IdrisiRasterBand() */
1264 : /************************************************************************/
1265 :
1266 26 : IdrisiRasterBand::~IdrisiRasterBand()
1267 : {
1268 26 : CPLFree( pabyScanLine );
1269 :
1270 26 : if( poDefaultRAT )
1271 : {
1272 0 : delete poDefaultRAT;
1273 : }
1274 26 : }
1275 :
1276 : /************************************************************************/
1277 : /* IReadBlock() */
1278 : /************************************************************************/
1279 :
1280 : CPLErr IdrisiRasterBand::IReadBlock( int nBlockXOff,
1281 : int nBlockYOff,
1282 162 : void *pImage )
1283 : {
1284 162 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1285 :
1286 162 : if( VSIFSeekL( poGDS->fp,
1287 : vsi_l_offset(nRecordSize) * nBlockYOff, SEEK_SET ) < 0 )
1288 : {
1289 : CPLError( CE_Failure, CPLE_FileIO,
1290 : "Can't seek(%s) block with X offset %d and Y offset %d.\n%s",
1291 0 : poGDS->pszFilename, nBlockXOff, nBlockYOff, VSIStrerror( errno ) );
1292 0 : return CE_Failure;
1293 : }
1294 :
1295 162 : if( (int) VSIFReadL( pabyScanLine, 1, nRecordSize, poGDS->fp ) < nRecordSize )
1296 : {
1297 : CPLError( CE_Failure, CPLE_FileIO,
1298 : "Can't read(%s) block with X offset %d and Y offset %d.\n%s",
1299 0 : poGDS->pszFilename, nBlockXOff, nBlockYOff, VSIStrerror( errno ) );
1300 0 : return CE_Failure;
1301 : }
1302 :
1303 162 : if( poGDS->nBands == 3 )
1304 : {
1305 : int i, j;
1306 5100 : for( i = 0, j = ( 3 - nBand ); i < nBlockXSize; i++, j += 3 )
1307 : {
1308 5000 : ( (GByte*) pImage )[i] = pabyScanLine[j];
1309 : }
1310 : }
1311 : else
1312 : {
1313 62 : memcpy( pImage, pabyScanLine, nRecordSize );
1314 : }
1315 :
1316 : #ifdef CPL_MSB
1317 : if( eDataType == GDT_Float32 )
1318 : GDALSwapWords( pImage, 4, nBlockXSize * nBlockYSize, 4 );
1319 : #endif
1320 :
1321 162 : return CE_None;
1322 : }
1323 :
1324 : /************************************************************************/
1325 : /* IWriteBlock() */
1326 : /************************************************************************/
1327 :
1328 : CPLErr IdrisiRasterBand::IWriteBlock( int nBlockXOff,
1329 : int nBlockYOff,
1330 170 : void *pImage )
1331 : {
1332 170 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1333 :
1334 : #ifdef CPL_MSB
1335 : // Swap in input buffer if needed.
1336 : if( eDataType == GDT_Float32 )
1337 : GDALSwapWords( pImage, 4, nBlockXSize * nBlockYSize, 4 );
1338 : #endif
1339 :
1340 170 : if( poGDS->nBands == 1 )
1341 : {
1342 20 : memcpy( pabyScanLine, pImage, nRecordSize );
1343 : }
1344 : else
1345 : {
1346 150 : if( nBand > 1 )
1347 : {
1348 : VSIFSeekL( poGDS->fp,
1349 100 : vsi_l_offset(nRecordSize) * nBlockYOff, SEEK_SET );
1350 100 : VSIFReadL( pabyScanLine, 1, nRecordSize, poGDS->fp );
1351 : }
1352 : int i, j;
1353 7650 : for( i = 0, j = ( 3 - nBand ); i < nBlockXSize; i++, j += 3 )
1354 : {
1355 7500 : pabyScanLine[j] = ( (GByte *) pImage )[i];
1356 : }
1357 : }
1358 :
1359 : #ifdef CPL_MSB
1360 : // Swap input buffer back to original form.
1361 : if( eDataType == GDT_Float32 )
1362 : GDALSwapWords( pImage, 4, nBlockXSize * nBlockYSize, 4 );
1363 : #endif
1364 :
1365 170 : VSIFSeekL( poGDS->fp, vsi_l_offset(nRecordSize) * nBlockYOff, SEEK_SET );
1366 :
1367 170 : if( (int) VSIFWriteL( pabyScanLine, 1, nRecordSize, poGDS->fp ) < nRecordSize )
1368 : {
1369 : CPLError( CE_Failure, CPLE_FileIO,
1370 : "Can't write(%s) block with X offset %d and Y offset %d.\n%s",
1371 0 : poGDS->pszFilename, nBlockXOff, nBlockYOff, VSIStrerror( errno ) );
1372 0 : return CE_Failure;
1373 : }
1374 :
1375 : // --------------------------------------------------------------------
1376 : // Store the first value as mimunm and maximum
1377 : // --------------------------------------------------------------------
1378 :
1379 170 : if( bFirstVal )
1380 : {
1381 4 : switch( eDataType )
1382 : {
1383 : case GDT_Float32:
1384 1 : dfMinimum = dfMaximum = (float) ((float*) pabyScanLine)[0];
1385 1 : break;
1386 : case GDT_Int16:
1387 0 : dfMinimum = dfMaximum = (float) ((GInt16*) pabyScanLine)[0];
1388 0 : break;
1389 : default:
1390 : dfMinimum = dfMaximum = (float) ((GByte*) pabyScanLine)[
1391 3 : poGDS->nBands == 1 ? 0 : 3 - this->nBand ];
1392 : }
1393 4 : bFirstVal = false;
1394 : }
1395 :
1396 : // --------------------------------------------------------------------
1397 : // Search for the minimum and maximum values
1398 : // --------------------------------------------------------------------
1399 :
1400 : int i, j;
1401 :
1402 170 : if( eDataType == GDT_Float32 )
1403 : {
1404 420 : for( i = 0; i < nBlockXSize; i++ )
1405 : {
1406 400 : dfMinimum = MIN( dfMinimum, ((float*) pabyScanLine)[i] );
1407 400 : dfMaximum = MAX( dfMaximum, ((float*) pabyScanLine)[i] );
1408 : }
1409 : }
1410 150 : else if( eDataType == GDT_Int16 )
1411 : {
1412 0 : for( i = 0; i < nBlockXSize; i++ )
1413 : {
1414 0 : dfMinimum = MIN( dfMinimum, ((GInt16*) pabyScanLine)[i] );
1415 0 : dfMaximum = MAX( dfMaximum, ((GInt16*) pabyScanLine)[i] );
1416 : }
1417 : }
1418 150 : else if( poGDS->nBands == 1 )
1419 : {
1420 0 : for( i = 0; i < nBlockXSize; i++ )
1421 : {
1422 0 : dfMinimum = MIN( dfMinimum, ((GByte*) pabyScanLine)[i] );
1423 0 : dfMaximum = MAX( dfMaximum, ((GByte*) pabyScanLine)[i] );
1424 : }
1425 : }
1426 : else
1427 : {
1428 7650 : for( i = 0, j = ( 3 - nBand ); i < nBlockXSize; i++, j += 3 )
1429 : {
1430 7500 : dfMinimum = MIN( dfMinimum, ((GByte*) pabyScanLine)[j] );
1431 7500 : dfMaximum = MAX( dfMaximum, ((GByte*) pabyScanLine)[j] );
1432 : }
1433 : }
1434 :
1435 170 : return CE_None;
1436 : }
1437 :
1438 : /************************************************************************/
1439 : /* GetMinimum() */
1440 : /************************************************************************/
1441 :
1442 0 : double IdrisiRasterBand::GetMinimum( int *pbSuccess )
1443 : {
1444 0 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1445 :
1446 : double adfMinValue[3];
1447 : sscanf( CSLFetchNameValue( poGDS->papszRDC, rdcMIN_VALUE ), "%lf %lf %lf",
1448 0 : &adfMinValue[0], &adfMinValue[1], &adfMinValue[2] );
1449 :
1450 0 : if( pbSuccess )
1451 : {
1452 0 : *pbSuccess = true;
1453 : }
1454 :
1455 0 : return adfMinValue[this->nBand - 1];
1456 : }
1457 :
1458 : /************************************************************************/
1459 : /* GetMaximum() */
1460 : /************************************************************************/
1461 :
1462 0 : double IdrisiRasterBand::GetMaximum( int *pbSuccess )
1463 : {
1464 0 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1465 :
1466 : double adfMaxValue[3];
1467 : sscanf( CSLFetchNameValue( poGDS->papszRDC, rdcMAX_VALUE ), "%lf %lf %lf",
1468 0 : &adfMaxValue[0], &adfMaxValue[1], &adfMaxValue[2] );
1469 :
1470 0 : if( pbSuccess )
1471 : {
1472 0 : *pbSuccess = true;
1473 : }
1474 :
1475 0 : return adfMaxValue[this->nBand - 1];
1476 : }
1477 :
1478 : /************************************************************************/
1479 : /* GetNoDataValue() */
1480 : /************************************************************************/
1481 :
1482 4 : double IdrisiRasterBand::GetNoDataValue( int *pbSuccess )
1483 : {
1484 4 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1485 :
1486 : double dfNoData;
1487 : const char *pszFlagDefn;
1488 :
1489 4 : if( CSLFetchNameValue( poGDS->papszRDC, rdcFLAG_DEFN ) != NULL )
1490 4 : pszFlagDefn = CSLFetchNameValue( poGDS->papszRDC, rdcFLAG_DEFN );
1491 0 : else if( CSLFetchNameValue( poGDS->papszRDC, rdcFLAG_DEFN2 ) != NULL )
1492 0 : pszFlagDefn = CSLFetchNameValue( poGDS->papszRDC, rdcFLAG_DEFN2 );
1493 : else
1494 0 : pszFlagDefn = CPLStrdup( "none" );
1495 :
1496 : // ------------------------------------------------------------------------
1497 : // If Flag_Def is not "none", Flag_Value means "background"
1498 : // or "missing data"
1499 : // ------------------------------------------------------------------------
1500 :
1501 4 : if( ! EQUAL( pszFlagDefn, "none" ) )
1502 : {
1503 0 : dfNoData = atof_nz( CSLFetchNameValue( poGDS->papszRDC, rdcFLAG_VALUE ) );
1504 0 : if( pbSuccess )
1505 0 : *pbSuccess = TRUE;
1506 : }
1507 : else
1508 : {
1509 4 : dfNoData = -9999.0; /* this value should be ignored */
1510 4 : if( pbSuccess )
1511 4 : *pbSuccess = FALSE;
1512 : }
1513 :
1514 4 : return dfNoData;
1515 : }
1516 :
1517 : /************************************************************************/
1518 : /* SetNoDataValue() */
1519 : /************************************************************************/
1520 :
1521 0 : CPLErr IdrisiRasterBand::SetNoDataValue( double dfNoDataValue )
1522 : {
1523 0 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1524 :
1525 : poGDS->papszRDC =
1526 0 : CSLSetNameValue( poGDS->papszRDC, rdcFLAG_VALUE, CPLSPrintf( "%.7g", dfNoDataValue ) );
1527 : poGDS->papszRDC =
1528 0 : CSLSetNameValue( poGDS->papszRDC, rdcFLAG_DEFN, "missing data" );
1529 :
1530 0 : return CE_None;
1531 : }
1532 :
1533 : /************************************************************************/
1534 : /* GetColorInterpretation() */
1535 : /************************************************************************/
1536 :
1537 0 : GDALColorInterp IdrisiRasterBand::GetColorInterpretation()
1538 : {
1539 0 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1540 :
1541 0 : if( poGDS->nBands == 3 )
1542 : {
1543 0 : switch( nBand )
1544 : {
1545 0 : case 1: return GCI_BlueBand;
1546 0 : case 2: return GCI_GreenBand;
1547 0 : case 3: return GCI_RedBand;
1548 : }
1549 : }
1550 0 : else if( poGDS->poColorTable->GetColorEntryCount() > 0 )
1551 : {
1552 0 : return GCI_PaletteIndex;
1553 : }
1554 0 : return GCI_GrayIndex;
1555 : }
1556 :
1557 : /************************************************************************/
1558 : /* GetCategoryNames() */
1559 : /************************************************************************/
1560 :
1561 0 : char **IdrisiRasterBand::GetCategoryNames()
1562 : {
1563 0 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1564 :
1565 0 : return poGDS->papszCategories;
1566 : }
1567 :
1568 : /************************************************************************/
1569 : /* SetCategoryNames() */
1570 : /************************************************************************/
1571 :
1572 0 : CPLErr IdrisiRasterBand::SetCategoryNames( char **papszCategoryNames )
1573 : {
1574 0 : int nCatCount = CSLCount( papszCategoryNames );
1575 :
1576 0 : if( nCatCount == 0 )
1577 0 : return CE_None;
1578 :
1579 0 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1580 :
1581 0 : CSLDestroy( poGDS->papszCategories );
1582 0 : poGDS->papszCategories = CSLDuplicate( papszCategoryNames );
1583 :
1584 : // ------------------------------------------------------
1585 : // Search for the "Legend cats : N" line
1586 : // ------------------------------------------------------
1587 :
1588 0 : int i, nLine = -1;
1589 0 : for( i = 0;( i < CSLCount( poGDS->papszRDC ) ) &&( nLine == -1 ); i++ )
1590 0 : if( EQUALN( poGDS->papszRDC[i], rdcLEGEND_CATS, 12 ) )
1591 0 : nLine = i;
1592 :
1593 0 : if( nLine < 0 )
1594 0 : return CE_None;
1595 :
1596 0 : int nCount = atoi_nz( CSLFetchNameValue( poGDS->papszRDC, rdcLEGEND_CATS ) );
1597 :
1598 : // ------------------------------------------------------
1599 : // Delte old instance of the categoty names
1600 : // ------------------------------------------------------
1601 :
1602 0 : if( nCount > 0 )
1603 0 : poGDS->papszRDC = CSLRemoveStrings( poGDS->papszRDC, nLine + 1, nCount, NULL );
1604 :
1605 0 : nCount = 0;
1606 :
1607 0 : for( i = 0; i < nCatCount; i++ )
1608 : {
1609 0 : if( ( strlen( papszCategoryNames[i] ) > 0 ) )
1610 : {
1611 : poGDS->papszRDC = CSLInsertString( poGDS->papszRDC,( nLine + nCount + 1 ),
1612 0 : CPLSPrintf( "%s:%s", CPLSPrintf( rdcCODE_N, i ), papszCategoryNames[i] ) );
1613 0 : nCount++;
1614 : }
1615 : }
1616 :
1617 0 : CSLSetNameValue( poGDS->papszRDC, rdcLEGEND_CATS, CPLSPrintf( "%d", nCount ) );
1618 :
1619 0 : return CE_None;
1620 : }
1621 :
1622 : /************************************************************************/
1623 : /* GetColorTable() */
1624 : /************************************************************************/
1625 :
1626 0 : GDALColorTable *IdrisiRasterBand::GetColorTable()
1627 : {
1628 0 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1629 :
1630 0 : if( poGDS->poColorTable->GetColorEntryCount() == 0 )
1631 : {
1632 0 : return NULL;
1633 : }
1634 : else
1635 : {
1636 0 : return poGDS->poColorTable;
1637 : }
1638 : }
1639 :
1640 : /************************************************************************/
1641 : /* SetColorTable() */
1642 : /************************************************************************/
1643 :
1644 0 : CPLErr IdrisiRasterBand::SetColorTable( GDALColorTable *poColorTable )
1645 : {
1646 0 : if( poColorTable == NULL )
1647 : {
1648 0 : return CE_None;
1649 : }
1650 :
1651 0 : if( poColorTable->GetColorEntryCount() == 0 )
1652 : {
1653 0 : return CE_None;
1654 : }
1655 :
1656 0 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1657 :
1658 0 : delete poGDS->poColorTable;
1659 :
1660 0 : poGDS->poColorTable = poColorTable->Clone();
1661 :
1662 : const char *pszSMPFilename;
1663 0 : pszSMPFilename = CPLResetExtension( poGDS->pszFilename, extSMP );
1664 : FILE *fpSMP;
1665 :
1666 0 : if( ( fpSMP = VSIFOpenL( pszSMPFilename, "w" ) ) != NULL )
1667 : {
1668 0 : VSIFWriteL( "[Idrisi]", 8, 1, fpSMP );
1669 0 : GByte nPlatform = 1; VSIFWriteL( &nPlatform, 1, 1, fpSMP );
1670 0 : GByte nVersion = 11; VSIFWriteL( &nVersion, 1, 1, fpSMP );
1671 0 : GByte nDepth = 8; VSIFWriteL( &nDepth, 1, 1, fpSMP );
1672 0 : GByte nHeadSz = 18; VSIFWriteL( &nHeadSz, 1, 1, fpSMP );
1673 0 : GUInt16 nCount = 255; VSIFWriteL( &nCount, 2, 1, fpSMP );
1674 0 : GUInt16 nMix = 0; VSIFWriteL( &nMix, 2, 1, fpSMP );
1675 0 : GUInt16 nMax = 255; VSIFWriteL( &nMax, 2, 1, fpSMP );
1676 :
1677 : GDALColorEntry oEntry;
1678 : GByte aucRGB[3];
1679 : int i;
1680 :
1681 0 : for( i = 0; i < poColorTable->GetColorEntryCount(); i++ )
1682 : {
1683 0 : poColorTable->GetColorEntryAsRGB( i, &oEntry );
1684 0 : aucRGB[0] = (GByte) oEntry.c1;
1685 0 : aucRGB[1] = (GByte) oEntry.c2;
1686 0 : aucRGB[2] = (GByte) oEntry.c3;
1687 0 : VSIFWriteL( &aucRGB, 3, 1, fpSMP );
1688 : }
1689 : /* smp files always have 256 occurences */
1690 0 : for( i = poColorTable->GetColorEntryCount(); i <= 255; i++ )
1691 : {
1692 0 : poColorTable->GetColorEntryAsRGB( i, &oEntry );
1693 0 : aucRGB[0] = (GByte) 0;
1694 0 : aucRGB[1] = (GByte) 0;
1695 0 : aucRGB[2] = (GByte) 0;
1696 0 : VSIFWriteL( &aucRGB, 3, 1, fpSMP );
1697 : }
1698 0 : VSIFCloseL( fpSMP );
1699 : }
1700 :
1701 0 : return CE_None;
1702 : }
1703 :
1704 : /************************************************************************/
1705 : /* GetUnitType() */
1706 : /************************************************************************/
1707 :
1708 0 : const char *IdrisiRasterBand::GetUnitType()
1709 : {
1710 0 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1711 :
1712 0 : return poGDS->pszUnitType;
1713 : }
1714 :
1715 : /************************************************************************/
1716 : /* SetUnitType() */
1717 : /************************************************************************/
1718 :
1719 0 : CPLErr IdrisiRasterBand::SetUnitType( const char *pszUnitType )
1720 : {
1721 0 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1722 :
1723 0 : if( strlen( pszUnitType ) == 0 )
1724 : {
1725 : poGDS->papszRDC =
1726 0 : CSLSetNameValue( poGDS->papszRDC, rdcVALUE_UNITS, "unspecified" );
1727 : }
1728 : else
1729 : {
1730 : poGDS->papszRDC =
1731 0 : CSLSetNameValue( poGDS->papszRDC, rdcVALUE_UNITS, pszUnitType );
1732 : }
1733 :
1734 0 : return CE_None;
1735 : }
1736 :
1737 : /************************************************************************/
1738 : /* SetStatistics() */
1739 : /************************************************************************/
1740 :
1741 13 : CPLErr IdrisiRasterBand::SetStatistics( double dfMin, double dfMax, double dfMean, double dfStdDev )
1742 : {
1743 13 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1744 :
1745 13 : double adfMin[3] = {0.0, 0.0, 0.0};
1746 13 : double adfMax[3] = {0.0, 0.0, 0.0};
1747 :
1748 13 : sscanf( CSLFetchNameValue( poGDS->papszRDC, rdcMIN_VALUE ), "%lf %lf %lf", &adfMin[0], &adfMin[1], &adfMin[2] );
1749 13 : sscanf( CSLFetchNameValue( poGDS->papszRDC, rdcMAX_VALUE ), "%lf %lf %lf", &adfMax[0], &adfMax[1], &adfMax[2] );
1750 :
1751 13 : adfMin[nBand - 1] = dfMin;
1752 13 : adfMax[nBand - 1] = dfMax;
1753 :
1754 13 : if( poGDS->nBands == 3 )
1755 : {
1756 : poGDS->papszRDC =
1757 9 : CSLSetNameValue( poGDS->papszRDC, rdcMIN_VALUE, CPLSPrintf( "%.8g %.8g %.8g", adfMin[0], adfMin[1], adfMin[2] ) );
1758 : poGDS->papszRDC =
1759 9 : CSLSetNameValue( poGDS->papszRDC, rdcMAX_VALUE, CPLSPrintf( "%.8g %.8g %.8g", adfMax[0], adfMax[1], adfMax[2] ) );
1760 : poGDS->papszRDC =
1761 9 : CSLSetNameValue( poGDS->papszRDC, rdcDISPLAY_MIN, CPLSPrintf( "%.8g %.8g %.8g", adfMin[0], adfMin[1], adfMin[2] ) );
1762 : poGDS->papszRDC =
1763 9 : CSLSetNameValue( poGDS->papszRDC, rdcDISPLAY_MAX, CPLSPrintf( "%.8g %.8g %.8g", adfMax[0], adfMax[1], adfMax[2] ) );
1764 : }
1765 : else
1766 : {
1767 : poGDS->papszRDC =
1768 4 : CSLSetNameValue( poGDS->papszRDC, rdcMIN_VALUE, CPLSPrintf( "%.8g", adfMin[0] ) );
1769 : poGDS->papszRDC =
1770 4 : CSLSetNameValue( poGDS->papszRDC, rdcMAX_VALUE, CPLSPrintf( "%.8g", adfMax[0] ) );
1771 : poGDS->papszRDC =
1772 4 : CSLSetNameValue( poGDS->papszRDC, rdcDISPLAY_MIN, CPLSPrintf( "%.8g", adfMin[0] ) );
1773 : poGDS->papszRDC =
1774 4 : CSLSetNameValue( poGDS->papszRDC, rdcDISPLAY_MAX, CPLSPrintf( "%.8g", adfMax[0] ) );
1775 : }
1776 :
1777 13 : return GDALRasterBand::SetStatistics( dfMin, dfMax, dfMean, dfStdDev );
1778 : }
1779 :
1780 : /************************************************************************/
1781 : /* SetDefaultRAT() */
1782 : /************************************************************************/
1783 :
1784 0 : CPLErr IdrisiRasterBand::SetDefaultRAT( const GDALRasterAttributeTable *poRAT )
1785 : {
1786 0 : if( ! poRAT )
1787 : {
1788 0 : return CE_Failure;
1789 : }
1790 :
1791 : // ----------------------------------------------------------
1792 : // Get field indecies
1793 : // ----------------------------------------------------------
1794 :
1795 0 : int iValue = -1;
1796 0 : int iRed = poRAT->GetColOfUsage( GFU_Red );
1797 0 : int iGreen = poRAT->GetColOfUsage( GFU_Green );
1798 0 : int iBlue = poRAT->GetColOfUsage( GFU_Blue );
1799 0 : int iName = -1;
1800 : int i;
1801 :
1802 0 : GDALColorTable *poCT = NULL;
1803 0 : char **papszNames = NULL;
1804 :
1805 0 : int nFact = 1;
1806 :
1807 : // ----------------------------------------------------------
1808 : // Seek for "Value" field index (AGIS standards field name)
1809 : // ----------------------------------------------------------
1810 :
1811 0 : if( GetColorTable() == 0 || GetColorTable()->GetColorEntryCount() == 0 )
1812 : {
1813 0 : for( i = 0; i < poRAT->GetColumnCount(); i++ )
1814 : {
1815 0 : if EQUALN( "Value", poRAT->GetNameOfCol( i ), 5 )
1816 : {
1817 0 : iValue = i;
1818 0 : break;
1819 : }
1820 : }
1821 :
1822 0 : if( iRed != -1 && iGreen != -1 && iBlue != -1 )
1823 : {
1824 0 : poCT = new GDALColorTable();
1825 0 : nFact = poRAT->GetTypeOfCol( iRed ) == GFT_Real ? 255 : 1;
1826 : }
1827 : }
1828 :
1829 : // ----------------------------------------------------------
1830 : // Seek for Name field index
1831 : // ----------------------------------------------------------
1832 :
1833 0 : if( GetCategoryNames() == 0 || CSLCount( GetCategoryNames() ) == 0 )
1834 : {
1835 0 : iName = poRAT->GetColOfUsage( GFU_Name );
1836 0 : if( iName == -1 )
1837 : {
1838 0 : for( i = 0; i < poRAT->GetColumnCount(); i++ )
1839 : {
1840 0 : if EQUALN( "Class_Name", poRAT->GetNameOfCol( i ), 10 )
1841 : {
1842 0 : iName = i;
1843 0 : break;
1844 : }
1845 0 : else if EQUALN( "Categor", poRAT->GetNameOfCol( i ), 7 )
1846 : {
1847 0 : iName = i;
1848 0 : break;
1849 : }
1850 0 : else if EQUALN( "Name", poRAT->GetNameOfCol( i ), 4 )
1851 : {
1852 0 : iName = i;
1853 0 : break;
1854 : }
1855 : }
1856 : }
1857 :
1858 : /* if still can't find it use the first String column */
1859 :
1860 0 : if( iName == -1 )
1861 : {
1862 0 : for( i = 0; i < poRAT->GetColumnCount(); i++ )
1863 : {
1864 0 : if( poRAT->GetTypeOfCol( i ) == GFT_String )
1865 : {
1866 0 : iName = i;
1867 0 : break;
1868 : }
1869 : }
1870 : }
1871 :
1872 : // ----------------------------------------------------------
1873 : // Incomplete Attribute Table;
1874 : // ----------------------------------------------------------
1875 :
1876 0 : if( iName == -1 )
1877 : {
1878 0 : iName = iValue;
1879 : }
1880 : }
1881 :
1882 : // ----------------------------------------------------------
1883 : // Initialization
1884 : // ----------------------------------------------------------
1885 :
1886 0 : double dRed = 0.0;
1887 0 : double dGreen = 0.0;
1888 0 : double dBlue = 0.0;
1889 :
1890 : // ----------------------------------------------------------
1891 : // Load values
1892 : // ----------------------------------------------------------
1893 :
1894 : GDALColorEntry sColor;
1895 0 : int iEntry = 0;
1896 0 : int iOut = 0;
1897 0 : int nEntryCount = poRAT->GetRowCount();
1898 0 : int nValue = 0;
1899 :
1900 0 : if( iValue != -1 )
1901 : {
1902 0 : nValue = poRAT->GetValueAsInt( iEntry, iValue );
1903 : }
1904 :
1905 0 : for( iOut = 0; iOut < 65535 && ( iEntry < nEntryCount ); iOut++ )
1906 : {
1907 0 : if( iOut == nValue )
1908 : {
1909 0 : if( poCT )
1910 : {
1911 0 : dRed = poRAT->GetValueAsDouble( iEntry, iRed );
1912 0 : dGreen = poRAT->GetValueAsDouble( iEntry, iGreen );
1913 0 : dBlue = poRAT->GetValueAsDouble( iEntry, iBlue );
1914 0 : sColor.c1 = (short) ( dRed * nFact );
1915 0 : sColor.c2 = (short) ( dGreen * nFact );
1916 0 : sColor.c3 = (short) ( dBlue * nFact );
1917 0 : sColor.c4 = (short) ( 255 / nFact );
1918 0 : poCT->SetColorEntry( iEntry, &sColor );
1919 : }
1920 :
1921 0 : if( iName != -1 )
1922 : {
1923 : papszNames = CSLAddString( papszNames,
1924 0 : poRAT->GetValueAsString( iEntry, iName ) );
1925 : }
1926 :
1927 : /* Advance on the table */
1928 :
1929 0 : if( ( ++iEntry ) < nEntryCount )
1930 : {
1931 0 : if( iValue != -1 )
1932 0 : nValue = poRAT->GetValueAsInt( iEntry, iValue );
1933 : else
1934 0 : nValue = iEntry;
1935 : }
1936 : }
1937 0 : else if( iOut < nValue )
1938 : {
1939 0 : if( poCT )
1940 : {
1941 0 : sColor.c1 = (short) 0;
1942 0 : sColor.c2 = (short) 0;
1943 0 : sColor.c3 = (short) 0;
1944 0 : sColor.c4 = (short) 255;
1945 0 : poCT->SetColorEntry( iEntry, &sColor );
1946 : }
1947 :
1948 0 : if( iName != -1 )
1949 0 : papszNames = CSLAddString( papszNames, "" );
1950 : }
1951 : }
1952 :
1953 : // ----------------------------------------------------------
1954 : // Set Color Table
1955 : // ----------------------------------------------------------
1956 :
1957 0 : if( poCT )
1958 : {
1959 0 : SetColorTable( poCT );
1960 0 : delete poCT;
1961 : }
1962 :
1963 : // ----------------------------------------------------------
1964 : // Update Category Names
1965 : // ----------------------------------------------------------
1966 :
1967 0 : if( papszNames )
1968 : {
1969 0 : SetCategoryNames( papszNames );
1970 0 : CSLDestroy( papszNames );
1971 : }
1972 :
1973 : // ----------------------------------------------------------
1974 : // Update Attribute Table
1975 : // ----------------------------------------------------------
1976 :
1977 0 : if( poDefaultRAT )
1978 : {
1979 0 : delete poDefaultRAT;
1980 : }
1981 :
1982 0 : poDefaultRAT = poRAT->Clone();
1983 :
1984 0 : return CE_None;
1985 : }
1986 :
1987 : /************************************************************************/
1988 : /* GetDefaultRAT() */
1989 : /************************************************************************/
1990 :
1991 0 : const GDALRasterAttributeTable *IdrisiRasterBand::GetDefaultRAT()
1992 : {
1993 0 : IdrisiDataset *poGDS = (IdrisiDataset *) poDS;
1994 :
1995 0 : if( poGDS->papszCategories == NULL )
1996 : {
1997 0 : return NULL;
1998 : }
1999 :
2000 0 : bool bHasColorTable = poGDS->poColorTable->GetColorEntryCount() > 0;
2001 :
2002 : // ----------------------------------------------------------
2003 : // Create the bands Attribute Table
2004 : // ----------------------------------------------------------
2005 :
2006 0 : if( poDefaultRAT )
2007 : {
2008 0 : delete poDefaultRAT;
2009 : }
2010 :
2011 0 : poDefaultRAT = new GDALRasterAttributeTable();
2012 :
2013 : // ----------------------------------------------------------
2014 : // Create (Value, Red, Green, Blue, Alpha, Class_Name) fields
2015 : // ----------------------------------------------------------
2016 :
2017 0 : poDefaultRAT->CreateColumn( "Value", GFT_Integer, GFU_Generic );
2018 0 : poDefaultRAT->CreateColumn( "Value_1", GFT_Integer, GFU_MinMax );
2019 :
2020 0 : if( bHasColorTable )
2021 : {
2022 0 : poDefaultRAT->CreateColumn( "Red", GFT_Integer, GFU_Red );
2023 0 : poDefaultRAT->CreateColumn( "Green", GFT_Integer, GFU_Green );
2024 0 : poDefaultRAT->CreateColumn( "Blue", GFT_Integer, GFU_Blue );
2025 0 : poDefaultRAT->CreateColumn( "Alpha", GFT_Integer, GFU_Alpha );
2026 : }
2027 0 : poDefaultRAT->CreateColumn( "Class_name", GFT_String, GFU_Name );
2028 :
2029 : // ----------------------------------------------------------
2030 : // Loop throught the Category Names
2031 : // ----------------------------------------------------------
2032 :
2033 : GDALColorEntry sEntry;
2034 0 : int iName = poDefaultRAT->GetColOfUsage( GFU_Name );
2035 0 : int nEntryCount = CSLCount( poGDS->papszCategories );
2036 0 : int iEntry = 0;
2037 0 : int iRows = 0;
2038 :
2039 0 : for( iEntry = 0; iEntry < nEntryCount; iEntry++ )
2040 : {
2041 0 : if EQUAL( poGDS->papszCategories[iEntry], "" )
2042 : {
2043 0 : continue; // Eliminate the empty ones
2044 : }
2045 0 : poDefaultRAT->SetRowCount( poDefaultRAT->GetRowCount() + 1 );
2046 0 : poDefaultRAT->SetValue( iRows, 0, iEntry );
2047 0 : poDefaultRAT->SetValue( iRows, 1, iEntry );
2048 0 : if( bHasColorTable )
2049 : {
2050 0 : poGDS->poColorTable->GetColorEntryAsRGB( iEntry, &sEntry );
2051 0 : poDefaultRAT->SetValue( iRows, 2, sEntry.c1 );
2052 0 : poDefaultRAT->SetValue( iRows, 3, sEntry.c2 );
2053 0 : poDefaultRAT->SetValue( iRows, 4, sEntry.c3 );
2054 0 : poDefaultRAT->SetValue( iRows, 5, sEntry.c4 );
2055 : }
2056 0 : poDefaultRAT->SetValue( iRows, iName, poGDS->papszCategories[iEntry] );
2057 0 : iRows++;
2058 : }
2059 :
2060 0 : return poDefaultRAT;
2061 : }
2062 :
2063 : /************************************************************************/
2064 : /* GeoReference2Wkt() */
2065 : /************************************************************************/
2066 :
2067 : /***
2068 : * Converts Idrisi geographic reference information to OpenGIS WKT.
2069 : *
2070 : * The Idrisi metadata file contain two fields that describe the
2071 : * geographic reference, RefSystem and RefUnit.
2072 : *
2073 : * RefSystem can contains the world "plane" or the name of a georeference
2074 : * file <refsystem>.ref that details the geographic reference
2075 : * system( coordinate system and projection parameters ). RefUnits
2076 : * indicates the unit of the image bounds.
2077 : *
2078 : * The georeference files are generally located in the product installation
2079 : * folder $IDRISIDIR\Georef, but they are first looked for in the same
2080 : * folder as the data file.
2081 : *
2082 : * If a Reference system names can be recognized by a name convention
2083 : * it will be interpreted without the need to read the georeference file.
2084 : * That includes "latlong" and all the UTM and State Plane zones.
2085 : *
2086 : * RefSystem "latlong" means that the data is not project and the coordinate
2087 : * system is WGS84. RefSystem "plane" means that the there is no coordinate
2088 : * system but the it is possible to calculate areas and distance by looking
2089 : * at the RefUnits.
2090 : *
2091 : * If the environment variable IDRISIDIR is not set and the georeference file
2092 : * need to be read then the projection string will result as unknown.
2093 : ***/
2094 :
2095 : CPLErr IdrisiDataset::GeoReference2Wkt( const char *pszRefSystem,
2096 : const char *pszRefUnits,
2097 1 : char **pszProjString )
2098 : {
2099 1 : OGRSpatialReference oSRS;
2100 :
2101 : // ---------------------------------------------------------
2102 : // Plane
2103 : // ---------------------------------------------------------
2104 :
2105 1 : if( EQUAL( pszRefSystem, rstPLANE ) )
2106 : {
2107 0 : oSRS.SetLocalCS( "Plane" );
2108 0 : int nUnit = GetUnitIndex( pszRefUnits );
2109 0 : if( nUnit > -1 )
2110 : {
2111 0 : int nDeft = aoLinearUnitsConv[nUnit].nDefaultG;
2112 : oSRS.SetLinearUnits( aoLinearUnitsConv[nDeft].pszName,
2113 0 : aoLinearUnitsConv[nDeft].dfConv );
2114 : }
2115 0 : oSRS.exportToWkt( pszProjString );
2116 0 : return CE_None;
2117 : }
2118 :
2119 : // ---------------------------------------------------------
2120 : // Latlong
2121 : // ---------------------------------------------------------
2122 :
2123 1 : if( EQUAL( pszRefSystem, rstLATLONG ) ||
2124 : EQUAL( pszRefSystem, rstLATLONG2 ) )
2125 : {
2126 1 : oSRS.SetWellKnownGeogCS( "WGS84" );
2127 1 : oSRS.exportToWkt( pszProjString );
2128 1 : return CE_None;
2129 : }
2130 :
2131 : // ---------------------------------------------------------
2132 : // Prepare for scanning in lower case
2133 : // ---------------------------------------------------------
2134 :
2135 : char *pszRefSystemLower;
2136 0 : pszRefSystemLower = CPLStrdup( pszRefSystem );
2137 0 : CPLStrlwr( pszRefSystemLower );
2138 :
2139 : // ---------------------------------------------------------
2140 : // UTM naming convention( ex.: utm-30n )
2141 : // ---------------------------------------------------------
2142 :
2143 0 : if( EQUALN( pszRefSystem, rstUTM, 3 ) )
2144 : {
2145 : int nZone;
2146 : char cNorth;
2147 0 : sscanf( pszRefSystemLower, rstUTM, &nZone, &cNorth );
2148 0 : oSRS.SetWellKnownGeogCS( "WGS84" );
2149 0 : oSRS.SetUTM( nZone,( cNorth == 'n' ) );
2150 0 : oSRS.exportToWkt( pszProjString );
2151 0 : CPLFree( pszRefSystemLower );
2152 0 : return CE_None;
2153 : }
2154 :
2155 : // ---------------------------------------------------------
2156 : // State Plane naming convention( ex.: spc83ma1 )
2157 : // ---------------------------------------------------------
2158 :
2159 0 : if( EQUALN( pszRefSystem, rstSPC, 3 ) )
2160 : {
2161 : int nNAD;
2162 : int nZone;
2163 : char szState[3];
2164 0 : sscanf( pszRefSystemLower, rstSPC, &nNAD, szState, &nZone );
2165 0 : int nSPCode = GetStateCode( szState );
2166 0 : if( nSPCode != -1 )
2167 : {
2168 0 : nZone = ( nZone == 1 ? nSPCode : nSPCode + nZone - 1 );
2169 :
2170 0 : if( oSRS.SetStatePlane( nZone, ( nNAD == 83 ) ) != OGRERR_FAILURE )
2171 : {
2172 0 : oSRS.exportToWkt( pszProjString );
2173 0 : CPLFree( pszRefSystemLower );
2174 0 : return CE_None;
2175 : }
2176 :
2177 : // ----------------------------------------------------------
2178 : // If SetStatePlane fails, set GeoCS as NAD Datum and let it
2179 : // try to read the projection info from georeference file( * )
2180 : // ----------------------------------------------------------
2181 :
2182 0 : oSRS.SetWellKnownGeogCS( CPLSPrintf( "NAD%d", nNAD ) );
2183 : }
2184 : }
2185 :
2186 0 : CPLFree( pszRefSystemLower );
2187 0 : pszRefSystemLower = NULL;
2188 :
2189 : // ------------------------------------------------------------------
2190 : // Search for georeference file <RefSystem>.ref
2191 : // ------------------------------------------------------------------
2192 :
2193 : const char *pszFName = CPLSPrintf( "%s%c%s.ref",
2194 0 : CPLGetDirname( pszFilename ), PATHDELIM, pszRefSystem );
2195 :
2196 0 : if( ! FileExists( pszFName ) )
2197 : {
2198 : // ------------------------------------------------------------------
2199 : // Look at $IDRISIDIR\Georef\<RefSystem>.ref
2200 : // ------------------------------------------------------------------
2201 :
2202 0 : const char *pszIdrisiDir = CPLGetConfigOption( "IDRISIDIR", NULL );
2203 :
2204 0 : if( ( pszIdrisiDir ) != NULL )
2205 : {
2206 : pszFName = CPLSPrintf( "%s%cgeoref%c%s.ref",
2207 0 : pszIdrisiDir, PATHDELIM, PATHDELIM, pszRefSystem );
2208 : }
2209 : }
2210 :
2211 : // ------------------------------------------------------------------
2212 : // Cannot find georeference file
2213 : // ------------------------------------------------------------------
2214 :
2215 0 : if( ! FileExists( pszFName ) )
2216 : {
2217 : CPLDebug( "RST", "Cannot find Idrisi georeference file %s",
2218 0 : pszRefSystem );
2219 :
2220 0 : if( oSRS.IsGeographic() == FALSE ) /* see State Plane remarks( * ) */
2221 : {
2222 0 : oSRS.SetLocalCS( "Unknown" );
2223 0 : int nUnit = GetUnitIndex( pszRefUnits );
2224 0 : if( nUnit > -1 )
2225 : {
2226 0 : int nDeft = aoLinearUnitsConv[nUnit].nDefaultG;
2227 : oSRS.SetLinearUnits( aoLinearUnitsConv[nDeft].pszName,
2228 0 : aoLinearUnitsConv[nDeft].dfConv );
2229 : }
2230 : }
2231 0 : oSRS.exportToWkt( pszProjString );
2232 0 : return CE_Failure;
2233 : }
2234 :
2235 : // ------------------------------------------------------------------
2236 : // Read values from georeference file
2237 : // ------------------------------------------------------------------
2238 :
2239 0 : char **papszRef = CSLLoad( pszFName );
2240 0 : CSLSetNameValueSeparator( papszRef, ":" );
2241 :
2242 : char *pszGeorefName;
2243 :
2244 0 : if( EQUAL( CSLFetchNameValue( papszRef, refREF_SYSTEM ), "" ) == FALSE )
2245 : {
2246 0 : pszGeorefName = CPLStrdup( CSLFetchNameValue( papszRef, refREF_SYSTEM ) );
2247 : }
2248 : else
2249 : {
2250 0 : pszGeorefName = CPLStrdup( CSLFetchNameValue( papszRef, refREF_SYSTEM2 ) );
2251 : }
2252 0 : char *pszProjName = CPLStrdup( CSLFetchNameValue( papszRef, refPROJECTION ) );
2253 0 : char *pszDatum = CPLStrdup( CSLFetchNameValue( papszRef, refDATUM ) );
2254 0 : char *pszEllipsoid = CPLStrdup( CSLFetchNameValue( papszRef, refELLIPSOID ) );
2255 0 : double dfCenterLat = atof_nz( CSLFetchNameValue( papszRef, refORIGIN_LAT ) );
2256 0 : double dfCenterLong = atof_nz( CSLFetchNameValue( papszRef, refORIGIN_LONG ) );
2257 0 : double dfSemiMajor = atof_nz( CSLFetchNameValue( papszRef, refMAJOR_SAX ) );
2258 0 : double dfSemiMinor = atof_nz( CSLFetchNameValue( papszRef, refMINOR_SAX ) );
2259 0 : double dfFalseEasting = atof_nz( CSLFetchNameValue( papszRef, refORIGIN_X ) );
2260 0 : double dfFalseNorthing = atof_nz( CSLFetchNameValue( papszRef, refORIGIN_Y ) );
2261 0 : double dfStdP1 = atof_nz( CSLFetchNameValue( papszRef, refSTANDL_1 ) );
2262 0 : double dfStdP2 = atof_nz( CSLFetchNameValue( papszRef, refSTANDL_2 ) );
2263 : double dfScale;
2264 : double adfToWGS84[3];
2265 :
2266 : sscanf( CSLFetchNameValue( papszRef, refDELTA_WGS84 ), "%lf %lf %lf",
2267 0 : &adfToWGS84[0], &adfToWGS84[1], &adfToWGS84[2] );
2268 :
2269 0 : if( EQUAL( CSLFetchNameValue( papszRef, refSCALE_FAC ), "na" ) )
2270 0 : dfScale = 1.0;
2271 : else
2272 0 : dfScale = atof_nz( CSLFetchNameValue( papszRef, refSCALE_FAC ) );
2273 :
2274 0 : CSLDestroy( papszRef );
2275 :
2276 : // ----------------------------------------------------------------------
2277 : // Set the Geographic Coordinate System
2278 : // ----------------------------------------------------------------------
2279 :
2280 0 : if( oSRS.IsGeographic() == FALSE ) /* see State Plane remarks(*) */
2281 : {
2282 0 : int nEPSG = 0;
2283 :
2284 : // ----------------------------------------------------------------------
2285 : // Is it a WGS84 equivalent?
2286 : // ----------------------------------------------------------------------
2287 :
2288 0 : if( ( EQUALN( pszEllipsoid, "WGS", 3 ) ) &&( strstr( pszEllipsoid, "84" ) ) &&
2289 : ( EQUALN( pszDatum, "WGS", 3 ) ) &&( strstr( pszDatum, "84" ) ) &&
2290 : ( adfToWGS84[0] == 0.0 ) &&( adfToWGS84[1] == 0.0 ) &&( adfToWGS84[2] == 0.0 ) )
2291 : {
2292 0 : nEPSG = 4326;
2293 : }
2294 :
2295 : // ----------------------------------------------------------------------
2296 : // Match GCS's DATUM_NAME by using 'ApproxString' over Datum
2297 : // ----------------------------------------------------------------------
2298 :
2299 0 : if( nEPSG == 0 )
2300 : {
2301 0 : nEPSG = atoi_nz( CSVGetField( CSVFilename( "gcs.csv" ),
2302 : "DATUM_NAME", pszDatum, CC_ApproxString, "COORD_REF_SYS_CODE" ) );
2303 : }
2304 :
2305 : // ----------------------------------------------------------------------
2306 : // Match GCS's COORD_REF_SYS_NAME by using 'ApproxString' over Datum
2307 : // ----------------------------------------------------------------------
2308 :
2309 0 : if( nEPSG == 0 )
2310 : {
2311 0 : nEPSG = atoi_nz( CSVGetField( CSVFilename( "gcs.csv" ),
2312 : "COORD_REF_SYS_NAME", pszDatum, CC_ApproxString, "COORD_REF_SYS_CODE" ) );
2313 : }
2314 :
2315 0 : if( nEPSG != 0 )
2316 : {
2317 0 : oSRS.importFromEPSG( nEPSG );
2318 : }
2319 : else
2320 : {
2321 : // --------------------------------------------------
2322 : // Create GeogCS based on the georeference file info
2323 : // --------------------------------------------------
2324 :
2325 : oSRS.SetGeogCS( pszRefSystem,
2326 : pszDatum,
2327 : pszEllipsoid,
2328 : dfSemiMajor,
2329 0 : ( -1.0 /( dfSemiMinor / dfSemiMajor - 1.0 ) ) );
2330 : }
2331 :
2332 : // ----------------------------------------------------------------------
2333 : // Note: That will override EPSG info:
2334 : // ----------------------------------------------------------------------
2335 :
2336 0 : oSRS.SetTOWGS84( adfToWGS84[0], adfToWGS84[1], adfToWGS84[2] );
2337 : }
2338 :
2339 : // ----------------------------------------------------------------------
2340 : // If the georeference file tells that it is a non project system:
2341 : // ----------------------------------------------------------------------
2342 :
2343 0 : if( EQUAL( pszProjName, "none" ) )
2344 : {
2345 0 : oSRS.exportToWkt( pszProjString );
2346 :
2347 0 : CPLFree( pszGeorefName );
2348 0 : CPLFree( pszProjName );
2349 0 : CPLFree( pszDatum );
2350 0 : CPLFree( pszEllipsoid );
2351 :
2352 0 : return CE_None;
2353 : }
2354 :
2355 : // ----------------------------------------------------------------------
2356 : // Create Projection information based on georeference file info
2357 : // ----------------------------------------------------------------------
2358 :
2359 : // Idrisi user's Manual, Supported Projection:
2360 : //
2361 : // Mercator
2362 : // Transverse Mercator
2363 : // Gauss-Kruger
2364 : // Lambert Conformal Conic
2365 : // Plate Carree
2366 : // Hammer Aitoff
2367 : // Lambert North Polar Azimuthal Equal Area
2368 : // Lambert South Polar Azimuthal Equal Area
2369 : // Lambert Transverse Azimuthal Equal Area
2370 : // Lambert Oblique Polar Azimuthal Equal Area
2371 : // North Polar Stereographic
2372 : // South Polar Stereographic
2373 : // Transverse Stereographic
2374 : // Oblique Stereographic
2375 : // Albers Equal Area Conic
2376 : // Sinusoidal
2377 : //
2378 :
2379 0 : if( EQUAL( pszProjName, "Mercator" ) )
2380 : {
2381 0 : oSRS.SetMercator( dfCenterLat, dfCenterLong, dfScale, dfFalseEasting, dfFalseNorthing );
2382 : }
2383 0 : else if EQUAL( pszProjName, "Transverse Mercator" )
2384 : {
2385 0 : oSRS.SetTM( dfCenterLat, dfCenterLong, dfScale, dfFalseEasting, dfFalseNorthing );
2386 : }
2387 0 : else if EQUALN( pszProjName, "Gauss-Kruger", 9 )
2388 : {
2389 0 : oSRS.SetTM( dfCenterLat, dfCenterLong, dfScale, dfFalseEasting, dfFalseNorthing );
2390 : }
2391 0 : else if( EQUAL( pszProjName, "Lambert Conformal Conic" ) )
2392 : {
2393 0 : oSRS.SetLCC( dfStdP1, dfStdP2, dfCenterLat, dfCenterLong, dfFalseEasting, dfFalseNorthing );
2394 : }
2395 0 : else if( EQUALN( pszProjName, "Plate Carr" "\xE9" "e", 10 ) ) /* 'eacute' in ISO-8859-1 */
2396 : {
2397 0 : oSRS.SetEquirectangular( dfCenterLat, dfCenterLong, dfFalseEasting, dfFalseNorthing );
2398 : }
2399 0 : else if( EQUAL( pszProjName, "Hammer Aitoff" ) )
2400 : {
2401 0 : oSRS.SetProjection( pszProjName );
2402 0 : oSRS.SetProjParm( SRS_PP_LATITUDE_OF_ORIGIN, dfCenterLat );
2403 0 : oSRS.SetProjParm( SRS_PP_CENTRAL_MERIDIAN, dfCenterLong );
2404 0 : oSRS.SetProjParm( SRS_PP_FALSE_EASTING, dfFalseEasting );
2405 0 : oSRS.SetProjParm( SRS_PP_FALSE_NORTHING, dfFalseNorthing );
2406 : }
2407 0 : else if( EQUALN( pszProjName, "Lambert North Polar Azimuthal Equal Area", 15 ) ||
2408 : EQUALN( pszProjName, "Lambert South Polar Azimuthal Equal Area", 15 ) ||
2409 : EQUALN( pszProjName, "Lambert Transverse Azimuthal Equal Area", 15 ) ||
2410 : EQUALN( pszProjName, "Lambert Oblique Polar Azimuthal Equal Area", 15 ) )
2411 : {
2412 0 : oSRS.SetLAEA( dfCenterLat, dfCenterLong, dfFalseEasting, dfFalseNorthing );
2413 : }
2414 0 : else if( EQUALN( pszProjName, "North Polar Stereographic", 15 ) ||
2415 : EQUALN( pszProjName, "South Polar Stereographic", 15 ) )
2416 : {
2417 0 : oSRS.SetPS( dfCenterLat, dfCenterLong, dfScale, dfFalseEasting, dfFalseNorthing );
2418 : }
2419 0 : else if( EQUALN( pszProjName, "Transverse Stereographic", 15 ) )
2420 : {
2421 0 : oSRS.SetStereographic( dfCenterLat, dfCenterLong, dfScale, dfFalseEasting, dfFalseNorthing );
2422 : }
2423 0 : else if( EQUALN( pszProjName, "Oblique Stereographic", 15 ) )
2424 : {
2425 0 : oSRS.SetOS( dfCenterLat, dfCenterLong, dfScale, dfFalseEasting, dfFalseNorthing );
2426 : }
2427 0 : else if( EQUAL( pszProjName, "Alber's Equal Area Conic" ) ||
2428 : EQUAL( pszProjName, "Albers Equal Area Conic" ) )
2429 : {
2430 0 : oSRS.SetACEA( dfStdP1, dfStdP2, dfCenterLat, dfCenterLong, dfFalseEasting, dfFalseNorthing );
2431 : }
2432 0 : else if( EQUAL( pszProjName, "Sinusoidal" ) )
2433 : {
2434 0 : oSRS.SetSinusoidal( dfCenterLong, dfFalseEasting, dfFalseNorthing );
2435 : }
2436 : else
2437 : {
2438 : CPLError( CE_Warning, CPLE_NotSupported,
2439 : "Projection not listed on Idrisi User's Manual( v.15.0/2005 ).\n\t"
2440 : "[\"%s\" in georeference file \"%s\"]",
2441 0 : pszProjName, pszFName );
2442 0 : oSRS.Clear();
2443 0 : oSRS.exportToWkt( pszProjString );
2444 :
2445 0 : CPLFree( pszGeorefName );
2446 0 : CPLFree( pszProjName );
2447 0 : CPLFree( pszDatum );
2448 0 : CPLFree( pszEllipsoid );
2449 :
2450 0 : return CE_Warning;
2451 : }
2452 :
2453 : // ----------------------------------------------------------------------
2454 : // Set the Linear Units
2455 : // ----------------------------------------------------------------------
2456 :
2457 0 : int nUnit = GetUnitIndex( pszRefUnits );
2458 0 : if( nUnit > -1 )
2459 : {
2460 0 : int nDeft = aoLinearUnitsConv[nUnit].nDefaultG;
2461 : oSRS.SetLinearUnits( aoLinearUnitsConv[nDeft].pszName,
2462 0 : aoLinearUnitsConv[nDeft].dfConv );
2463 : }
2464 : else
2465 : {
2466 0 : oSRS.SetLinearUnits( "unknown", 1.0 );
2467 : }
2468 :
2469 : // ----------------------------------------------------------------------
2470 : // Name ProjCS with the name on the georeference file
2471 : // ----------------------------------------------------------------------
2472 :
2473 0 : oSRS.SetProjCS( pszGeorefName );
2474 :
2475 0 : oSRS.exportToWkt( pszProjString );
2476 :
2477 0 : CPLFree( pszGeorefName );
2478 0 : CPLFree( pszProjName );
2479 0 : CPLFree( pszDatum );
2480 0 : CPLFree( pszEllipsoid );
2481 :
2482 0 : return CE_None;
2483 : }
2484 :
2485 : /************************************************************************/
2486 : /* Wkt2GeoReference() */
2487 : /************************************************************************/
2488 :
2489 : /***
2490 : * Converts OpenGIS WKT to Idrisi geographic reference information.
2491 : *
2492 : * That function will fill up the two parameters RefSystem and RefUnit
2493 : * that goes into the Idrisi metadata. But it could also create
2494 : * a companying georeference file to the output if necessary.
2495 : *
2496 : * First it will try to identify the ProjString as Local, WGS84 or
2497 : * one of the Idrisi name convention reference systems
2498 : * otherwise, if the projection system is supported by Idrisi,
2499 : * it will create a companying georeference files.
2500 : ***/
2501 :
2502 : CPLErr IdrisiDataset::Wkt2GeoReference( const char *pszProjString,
2503 : char **pszRefSystem,
2504 5 : char **pszRefUnit )
2505 : {
2506 : // -----------------------------------------------------
2507 : // Plane with default "Meters"
2508 : // -----------------------------------------------------
2509 :
2510 5 : if EQUAL( pszProjString, "" )
2511 : {
2512 0 : *pszRefSystem = CPLStrdup( rstPLANE );
2513 0 : *pszRefUnit = CPLStrdup( rstMETER );
2514 0 : return CE_None;
2515 : }
2516 :
2517 5 : OGRSpatialReference oSRS;
2518 5 : oSRS.importFromWkt( (char **) &pszProjString );
2519 :
2520 : // -----------------------------------------------------
2521 : // Local => Plane + Linear Unit
2522 : // -----------------------------------------------------
2523 :
2524 5 : if( oSRS.IsLocal() )
2525 : {
2526 0 : *pszRefSystem = CPLStrdup( rstPLANE );
2527 : *pszRefUnit = GetUnitDefault( oSRS.GetAttrValue( "UNIT" ),
2528 0 : CPLSPrintf( "%f", oSRS.GetLinearUnits() ) );
2529 5 : return CE_None;
2530 : }
2531 :
2532 : // -----------------------------------------------------
2533 : // Test to identify WGS84 => Latlong + Angular Unit
2534 : // -----------------------------------------------------
2535 :
2536 5 : if( oSRS.IsGeographic() )
2537 : {
2538 5 : char *pszSpheroid = CPLStrdup( oSRS.GetAttrValue( "SPHEROID" ) );
2539 5 : char *pszAuthName = CPLStrdup( oSRS.GetAuthorityName( "GEOGCS" ) );
2540 5 : char *pszDatum = CPLStrdup( oSRS.GetAttrValue( "DATUM" ) );
2541 5 : int nGCSCode = -1;
2542 5 : if EQUAL( pszAuthName, "EPSG" )
2543 : {
2544 1 : nGCSCode = atoi( oSRS.GetAuthorityCode( "GEOGCS" ) );
2545 : }
2546 5 : if( ( nGCSCode == 4326 ) ||(
2547 : ( EQUALN( pszSpheroid, "WGS", 3 ) ) &&( strstr( pszSpheroid, "84" ) ) &&
2548 : ( EQUALN( pszDatum, "WGS", 3 ) ) &&( strstr( pszDatum, "84" ) ) ) )
2549 : {
2550 5 : *pszRefSystem = CPLStrdup( rstLATLONG );
2551 5 : *pszRefUnit = CPLStrdup( rstDEGREE );
2552 :
2553 5 : CPLFree( pszSpheroid );
2554 5 : CPLFree( pszAuthName );
2555 5 : CPLFree( pszDatum );
2556 :
2557 5 : return CE_None;
2558 : }
2559 :
2560 0 : CPLFree( pszSpheroid );
2561 0 : CPLFree( pszAuthName );
2562 0 : CPLFree( pszDatum );
2563 : }
2564 :
2565 : // -----------------------------------------------------
2566 : // Prepare to match some projections
2567 : // -----------------------------------------------------
2568 :
2569 0 : const char *pszProjName = oSRS.GetAttrValue( "PROJECTION" );
2570 :
2571 0 : if( pszProjName == NULL )
2572 : {
2573 0 : pszProjName = "";
2574 : }
2575 :
2576 : // -----------------------------------------------------
2577 : // Check for UTM zones
2578 : // -----------------------------------------------------
2579 :
2580 0 : if( EQUAL( pszProjName, SRS_PT_TRANSVERSE_MERCATOR ) )
2581 : {
2582 0 : int nZone = oSRS.GetUTMZone();
2583 :
2584 0 : if( ( nZone != 0 ) && ( EQUAL( oSRS.GetAttrValue( "DATUM" ), SRS_DN_WGS84 ) ) )
2585 : {
2586 0 : double dfNorth = oSRS.GetProjParm( SRS_PP_FALSE_NORTHING );
2587 0 : *pszRefSystem = CPLStrdup( CPLSPrintf( rstUTM, nZone,( dfNorth == 0.0 ? 'n' : 's' ) ) );
2588 0 : *pszRefUnit = CPLStrdup( rstMETER );
2589 0 : return CE_None;
2590 : }
2591 : }
2592 :
2593 : // -----------------------------------------------------
2594 : // Check for State Plane
2595 : // -----------------------------------------------------
2596 :
2597 : #ifndef GDAL_RST_PLUGIN
2598 :
2599 0 : if( EQUAL( pszProjName, SRS_PT_LAMBERT_CONFORMAL_CONIC_2SP ) ||
2600 : EQUAL( pszProjName, SRS_PT_TRANSVERSE_MERCATOR ) )
2601 : {
2602 : const char *pszPCSCode;
2603 0 : const char *pszID = CPLStrdup( oSRS.GetAuthorityCode( "PROJCS" ) );
2604 0 : if( strlen( pszID ) > 0 )
2605 : {
2606 : pszPCSCode = CPLStrdup( CSVGetField( CSVFilename( "stateplane.csv" ),
2607 0 : "EPSG_PCS_CODE", pszID, CC_Integer, "ID" ) );
2608 0 : if( strlen( pszPCSCode ) > 0 )
2609 : {
2610 0 : int nNADYear = 83;
2611 0 : int nZone = pszPCSCode[strlen( pszPCSCode ) - 1] - '0';
2612 0 : int nSPCode = atoi_nz( pszPCSCode );
2613 :
2614 0 : if( nZone == 0 )
2615 0 : nZone = 1;
2616 : else
2617 0 : nSPCode = nSPCode - nZone + 1;
2618 :
2619 0 : if( nSPCode > 10000 )
2620 : {
2621 0 : nNADYear = 27;
2622 0 : nSPCode -= 10000;
2623 : }
2624 0 : char *pszState = CPLStrdup( GetStateName( nSPCode ) );
2625 0 : if( ! EQUAL( pszState, "" ) )
2626 : {
2627 0 : *pszRefSystem = CPLStrdup( CPLSPrintf( rstSPC, nNADYear, pszState, nZone ) );
2628 : *pszRefUnit = GetUnitDefault( oSRS.GetAttrValue( "UNIT" ),
2629 0 : CPLSPrintf( "%f", oSRS.GetLinearUnits() ) );
2630 0 : return CE_None;
2631 : }
2632 : }
2633 : }
2634 : }
2635 :
2636 : #endif // GDAL_RST_PLUGIN
2637 :
2638 0 : const char *pszProjectionOut = NULL;
2639 :
2640 0 : if( oSRS.IsProjected() )
2641 : {
2642 : // ---------------------------------------------------------
2643 : // Check for supported projections
2644 : // ---------------------------------------------------------
2645 :
2646 0 : if EQUAL( pszProjName, SRS_PT_MERCATOR_1SP )
2647 : {
2648 0 : pszProjectionOut = "Mercator" ;
2649 : }
2650 0 : if EQUAL( pszProjName, SRS_PT_TRANSVERSE_MERCATOR )
2651 : {
2652 0 : pszProjectionOut = "Transverse Mercator" ;
2653 : }
2654 0 : else if EQUAL( pszProjName, SRS_PT_LAMBERT_CONFORMAL_CONIC_2SP )
2655 : {
2656 0 : pszProjectionOut = "Lambert Conformal Conic" ;
2657 : }
2658 0 : else if EQUAL( pszProjName, SRS_PT_EQUIRECTANGULAR )
2659 : {
2660 0 : pszProjectionOut = "Plate Carr" "\xE9" "e" ; /* 'eacute' in ISO-8859-1 */
2661 : }
2662 0 : else if EQUAL( pszProjName, SRS_PT_LAMBERT_AZIMUTHAL_EQUAL_AREA )
2663 : {
2664 0 : double dfCenterLat = oSRS.GetProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0, NULL );
2665 0 : if( dfCenterLat == 0.0 )
2666 0 : pszProjectionOut = "Lambert Transverse Azimuthal Equal Area" ;
2667 0 : else if( fabs( dfCenterLat ) == 90.0 )
2668 0 : pszProjectionOut = "Lambert Oblique Polar Azimuthal Equal Area" ;
2669 0 : else if( dfCenterLat > 0.0 )
2670 0 : pszProjectionOut = "Lambert North Oblique Azimuthal Equal Area" ;
2671 : else
2672 0 : pszProjectionOut = "Lambert South Oblique Azimuthal Equal Area" ;
2673 : }
2674 0 : else if EQUAL( pszProjName, SRS_PT_POLAR_STEREOGRAPHIC )
2675 : {
2676 0 : if( oSRS.GetProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0, NULL ) > 0 )
2677 0 : pszProjectionOut = "North Polar Stereographic" ;
2678 : else
2679 0 : pszProjectionOut = "South Polar Stereographic" ;
2680 : }
2681 0 : else if EQUAL( pszProjName, SRS_PT_STEREOGRAPHIC )
2682 : {
2683 0 : pszProjectionOut = "Transverse Stereographic" ;
2684 : }
2685 0 : else if EQUAL( pszProjName, SRS_PT_OBLIQUE_STEREOGRAPHIC )
2686 : {
2687 0 : pszProjectionOut = "Oblique Stereographic" ;
2688 : }
2689 0 : else if EQUAL( pszProjName, SRS_PT_SINUSOIDAL )
2690 : {
2691 0 : pszProjectionOut = "Sinusoidal" ;
2692 : }
2693 0 : else if EQUAL( pszProjName, SRS_PT_ALBERS_CONIC_EQUAL_AREA )
2694 : {
2695 0 : pszProjectionOut = "Alber's Equal Area Conic" ;
2696 : }
2697 :
2698 : // ---------------------------------------------------------
2699 : // Failure, Projection system not suppotted
2700 : // ---------------------------------------------------------
2701 :
2702 0 : if( pszProjectionOut == NULL )
2703 : {
2704 : CPLDebug( "RST",
2705 : "Not supported by RST driver: PROJECTION[\"%s\"]",
2706 0 : pszProjName );
2707 :
2708 0 : *pszRefSystem = CPLStrdup( rstPLANE );
2709 0 : *pszRefUnit = CPLStrdup( rstMETER );
2710 0 : return CE_Failure;
2711 : }
2712 : }
2713 : else
2714 : {
2715 0 : pszProjectionOut = "none" ;
2716 : }
2717 :
2718 : // ---------------------------------------------------------
2719 : // Prepare to write ref file
2720 : // ---------------------------------------------------------
2721 :
2722 0 : char *pszGeorefName = CPLStrdup( "Unknown" );
2723 0 : char *pszDatum = CPLStrdup( oSRS.GetAttrValue( "DATUM" ) );
2724 0 : char *pszEllipsoid = CPLStrdup( oSRS.GetAttrValue( "SPHEROID" ) );
2725 0 : double dfSemiMajor = oSRS.GetSemiMajor();
2726 0 : double dfSemiMinor = oSRS.GetSemiMinor();
2727 : double adfToWGS84[3];
2728 0 : oSRS.GetTOWGS84( adfToWGS84, 3 );
2729 :
2730 0 : double dfCenterLat = 0.0;
2731 0 : double dfCenterLong = 0.0;
2732 0 : double dfFalseNorthing = 0.0;
2733 0 : double dfFalseEasting = 0.0;
2734 0 : double dfScale = 1.0;
2735 0 : int nParameters = 0;
2736 0 : double dfStdP1 = 0.0;
2737 0 : double dfStdP2 = 0.0;
2738 0 : char *pszAngularUnit = CPLStrdup( oSRS.GetAttrValue( "GEOGCS|UNIT" ) );
2739 : char *pszLinearUnit;
2740 :
2741 0 : if( oSRS.IsProjected() )
2742 : {
2743 0 : CPLFree( pszGeorefName );
2744 0 : pszGeorefName = CPLStrdup( oSRS.GetAttrValue( "PROJCS" ) );
2745 0 : dfCenterLat = oSRS.GetProjParm( SRS_PP_LATITUDE_OF_ORIGIN, 0.0, NULL );
2746 0 : dfCenterLong = oSRS.GetProjParm( SRS_PP_CENTRAL_MERIDIAN, 0.0, NULL );
2747 0 : dfFalseNorthing = oSRS.GetProjParm( SRS_PP_FALSE_NORTHING, 0.0, NULL );
2748 0 : dfFalseEasting = oSRS.GetProjParm( SRS_PP_FALSE_EASTING, 0.0, NULL );
2749 0 : dfScale = oSRS.GetProjParm( SRS_PP_SCALE_FACTOR, 0.0, NULL );
2750 0 : dfStdP1 = oSRS.GetProjParm( SRS_PP_STANDARD_PARALLEL_1, -0.1, NULL );
2751 0 : dfStdP2 = oSRS.GetProjParm( SRS_PP_STANDARD_PARALLEL_2, -0.1, NULL );
2752 0 : if( dfStdP1 != -0.1 )
2753 : {
2754 0 : nParameters = 1;
2755 0 : if( dfStdP2 != -0.1 )
2756 0 : nParameters = 2;
2757 : }
2758 : pszLinearUnit = GetUnitDefault( oSRS.GetAttrValue( "PROJCS|UNIT" ),
2759 0 : CPLSPrintf( "%f", oSRS.GetLinearUnits() ) ) ;
2760 : }
2761 : else
2762 : {
2763 0 : pszLinearUnit = GetUnitDefault( pszAngularUnit );
2764 : }
2765 :
2766 : // ---------------------------------------------------------
2767 : // Create a compaining georeference file for this dataset
2768 : // ---------------------------------------------------------
2769 :
2770 0 : char **papszRef = NULL;
2771 0 : papszRef = CSLAddNameValue( papszRef, refREF_SYSTEM, pszGeorefName );
2772 0 : papszRef = CSLAddNameValue( papszRef, refPROJECTION, pszProjectionOut );
2773 0 : papszRef = CSLAddNameValue( papszRef, refDATUM, pszDatum );
2774 : papszRef = CSLAddNameValue( papszRef, refDELTA_WGS84, CPLSPrintf( "%.3g %.3g %.3g",
2775 0 : adfToWGS84[0], adfToWGS84[1], adfToWGS84[2] ) );
2776 0 : papszRef = CSLAddNameValue( papszRef, refELLIPSOID, pszEllipsoid );
2777 0 : papszRef = CSLAddNameValue( papszRef, refMAJOR_SAX, CPLSPrintf( "%.3f", dfSemiMajor ) );
2778 0 : papszRef = CSLAddNameValue( papszRef, refMINOR_SAX, CPLSPrintf( "%.3f", dfSemiMinor ) );
2779 0 : papszRef = CSLAddNameValue( papszRef, refORIGIN_LONG, CPLSPrintf( "%.9g", dfCenterLong ) );
2780 0 : papszRef = CSLAddNameValue( papszRef, refORIGIN_LAT, CPLSPrintf( "%.9g", dfCenterLat ) );
2781 0 : papszRef = CSLAddNameValue( papszRef, refORIGIN_X, CPLSPrintf( "%.9g", dfFalseEasting ) );
2782 0 : papszRef = CSLAddNameValue( papszRef, refORIGIN_Y, CPLSPrintf( "%.9g", dfFalseNorthing ) );
2783 0 : papszRef = CSLAddNameValue( papszRef, refSCALE_FAC, CPLSPrintf( "%.9g", dfScale ) );
2784 0 : papszRef = CSLAddNameValue( papszRef, refUNITS, pszLinearUnit );
2785 0 : papszRef = CSLAddNameValue( papszRef, refPARAMETERS, CPLSPrintf( "%1d", nParameters ) );
2786 0 : if( nParameters > 0 )
2787 0 : papszRef = CSLAddNameValue( papszRef, refSTANDL_1, CPLSPrintf( "%.9g", dfStdP1 ) );
2788 0 : if( nParameters > 1 )
2789 0 : papszRef = CSLAddNameValue( papszRef, refSTANDL_2, CPLSPrintf( "%.9g", dfStdP2 ) );
2790 0 : CSLSetNameValueSeparator( papszRef, ": " );
2791 0 : SaveAsCRLF( papszRef, CPLResetExtension( pszFilename, extREF ) );
2792 0 : CSLDestroy( papszRef );
2793 :
2794 0 : *pszRefSystem = CPLStrdup( CPLGetBasename( pszFilename ) );
2795 0 : *pszRefUnit = CPLStrdup( pszLinearUnit );
2796 :
2797 0 : CPLFree( pszGeorefName );
2798 0 : CPLFree( pszDatum );
2799 0 : CPLFree( pszEllipsoid );
2800 0 : CPLFree( pszLinearUnit );
2801 0 : CPLFree( pszAngularUnit );
2802 :
2803 0 : return CE_None;
2804 : }
2805 :
2806 : /************************************************************************/
2807 : /* FileExists() */
2808 : /************************************************************************/
2809 :
2810 26 : bool FileExists( const char *pszPath )
2811 : {
2812 : VSIStatBuf sStat;
2813 :
2814 26 : return (bool) ( CPLStat( pszPath, &sStat ) == 0 );
2815 : }
2816 :
2817 : /************************************************************************/
2818 : /* GetStateCode() */
2819 : /************************************************************************/
2820 :
2821 0 : int GetStateCode( const char *pszState )
2822 : {
2823 : unsigned int i;
2824 :
2825 0 : for( i = 0; i < US_STATE_COUNT; i++ )
2826 : {
2827 0 : if EQUAL( pszState, aoUSStateTable[i].pszName )
2828 : {
2829 0 : return aoUSStateTable[i].nCode;
2830 : }
2831 : }
2832 0 : return -1;
2833 : }
2834 :
2835 : /************************************************************************/
2836 : /* GetStateCode() */
2837 : /************************************************************************/
2838 :
2839 0 : const char *GetStateName( int nCode )
2840 : {
2841 : unsigned int i;
2842 :
2843 0 : for( i = 0; i < US_STATE_COUNT; i++ )
2844 : {
2845 0 : if( nCode == aoUSStateTable[i].nCode )
2846 : {
2847 0 : return aoUSStateTable[i].pszName;
2848 : }
2849 : }
2850 0 : return NULL;
2851 : }
2852 :
2853 : /************************************************************************/
2854 : /* GetUnitIndex() */
2855 : /************************************************************************/
2856 :
2857 0 : int GetUnitIndex( const char *pszUnitName )
2858 : {
2859 : int i;
2860 :
2861 0 : for( i = 0; i < (int) LINEAR_UNITS_COUNT; i++ )
2862 : {
2863 0 : if EQUAL( pszUnitName, aoLinearUnitsConv[i].pszName )
2864 : {
2865 0 : return i;
2866 : }
2867 : }
2868 0 : return -1;
2869 : }
2870 :
2871 : /************************************************************************/
2872 : /* GetToMeterIndex() */
2873 : /************************************************************************/
2874 :
2875 0 : int GetToMeterIndex( const char *pszToMeter )
2876 : {
2877 0 : double dfToMeter = atof_nz(pszToMeter);
2878 :
2879 0 : if( dfToMeter != 0.0 )
2880 : {
2881 : int i;
2882 :
2883 0 : for( i = 0; i < (int) LINEAR_UNITS_COUNT; i++ )
2884 : {
2885 0 : if ( ABS( aoLinearUnitsConv[i].dfConv - dfToMeter ) < 0.00001 )
2886 : {
2887 0 : return i;
2888 : }
2889 : }
2890 : }
2891 :
2892 0 : return -1;
2893 : }
2894 :
2895 : /************************************************************************/
2896 : /* GetUnitDefault() */
2897 : /************************************************************************/
2898 :
2899 0 : char *GetUnitDefault( const char *pszUnitName, const char *pszToMeter )
2900 : {
2901 0 : int nIndex = GetUnitIndex( pszUnitName );
2902 :
2903 0 : if( nIndex == -1 && pszToMeter != NULL )
2904 : {
2905 0 : nIndex = GetToMeterIndex( pszToMeter );
2906 : }
2907 :
2908 0 : if( nIndex == -1 )
2909 : {
2910 0 : return CPLStrdup( "Unknown" );
2911 : }
2912 :
2913 0 : return CPLStrdup( aoLinearUnitsConv[aoLinearUnitsConv[nIndex].nDefaultI].pszName );
2914 : }
2915 :
2916 : /************************************************************************/
2917 : /* CSLSaveCRLF() */
2918 : /************************************************************************/
2919 :
2920 : /***
2921 : * Write a stringlist to a CR + LF terminated text file.
2922 : *
2923 : * Returns the number of lines written, or 0 if the file could not
2924 : * be written.
2925 : */
2926 :
2927 16 : int SaveAsCRLF(char **papszStrList, const char *pszFname)
2928 : {
2929 : FILE *fp;
2930 16 : int nLines = 0;
2931 :
2932 16 : if (papszStrList)
2933 : {
2934 16 : if ((fp = VSIFOpenL(pszFname, "wt")) != NULL)
2935 : {
2936 448 : while(*papszStrList != NULL)
2937 : {
2938 416 : if( VSIFPrintfL( fp, "%s\r\n", *papszStrList ) < 1 )
2939 : {
2940 : CPLError( CE_Failure, CPLE_FileIO,
2941 : "CSLSaveCRLF(\"%s\") failed: unable to write to output file.",
2942 0 : pszFname );
2943 0 : break;
2944 : }
2945 :
2946 416 : nLines++;
2947 416 : papszStrList++;
2948 : }
2949 :
2950 16 : VSIFCloseL(fp);
2951 : }
2952 : else
2953 : {
2954 : CPLError( CE_Failure, CPLE_OpenFailed,
2955 : "CSLSaveCRLF(\"%s\") failed: unable to open output file.",
2956 0 : pszFname );
2957 : }
2958 : }
2959 :
2960 16 : return nLines;
2961 : }
2962 :
2963 : /************************************************************************/
2964 : /* GDALRegister_IDRISI() */
2965 : /************************************************************************/
2966 409 : void GDALRegister_IDRISI()
2967 : {
2968 : GDALDriver *poDriver;
2969 :
2970 409 : if( GDALGetDriverByName( "RST" ) == NULL )
2971 : {
2972 392 : poDriver = new GDALDriver();
2973 :
2974 392 : poDriver->SetDescription( "RST" );
2975 392 : poDriver->SetMetadataItem( GDAL_DMD_LONGNAME, rstVERSION );
2976 392 : poDriver->SetMetadataItem( GDAL_DMD_HELPTOPIC, "frmt_Idrisi.html" );
2977 392 : poDriver->SetMetadataItem( GDAL_DMD_EXTENSION, extRST );
2978 392 : poDriver->SetMetadataItem( GDAL_DMD_CREATIONDATATYPES, "Byte Int16 Float32" );
2979 392 : poDriver->pfnOpen = IdrisiDataset::Open;
2980 392 : poDriver->pfnCreate = IdrisiDataset::Create;
2981 392 : poDriver->pfnCreateCopy = IdrisiDataset::CreateCopy;
2982 :
2983 392 : GetGDALDriverManager()->RegisterDriver( poDriver );
2984 : }
2985 409 : }
2986 :
|