1 : /******************************************************************************
2 : * $Id: zmapdataset.cpp 22671 2011-07-08 18:08:12Z rouault $
3 : *
4 : * Project: ZMap driver
5 : * Purpose: GDALDataset driver for ZMap dataset.
6 : * Author: Even Rouault, <even dot rouault at mines dash paris dot org>
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 2011, Even Rouault
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 "cpl_vsi_virtual.h"
31 : #include "cpl_string.h"
32 : #include "gdal_pam.h"
33 :
34 : CPL_CVSID("$Id: zmapdataset.cpp 22671 2011-07-08 18:08:12Z rouault $");
35 :
36 : CPL_C_START
37 : void GDALRegister_ZMap(void);
38 : CPL_C_END
39 :
40 : /************************************************************************/
41 : /* ==================================================================== */
42 : /* ZMapDataset */
43 : /* ==================================================================== */
44 : /************************************************************************/
45 :
46 : class ZMapRasterBand;
47 :
48 : class ZMapDataset : public GDALPamDataset
49 : {
50 : friend class ZMapRasterBand;
51 :
52 : VSILFILE *fp;
53 : int nValuesPerLine;
54 : int nFieldSize;
55 : int nDecimalCount;
56 : int nColNum;
57 : double dfNoDataValue;
58 : vsi_l_offset nDataStartOff;
59 : double adfGeoTransform[6];
60 :
61 : public:
62 : ZMapDataset();
63 : virtual ~ZMapDataset();
64 :
65 : virtual CPLErr GetGeoTransform( double * );
66 :
67 : static GDALDataset *Open( GDALOpenInfo * );
68 : static int Identify( GDALOpenInfo * );
69 : static GDALDataset *CreateCopy( const char * pszFilename, GDALDataset *poSrcDS,
70 : int bStrict, char ** papszOptions,
71 : GDALProgressFunc pfnProgress, void * pProgressData );
72 : };
73 :
74 : /************************************************************************/
75 : /* ==================================================================== */
76 : /* ZMapRasterBand */
77 : /* ==================================================================== */
78 : /************************************************************************/
79 :
80 : class ZMapRasterBand : public GDALPamRasterBand
81 14 : {
82 : friend class ZMapDataset;
83 :
84 : public:
85 :
86 : ZMapRasterBand( ZMapDataset * );
87 :
88 : virtual CPLErr IReadBlock( int, int, void * );
89 :
90 : virtual double GetNoDataValue( int *pbSuccess = NULL );
91 : };
92 :
93 :
94 : /************************************************************************/
95 : /* ZMapRasterBand() */
96 : /************************************************************************/
97 :
98 14 : ZMapRasterBand::ZMapRasterBand( ZMapDataset *poDS )
99 :
100 : {
101 14 : this->poDS = poDS;
102 14 : this->nBand = nBand;
103 :
104 14 : eDataType = GDT_Float64;
105 :
106 14 : nBlockXSize = 1;
107 14 : nBlockYSize = poDS->GetRasterYSize();
108 14 : }
109 :
110 : /************************************************************************/
111 : /* IReadBlock() */
112 : /************************************************************************/
113 :
114 40 : CPLErr ZMapRasterBand::IReadBlock( int nBlockXOff, int nBlockYOff,
115 : void * pImage )
116 :
117 : {
118 : int i;
119 40 : ZMapDataset *poGDS = (ZMapDataset *) poDS;
120 :
121 40 : if (poGDS->fp == NULL)
122 0 : return CE_Failure;
123 :
124 40 : if (nBlockXOff < poGDS->nColNum + 1)
125 : {
126 0 : VSIFSeekL(poGDS->fp, poGDS->nDataStartOff, SEEK_SET);
127 0 : poGDS->nColNum = -1;
128 : }
129 :
130 40 : if (nBlockXOff > poGDS->nColNum + 1)
131 : {
132 0 : for(i=poGDS->nColNum + 1;i<nBlockXOff;i++)
133 : {
134 0 : if (IReadBlock(i,0,pImage) != CE_None)
135 0 : return CE_Failure;
136 : }
137 : }
138 :
139 : char* pszLine;
140 40 : i = 0;
141 40 : double dfExp = pow(10.0, poGDS->nDecimalCount);
142 280 : while(i<nRasterYSize)
143 : {
144 200 : pszLine = (char*)CPLReadLineL(poGDS->fp);
145 200 : if (pszLine == NULL)
146 0 : return CE_Failure;
147 200 : int nExpected = nRasterYSize - i;
148 200 : if (nExpected > poGDS->nValuesPerLine)
149 160 : nExpected = poGDS->nValuesPerLine;
150 200 : if ((int)strlen(pszLine) != nExpected * poGDS->nFieldSize)
151 0 : return CE_Failure;
152 :
153 1000 : for(int j=0;j<nExpected;j++)
154 : {
155 800 : char* pszValue = pszLine + j * poGDS->nFieldSize;
156 800 : char chSaved = pszValue[poGDS->nFieldSize];
157 800 : pszValue[poGDS->nFieldSize] = 0;
158 800 : if (strchr(pszValue, '.') != NULL)
159 800 : ((double*)pImage)[i+j] = CPLAtofM(pszValue);
160 : else
161 0 : ((double*)pImage)[i+j] = atoi(pszValue) * dfExp;
162 800 : pszValue[poGDS->nFieldSize] = chSaved;
163 : }
164 :
165 200 : i += nExpected;
166 : }
167 :
168 40 : poGDS->nColNum ++;
169 :
170 40 : return CE_None;
171 : }
172 :
173 : /************************************************************************/
174 : /* GetNoDataValue() */
175 : /************************************************************************/
176 :
177 2 : double ZMapRasterBand::GetNoDataValue( int *pbSuccess )
178 : {
179 2 : ZMapDataset *poGDS = (ZMapDataset *) poDS;
180 :
181 2 : if (pbSuccess)
182 2 : *pbSuccess = TRUE;
183 :
184 2 : return poGDS->dfNoDataValue;
185 : }
186 :
187 : /************************************************************************/
188 : /* ~ZMapDataset() */
189 : /************************************************************************/
190 :
191 14 : ZMapDataset::ZMapDataset()
192 : {
193 14 : fp = NULL;
194 14 : nDataStartOff = 0;
195 14 : nColNum = -1;
196 14 : nValuesPerLine = 0;
197 14 : nFieldSize = 0;
198 14 : nDecimalCount = 0;
199 14 : dfNoDataValue = 0.0;
200 14 : adfGeoTransform[0] = 0;
201 14 : adfGeoTransform[1] = 1;
202 14 : adfGeoTransform[2] = 0;
203 14 : adfGeoTransform[3] = 0;
204 14 : adfGeoTransform[4] = 0;
205 14 : adfGeoTransform[5] = 1;
206 14 : }
207 :
208 : /************************************************************************/
209 : /* ~ZMapDataset() */
210 : /************************************************************************/
211 :
212 14 : ZMapDataset::~ZMapDataset()
213 :
214 : {
215 14 : FlushCache();
216 14 : if (fp)
217 14 : VSIFCloseL(fp);
218 14 : }
219 :
220 : /************************************************************************/
221 : /* Identify() */
222 : /************************************************************************/
223 :
224 10491 : int ZMapDataset::Identify( GDALOpenInfo * poOpenInfo )
225 : {
226 : int i;
227 :
228 10491 : if (poOpenInfo->nHeaderBytes == 0)
229 10328 : return FALSE;
230 :
231 : /* -------------------------------------------------------------------- */
232 : /* Chech that it looks roughly as a ZMap dataset */
233 : /* -------------------------------------------------------------------- */
234 163 : const char* pszData = (const char*)poOpenInfo->pabyHeader;
235 :
236 : /* Skip comments line at the beginning */
237 163 : i=0;
238 163 : if (pszData[i] == '!')
239 : {
240 14 : i++;
241 280 : for(;i<poOpenInfo->nHeaderBytes;i++)
242 : {
243 280 : char ch = pszData[i];
244 280 : if (ch == 13 || ch == 10)
245 : {
246 42 : i++;
247 42 : if (ch == 13 && pszData[i] == 10)
248 0 : i++;
249 42 : if (pszData[i] != '!')
250 14 : break;
251 : }
252 : }
253 : }
254 :
255 163 : if (pszData[i] != '@')
256 148 : return FALSE;
257 15 : i++;
258 :
259 15 : char** papszTokens = CSLTokenizeString2( pszData+i, ",", 0 );
260 15 : if (CSLCount(papszTokens) < 3)
261 : {
262 1 : CSLDestroy(papszTokens);
263 1 : return FALSE;
264 : }
265 :
266 14 : const char* pszToken = papszTokens[1];
267 42 : while (*pszToken == ' ')
268 14 : pszToken ++;
269 :
270 14 : if (strncmp(pszToken, "GRID", 4) != 0)
271 : {
272 0 : CSLDestroy(papszTokens);
273 0 : return FALSE;
274 : }
275 :
276 14 : CSLDestroy(papszTokens);
277 14 : return TRUE;
278 : }
279 :
280 : /************************************************************************/
281 : /* Open() */
282 : /************************************************************************/
283 :
284 1253 : GDALDataset *ZMapDataset::Open( GDALOpenInfo * poOpenInfo )
285 :
286 : {
287 1253 : if (!Identify(poOpenInfo))
288 1239 : return NULL;
289 :
290 : /* -------------------------------------------------------------------- */
291 : /* Find dataset characteristics */
292 : /* -------------------------------------------------------------------- */
293 14 : VSILFILE* fp = VSIFOpenL(poOpenInfo->pszFilename, "rb");
294 14 : if (fp == NULL)
295 0 : return NULL;
296 :
297 : const char* pszLine;
298 :
299 70 : while((pszLine = CPLReadLine2L(fp, 100, NULL)) != NULL)
300 : {
301 56 : if (*pszLine == '!')
302 : {
303 42 : continue;
304 : }
305 : else
306 14 : break;
307 : }
308 14 : if (pszLine == NULL)
309 : {
310 0 : VSIFCloseL(fp);
311 0 : return NULL;
312 : }
313 :
314 : /* Parse first header line */
315 14 : char** papszTokens = CSLTokenizeString2( pszLine, ",", 0 );
316 14 : if (CSLCount(papszTokens) != 3)
317 : {
318 0 : CSLDestroy(papszTokens);
319 0 : VSIFCloseL(fp);
320 0 : return NULL;
321 : }
322 :
323 14 : int nValuesPerLine = atoi(papszTokens[2]);
324 14 : if (nValuesPerLine <= 0)
325 : {
326 0 : CSLDestroy(papszTokens);
327 0 : VSIFCloseL(fp);
328 0 : return NULL;
329 : }
330 :
331 14 : CSLDestroy(papszTokens);
332 14 : papszTokens = NULL;
333 :
334 : /* Parse second header line */
335 14 : pszLine = CPLReadLine2L(fp, 100, NULL);
336 14 : if (pszLine == NULL)
337 : {
338 0 : VSIFCloseL(fp);
339 0 : return NULL;
340 : }
341 14 : papszTokens = CSLTokenizeString2( pszLine, ",", 0 );
342 14 : if (CSLCount(papszTokens) != 5)
343 : {
344 0 : CSLDestroy(papszTokens);
345 0 : VSIFCloseL(fp);
346 0 : return NULL;
347 : }
348 :
349 14 : int nFieldSize = atoi(papszTokens[0]);
350 14 : double dfNoDataValue = CPLAtofM(papszTokens[1]);
351 14 : int nDecimalCount = atoi(papszTokens[3]);
352 14 : int nColumnNumber = atoi(papszTokens[4]);
353 :
354 14 : CSLDestroy(papszTokens);
355 14 : papszTokens = NULL;
356 :
357 14 : if (nFieldSize <= 0 || nFieldSize >= 40 ||
358 : nDecimalCount <= 0 || nDecimalCount >= nFieldSize ||
359 : nColumnNumber != 1)
360 : {
361 : CPLDebug("ZMap", "nFieldSize=%d, nDecimalCount=%d, nColumnNumber=%d",
362 0 : nFieldSize, nDecimalCount, nColumnNumber);
363 0 : VSIFCloseL(fp);
364 0 : return NULL;
365 : }
366 :
367 : /* Parse third header line */
368 14 : pszLine = CPLReadLine2L(fp, 100, NULL);
369 14 : if (pszLine == NULL)
370 : {
371 0 : VSIFCloseL(fp);
372 0 : return NULL;
373 : }
374 14 : papszTokens = CSLTokenizeString2( pszLine, ",", 0 );
375 14 : if (CSLCount(papszTokens) != 6)
376 : {
377 0 : CSLDestroy(papszTokens);
378 0 : VSIFCloseL(fp);
379 0 : return NULL;
380 : }
381 :
382 14 : int nRows = atoi(papszTokens[0]);
383 14 : int nCols = atoi(papszTokens[1]);
384 14 : double dfMinX = CPLAtofM(papszTokens[2]);
385 14 : double dfMaxX = CPLAtofM(papszTokens[3]);
386 14 : double dfMinY = CPLAtofM(papszTokens[4]);
387 14 : double dfMaxY = CPLAtofM(papszTokens[5]);
388 :
389 14 : CSLDestroy(papszTokens);
390 14 : papszTokens = NULL;
391 :
392 14 : if (!GDALCheckDatasetDimensions(nCols, nRows) ||
393 : nCols == 1 || nRows == 1)
394 : {
395 0 : VSIFCloseL(fp);
396 0 : return NULL;
397 : }
398 :
399 : /* Ignore fourth header line */
400 14 : pszLine = CPLReadLine2L(fp, 100, NULL);
401 14 : if (pszLine == NULL)
402 : {
403 0 : VSIFCloseL(fp);
404 0 : return NULL;
405 : }
406 :
407 : /* Check fifth header line */
408 14 : pszLine = CPLReadLine2L(fp, 100, NULL);
409 14 : if (pszLine == NULL || pszLine[0] != '@')
410 : {
411 0 : VSIFCloseL(fp);
412 0 : return NULL;
413 : }
414 :
415 : /* -------------------------------------------------------------------- */
416 : /* Create a corresponding GDALDataset. */
417 : /* -------------------------------------------------------------------- */
418 : ZMapDataset *poDS;
419 :
420 14 : poDS = new ZMapDataset();
421 14 : poDS->fp = fp;
422 14 : poDS->nDataStartOff = VSIFTellL(fp);
423 14 : poDS->nValuesPerLine = nValuesPerLine;
424 14 : poDS->nFieldSize = nFieldSize;
425 14 : poDS->nDecimalCount = nDecimalCount;
426 14 : poDS->nRasterXSize = nCols;
427 14 : poDS->nRasterYSize = nRows;
428 14 : poDS->dfNoDataValue = dfNoDataValue;
429 :
430 14 : if (CSLTestBoolean(CPLGetConfigOption("ZMAP_PIXEL_IS_POINT", "FALSE")))
431 : {
432 0 : double dfStepX = (dfMaxX - dfMinX) / (nCols - 1);
433 0 : double dfStepY = (dfMaxY - dfMinY) / (nRows - 1);
434 :
435 0 : poDS->adfGeoTransform[0] = dfMinX - dfStepX / 2;
436 0 : poDS->adfGeoTransform[1] = dfStepX;
437 0 : poDS->adfGeoTransform[3] = dfMaxY + dfStepY / 2;
438 0 : poDS->adfGeoTransform[5] = -dfStepY;
439 : }
440 : else
441 : {
442 14 : double dfStepX = (dfMaxX - dfMinX) / nCols ;
443 14 : double dfStepY = (dfMaxY - dfMinY) / nRows;
444 :
445 14 : poDS->adfGeoTransform[0] = dfMinX;
446 14 : poDS->adfGeoTransform[1] = dfStepX;
447 14 : poDS->adfGeoTransform[3] = dfMaxY;
448 14 : poDS->adfGeoTransform[5] = -dfStepY;
449 : }
450 :
451 : /* -------------------------------------------------------------------- */
452 : /* Create band information objects. */
453 : /* -------------------------------------------------------------------- */
454 14 : poDS->nBands = 1;
455 14 : poDS->SetBand( 1, new ZMapRasterBand( poDS ) );
456 :
457 : /* -------------------------------------------------------------------- */
458 : /* Initialize any PAM information. */
459 : /* -------------------------------------------------------------------- */
460 14 : poDS->SetDescription( poOpenInfo->pszFilename );
461 14 : poDS->TryLoadXML();
462 :
463 : /* -------------------------------------------------------------------- */
464 : /* Support overviews. */
465 : /* -------------------------------------------------------------------- */
466 14 : poDS->oOvManager.Initialize( poDS, poOpenInfo->pszFilename );
467 14 : return( poDS );
468 : }
469 :
470 :
471 : /************************************************************************/
472 : /* WriteRightJustified() */
473 : /************************************************************************/
474 :
475 1632 : static void WriteRightJustified(VSILFILE* fp, const char *pszValue, int nWidth)
476 : {
477 1632 : int nLen = strlen(pszValue);
478 1632 : CPLAssert(nLen <= nWidth);
479 : int i;
480 16072 : for(i=0;i<nWidth -nLen;i++)
481 14440 : VSIFWriteL(" ", 1, 1, fp);
482 1632 : VSIFWriteL(pszValue, 1, nLen, fp);
483 1632 : }
484 :
485 60 : static void WriteRightJustified(VSILFILE* fp, int nValue, int nWidth)
486 : {
487 60 : CPLString osValue(CPLSPrintf("%d", nValue));
488 60 : WriteRightJustified(fp, osValue.c_str(), nWidth);
489 60 : }
490 :
491 1560 : static void WriteRightJustified(VSILFILE* fp, double dfValue, int nWidth,
492 : int nDecimals = -1)
493 : {
494 : char szFormat[32];
495 1560 : if (nDecimals >= 0)
496 1548 : sprintf(szFormat, "%%.%df", nDecimals);
497 : else
498 12 : sprintf(szFormat, "%%g");
499 1560 : char* pszValue = (char*)CPLSPrintf(szFormat, dfValue);
500 1560 : char* pszE = strchr(pszValue, 'e');
501 1560 : if (pszE)
502 12 : *pszE = 'E';
503 :
504 1560 : if ((int)strlen(pszValue) > nWidth)
505 : {
506 2 : sprintf(szFormat, "%%.%dg", nDecimals);
507 2 : pszValue = (char*)CPLSPrintf(szFormat, dfValue);
508 2 : pszE = strchr(pszValue, 'e');
509 2 : if (pszE)
510 0 : *pszE = 'E';
511 : }
512 :
513 1560 : CPLString osValue(pszValue);
514 1560 : WriteRightJustified(fp, osValue.c_str(), nWidth);
515 1560 : }
516 :
517 : /************************************************************************/
518 : /* CreateCopy() */
519 : /************************************************************************/
520 :
521 19 : GDALDataset* ZMapDataset::CreateCopy( const char * pszFilename,
522 : GDALDataset *poSrcDS,
523 : int bStrict, char ** papszOptions,
524 : GDALProgressFunc pfnProgress,
525 : void * pProgressData )
526 : {
527 : /* -------------------------------------------------------------------- */
528 : /* Some some rudimentary checks */
529 : /* -------------------------------------------------------------------- */
530 19 : int nBands = poSrcDS->GetRasterCount();
531 19 : if (nBands == 0)
532 : {
533 : CPLError( CE_Failure, CPLE_NotSupported,
534 1 : "ZMap driver does not support source dataset with zero band.\n");
535 1 : return NULL;
536 : }
537 :
538 18 : if (nBands != 1)
539 : {
540 : CPLError( (bStrict) ? CE_Failure : CE_Warning, CPLE_NotSupported,
541 4 : "ZMap driver only uses the first band of the dataset.\n");
542 4 : if (bStrict)
543 4 : return NULL;
544 : }
545 :
546 14 : if( pfnProgress && !pfnProgress( 0.0, NULL, pProgressData ) )
547 0 : return NULL;
548 :
549 : /* -------------------------------------------------------------------- */
550 : /* Get source dataset info */
551 : /* -------------------------------------------------------------------- */
552 :
553 14 : int nXSize = poSrcDS->GetRasterXSize();
554 14 : int nYSize = poSrcDS->GetRasterYSize();
555 14 : if (nXSize == 1 || nYSize == 1)
556 : {
557 0 : return NULL;
558 : }
559 :
560 : double adfGeoTransform[6];
561 14 : poSrcDS->GetGeoTransform(adfGeoTransform);
562 14 : if (adfGeoTransform[2] != 0 || adfGeoTransform[4] != 0)
563 : {
564 : CPLError( CE_Failure, CPLE_NotSupported,
565 0 : "ZMap driver does not support CreateCopy() from skewed or rotated dataset.\n");
566 0 : return NULL;
567 : }
568 :
569 : /* -------------------------------------------------------------------- */
570 : /* Create target file */
571 : /* -------------------------------------------------------------------- */
572 :
573 14 : VSILFILE* fp = VSIFOpenL(pszFilename, "wb");
574 14 : if (fp == NULL)
575 : {
576 : CPLError( CE_Failure, CPLE_AppDefined,
577 2 : "Cannot create %s", pszFilename );
578 2 : return NULL;
579 : }
580 :
581 12 : int nFieldSize = 20;
582 12 : int nValuesPerLine = 4;
583 12 : int nDecimalCount = 7;
584 :
585 12 : int bHasNoDataValue = FALSE;
586 : double dfNoDataValue =
587 12 : poSrcDS->GetRasterBand(1)->GetNoDataValue(&bHasNoDataValue);
588 12 : if (!bHasNoDataValue)
589 12 : dfNoDataValue = 1.e30;
590 :
591 12 : VSIFPrintfL(fp, "!\n");
592 12 : VSIFPrintfL(fp, "! Created by GDAL.\n");
593 12 : VSIFPrintfL(fp, "!\n");
594 12 : VSIFPrintfL(fp, "@GRID FILE, GRID, %d\n", nValuesPerLine);
595 :
596 12 : WriteRightJustified(fp, nFieldSize, 10);
597 12 : VSIFPrintfL(fp, ",");
598 12 : WriteRightJustified(fp, dfNoDataValue, 10);
599 12 : VSIFPrintfL(fp, ",");
600 12 : WriteRightJustified(fp, "", 10);
601 12 : VSIFPrintfL(fp, ",");
602 12 : WriteRightJustified(fp, nDecimalCount, 10);
603 12 : VSIFPrintfL(fp, ",");
604 12 : WriteRightJustified(fp, 1, 10);
605 12 : VSIFPrintfL(fp, "\n");
606 :
607 12 : WriteRightJustified(fp, nYSize, 10);
608 12 : VSIFPrintfL(fp, ",");
609 12 : WriteRightJustified(fp, nXSize, 10);
610 12 : VSIFPrintfL(fp, ",");
611 :
612 12 : if (CSLTestBoolean(CPLGetConfigOption("ZMAP_PIXEL_IS_POINT", "FALSE")))
613 : {
614 0 : WriteRightJustified(fp, adfGeoTransform[0] + adfGeoTransform[1] / 2, 14, 7);
615 0 : VSIFPrintfL(fp, ",");
616 0 : WriteRightJustified(fp, adfGeoTransform[0] + adfGeoTransform[1] * nXSize -
617 0 : adfGeoTransform[1] / 2, 14, 7);
618 0 : VSIFPrintfL(fp, ",");
619 0 : WriteRightJustified(fp, adfGeoTransform[3] + adfGeoTransform[5] * nYSize -
620 0 : adfGeoTransform[5] / 2, 14, 7);
621 0 : VSIFPrintfL(fp, ",");
622 0 : WriteRightJustified(fp, adfGeoTransform[3] + adfGeoTransform[5] / 2, 14, 7);
623 : }
624 : else
625 : {
626 12 : WriteRightJustified(fp, adfGeoTransform[0], 14, 7);
627 12 : VSIFPrintfL(fp, ",");
628 12 : WriteRightJustified(fp, adfGeoTransform[0] + adfGeoTransform[1] * nXSize, 14, 7);
629 12 : VSIFPrintfL(fp, ",");
630 12 : WriteRightJustified(fp, adfGeoTransform[3] + adfGeoTransform[5] * nYSize, 14, 7);
631 12 : VSIFPrintfL(fp, ",");
632 12 : WriteRightJustified(fp, adfGeoTransform[3], 14, 7);
633 : }
634 :
635 12 : VSIFPrintfL(fp, "\n");
636 :
637 12 : VSIFPrintfL(fp, "0.0, 0.0, 0.0\n");
638 12 : VSIFPrintfL(fp, "@\n");
639 :
640 : /* -------------------------------------------------------------------- */
641 : /* Copy imagery */
642 : /* -------------------------------------------------------------------- */
643 12 : double* padfLineBuffer = (double*) CPLMalloc(nYSize * sizeof(double));
644 : int i, j;
645 12 : CPLErr eErr = CE_None;
646 142 : for(i=0;i<nXSize && eErr == CE_None;i++)
647 : {
648 : eErr = poSrcDS->GetRasterBand(1)->RasterIO(
649 : GF_Read, i, 0, 1, nYSize,
650 : padfLineBuffer, 1, nYSize,
651 130 : GDT_Float64, 0, 0);
652 130 : if (eErr != CE_None)
653 0 : break;
654 130 : int bEOLPrinted = FALSE;
655 1630 : for(j=0;j<nYSize;j++)
656 : {
657 1500 : WriteRightJustified(fp, padfLineBuffer[j], nFieldSize, nDecimalCount);
658 1500 : if (((j + 1) % nValuesPerLine) == 0)
659 : {
660 320 : bEOLPrinted = TRUE;
661 320 : VSIFPrintfL(fp, "\n");
662 : }
663 : else
664 1180 : bEOLPrinted = FALSE;
665 : }
666 130 : if (!bEOLPrinted)
667 110 : VSIFPrintfL(fp, "\n");
668 :
669 130 : if (!pfnProgress( (j+1) * 1.0 / nYSize, NULL, pProgressData))
670 0 : break;
671 : }
672 12 : CPLFree(padfLineBuffer);
673 12 : VSIFCloseL(fp);
674 :
675 12 : if (eErr != CE_None)
676 0 : return NULL;
677 :
678 12 : return (GDALDataset*) GDALOpen(pszFilename, GA_ReadOnly);
679 : }
680 :
681 : /************************************************************************/
682 : /* GetGeoTransform() */
683 : /************************************************************************/
684 :
685 1 : CPLErr ZMapDataset::GetGeoTransform( double * padfTransform )
686 :
687 : {
688 1 : memcpy(padfTransform, adfGeoTransform, 6 * sizeof(double));
689 :
690 1 : return( CE_None );
691 : }
692 :
693 : /************************************************************************/
694 : /* GDALRegister_ZMap() */
695 : /************************************************************************/
696 :
697 558 : void GDALRegister_ZMap()
698 :
699 : {
700 : GDALDriver *poDriver;
701 :
702 558 : if( GDALGetDriverByName( "ZMap" ) == NULL )
703 : {
704 537 : poDriver = new GDALDriver();
705 :
706 537 : poDriver->SetDescription( "ZMap" );
707 : poDriver->SetMetadataItem( GDAL_DMD_LONGNAME,
708 537 : "ZMap Plus Grid" );
709 : poDriver->SetMetadataItem( GDAL_DMD_HELPTOPIC,
710 537 : "frmt_various.html#ZMap" );
711 537 : poDriver->SetMetadataItem( GDAL_DMD_EXTENSION, "dat" );
712 :
713 537 : poDriver->SetMetadataItem( GDAL_DCAP_VIRTUALIO, "YES" );
714 :
715 537 : poDriver->pfnOpen = ZMapDataset::Open;
716 537 : poDriver->pfnIdentify = ZMapDataset::Identify;
717 537 : poDriver->pfnCreateCopy = ZMapDataset::CreateCopy;
718 :
719 537 : GetGDALDriverManager()->RegisterDriver( poDriver );
720 : }
721 558 : }
722 :
|
