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 28 : {
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 28 : ZMapRasterBand::ZMapRasterBand( ZMapDataset *poDS )
99 :
100 : {
101 28 : this->poDS = poDS;
102 28 : this->nBand = nBand;
103 :
104 28 : eDataType = GDT_Float64;
105 :
106 28 : nBlockXSize = 1;
107 28 : nBlockYSize = poDS->GetRasterYSize();
108 28 : }
109 :
110 : /************************************************************************/
111 : /* IReadBlock() */
112 : /************************************************************************/
113 :
114 80 : CPLErr ZMapRasterBand::IReadBlock( int nBlockXOff, int nBlockYOff,
115 : void * pImage )
116 :
117 : {
118 : int i;
119 80 : ZMapDataset *poGDS = (ZMapDataset *) poDS;
120 :
121 80 : if (poGDS->fp == NULL)
122 0 : return CE_Failure;
123 :
124 80 : if (nBlockXOff < poGDS->nColNum + 1)
125 : {
126 0 : VSIFSeekL(poGDS->fp, poGDS->nDataStartOff, SEEK_SET);
127 0 : poGDS->nColNum = -1;
128 : }
129 :
130 80 : 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 80 : i = 0;
141 80 : double dfExp = pow(10.0, poGDS->nDecimalCount);
142 560 : while(i<nRasterYSize)
143 : {
144 400 : pszLine = (char*)CPLReadLineL(poGDS->fp);
145 400 : if (pszLine == NULL)
146 0 : return CE_Failure;
147 400 : int nExpected = nRasterYSize - i;
148 400 : if (nExpected > poGDS->nValuesPerLine)
149 320 : nExpected = poGDS->nValuesPerLine;
150 400 : if ((int)strlen(pszLine) != nExpected * poGDS->nFieldSize)
151 0 : return CE_Failure;
152 :
153 2000 : for(int j=0;j<nExpected;j++)
154 : {
155 1600 : char* pszValue = pszLine + j * poGDS->nFieldSize;
156 1600 : char chSaved = pszValue[poGDS->nFieldSize];
157 1600 : pszValue[poGDS->nFieldSize] = 0;
158 1600 : if (strchr(pszValue, '.') != NULL)
159 1600 : ((double*)pImage)[i+j] = CPLAtofM(pszValue);
160 : else
161 0 : ((double*)pImage)[i+j] = atoi(pszValue) * dfExp;
162 1600 : pszValue[poGDS->nFieldSize] = chSaved;
163 : }
164 :
165 400 : i += nExpected;
166 : }
167 :
168 80 : poGDS->nColNum ++;
169 :
170 80 : return CE_None;
171 : }
172 :
173 : /************************************************************************/
174 : /* GetNoDataValue() */
175 : /************************************************************************/
176 :
177 4 : double ZMapRasterBand::GetNoDataValue( int *pbSuccess )
178 : {
179 4 : ZMapDataset *poGDS = (ZMapDataset *) poDS;
180 :
181 4 : if (pbSuccess)
182 4 : *pbSuccess = TRUE;
183 :
184 4 : return poGDS->dfNoDataValue;
185 : }
186 :
187 : /************************************************************************/
188 : /* ~ZMapDataset() */
189 : /************************************************************************/
190 :
191 28 : ZMapDataset::ZMapDataset()
192 : {
193 28 : fp = NULL;
194 28 : nDataStartOff = 0;
195 28 : nColNum = -1;
196 28 : nValuesPerLine = 0;
197 28 : nFieldSize = 0;
198 28 : nDecimalCount = 0;
199 28 : dfNoDataValue = 0.0;
200 28 : adfGeoTransform[0] = 0;
201 28 : adfGeoTransform[1] = 1;
202 28 : adfGeoTransform[2] = 0;
203 28 : adfGeoTransform[3] = 0;
204 28 : adfGeoTransform[4] = 0;
205 28 : adfGeoTransform[5] = 1;
206 28 : }
207 :
208 : /************************************************************************/
209 : /* ~ZMapDataset() */
210 : /************************************************************************/
211 :
212 28 : ZMapDataset::~ZMapDataset()
213 :
214 : {
215 28 : FlushCache();
216 28 : if (fp)
217 28 : VSIFCloseL(fp);
218 28 : }
219 :
220 : /************************************************************************/
221 : /* Identify() */
222 : /************************************************************************/
223 :
224 21484 : int ZMapDataset::Identify( GDALOpenInfo * poOpenInfo )
225 : {
226 : int i;
227 :
228 21484 : if (poOpenInfo->nHeaderBytes == 0)
229 21138 : return FALSE;
230 :
231 : /* -------------------------------------------------------------------- */
232 : /* Chech that it looks roughly as a ZMap dataset */
233 : /* -------------------------------------------------------------------- */
234 346 : const char* pszData = (const char*)poOpenInfo->pabyHeader;
235 :
236 : /* Skip comments line at the beginning */
237 346 : i=0;
238 346 : if (pszData[i] == '!')
239 : {
240 28 : i++;
241 560 : for(;i<poOpenInfo->nHeaderBytes;i++)
242 : {
243 560 : char ch = pszData[i];
244 560 : if (ch == 13 || ch == 10)
245 : {
246 84 : i++;
247 84 : if (ch == 13 && pszData[i] == 10)
248 0 : i++;
249 84 : if (pszData[i] != '!')
250 28 : break;
251 : }
252 : }
253 : }
254 :
255 346 : if (pszData[i] != '@')
256 316 : return FALSE;
257 30 : i++;
258 :
259 30 : char** papszTokens = CSLTokenizeString2( pszData+i, ",", 0 );
260 30 : if (CSLCount(papszTokens) < 3)
261 : {
262 2 : CSLDestroy(papszTokens);
263 2 : return FALSE;
264 : }
265 :
266 28 : const char* pszToken = papszTokens[1];
267 84 : while (*pszToken == ' ')
268 28 : pszToken ++;
269 :
270 28 : if (strncmp(pszToken, "GRID", 4) != 0)
271 : {
272 0 : CSLDestroy(papszTokens);
273 0 : return FALSE;
274 : }
275 :
276 28 : CSLDestroy(papszTokens);
277 28 : return TRUE;
278 : }
279 :
280 : /************************************************************************/
281 : /* Open() */
282 : /************************************************************************/
283 :
284 2590 : GDALDataset *ZMapDataset::Open( GDALOpenInfo * poOpenInfo )
285 :
286 : {
287 2590 : if (!Identify(poOpenInfo))
288 2562 : return NULL;
289 :
290 : /* -------------------------------------------------------------------- */
291 : /* Find dataset characteristics */
292 : /* -------------------------------------------------------------------- */
293 28 : VSILFILE* fp = VSIFOpenL(poOpenInfo->pszFilename, "rb");
294 28 : if (fp == NULL)
295 0 : return NULL;
296 :
297 : const char* pszLine;
298 :
299 140 : while((pszLine = CPLReadLine2L(fp, 100, NULL)) != NULL)
300 : {
301 112 : if (*pszLine == '!')
302 : {
303 84 : continue;
304 : }
305 : else
306 28 : break;
307 : }
308 28 : if (pszLine == NULL)
309 : {
310 0 : VSIFCloseL(fp);
311 0 : return NULL;
312 : }
313 :
314 : /* Parse first header line */
315 28 : char** papszTokens = CSLTokenizeString2( pszLine, ",", 0 );
316 28 : if (CSLCount(papszTokens) != 3)
317 : {
318 0 : CSLDestroy(papszTokens);
319 0 : VSIFCloseL(fp);
320 0 : return NULL;
321 : }
322 :
323 28 : int nValuesPerLine = atoi(papszTokens[2]);
324 28 : if (nValuesPerLine <= 0)
325 : {
326 0 : CSLDestroy(papszTokens);
327 0 : VSIFCloseL(fp);
328 0 : return NULL;
329 : }
330 :
331 28 : CSLDestroy(papszTokens);
332 28 : papszTokens = NULL;
333 :
334 : /* Parse second header line */
335 28 : pszLine = CPLReadLine2L(fp, 100, NULL);
336 28 : if (pszLine == NULL)
337 : {
338 0 : VSIFCloseL(fp);
339 0 : return NULL;
340 : }
341 28 : papszTokens = CSLTokenizeString2( pszLine, ",", 0 );
342 28 : if (CSLCount(papszTokens) != 5)
343 : {
344 0 : CSLDestroy(papszTokens);
345 0 : VSIFCloseL(fp);
346 0 : return NULL;
347 : }
348 :
349 28 : int nFieldSize = atoi(papszTokens[0]);
350 28 : double dfNoDataValue = CPLAtofM(papszTokens[1]);
351 28 : int nDecimalCount = atoi(papszTokens[3]);
352 28 : int nColumnNumber = atoi(papszTokens[4]);
353 :
354 28 : CSLDestroy(papszTokens);
355 28 : papszTokens = NULL;
356 :
357 28 : 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 28 : pszLine = CPLReadLine2L(fp, 100, NULL);
369 28 : if (pszLine == NULL)
370 : {
371 0 : VSIFCloseL(fp);
372 0 : return NULL;
373 : }
374 28 : papszTokens = CSLTokenizeString2( pszLine, ",", 0 );
375 28 : if (CSLCount(papszTokens) != 6)
376 : {
377 0 : CSLDestroy(papszTokens);
378 0 : VSIFCloseL(fp);
379 0 : return NULL;
380 : }
381 :
382 28 : int nRows = atoi(papszTokens[0]);
383 28 : int nCols = atoi(papszTokens[1]);
384 28 : double dfMinX = CPLAtofM(papszTokens[2]);
385 28 : double dfMaxX = CPLAtofM(papszTokens[3]);
386 28 : double dfMinY = CPLAtofM(papszTokens[4]);
387 28 : double dfMaxY = CPLAtofM(papszTokens[5]);
388 :
389 28 : CSLDestroy(papszTokens);
390 28 : papszTokens = NULL;
391 :
392 28 : 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 28 : pszLine = CPLReadLine2L(fp, 100, NULL);
401 28 : if (pszLine == NULL)
402 : {
403 0 : VSIFCloseL(fp);
404 0 : return NULL;
405 : }
406 :
407 : /* Check fifth header line */
408 28 : pszLine = CPLReadLine2L(fp, 100, NULL);
409 28 : 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 28 : poDS = new ZMapDataset();
421 28 : poDS->fp = fp;
422 28 : poDS->nDataStartOff = VSIFTellL(fp);
423 28 : poDS->nValuesPerLine = nValuesPerLine;
424 28 : poDS->nFieldSize = nFieldSize;
425 28 : poDS->nDecimalCount = nDecimalCount;
426 28 : poDS->nRasterXSize = nCols;
427 28 : poDS->nRasterYSize = nRows;
428 28 : poDS->dfNoDataValue = dfNoDataValue;
429 :
430 28 : 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 28 : double dfStepX = (dfMaxX - dfMinX) / nCols ;
443 28 : double dfStepY = (dfMaxY - dfMinY) / nRows;
444 :
445 28 : poDS->adfGeoTransform[0] = dfMinX;
446 28 : poDS->adfGeoTransform[1] = dfStepX;
447 28 : poDS->adfGeoTransform[3] = dfMaxY;
448 28 : poDS->adfGeoTransform[5] = -dfStepY;
449 : }
450 :
451 : /* -------------------------------------------------------------------- */
452 : /* Create band information objects. */
453 : /* -------------------------------------------------------------------- */
454 28 : poDS->nBands = 1;
455 28 : poDS->SetBand( 1, new ZMapRasterBand( poDS ) );
456 :
457 : /* -------------------------------------------------------------------- */
458 : /* Initialize any PAM information. */
459 : /* -------------------------------------------------------------------- */
460 28 : poDS->SetDescription( poOpenInfo->pszFilename );
461 28 : poDS->TryLoadXML();
462 :
463 : /* -------------------------------------------------------------------- */
464 : /* Support overviews. */
465 : /* -------------------------------------------------------------------- */
466 28 : poDS->oOvManager.Initialize( poDS, poOpenInfo->pszFilename );
467 28 : return( poDS );
468 : }
469 :
470 :
471 : /************************************************************************/
472 : /* WriteRightJustified() */
473 : /************************************************************************/
474 :
475 3264 : static void WriteRightJustified(VSILFILE* fp, const char *pszValue, int nWidth)
476 : {
477 3264 : int nLen = strlen(pszValue);
478 3264 : CPLAssert(nLen <= nWidth);
479 : int i;
480 32144 : for(i=0;i<nWidth -nLen;i++)
481 28880 : VSIFWriteL(" ", 1, 1, fp);
482 3264 : VSIFWriteL(pszValue, 1, nLen, fp);
483 3264 : }
484 :
485 120 : static void WriteRightJustified(VSILFILE* fp, int nValue, int nWidth)
486 : {
487 120 : CPLString osValue(CPLSPrintf("%d", nValue));
488 120 : WriteRightJustified(fp, osValue.c_str(), nWidth);
489 120 : }
490 :
491 3120 : static void WriteRightJustified(VSILFILE* fp, double dfValue, int nWidth,
492 : int nDecimals = -1)
493 : {
494 : char szFormat[32];
495 3120 : if (nDecimals >= 0)
496 3096 : sprintf(szFormat, "%%.%df", nDecimals);
497 : else
498 24 : sprintf(szFormat, "%%g");
499 3120 : char* pszValue = (char*)CPLSPrintf(szFormat, dfValue);
500 3120 : char* pszE = strchr(pszValue, 'e');
501 3120 : if (pszE)
502 24 : *pszE = 'E';
503 :
504 3120 : if ((int)strlen(pszValue) > nWidth)
505 : {
506 4 : sprintf(szFormat, "%%.%dg", nDecimals);
507 4 : pszValue = (char*)CPLSPrintf(szFormat, dfValue);
508 4 : pszE = strchr(pszValue, 'e');
509 4 : if (pszE)
510 0 : *pszE = 'E';
511 : }
512 :
513 3120 : CPLString osValue(pszValue);
514 3120 : WriteRightJustified(fp, osValue.c_str(), nWidth);
515 3120 : }
516 :
517 : /************************************************************************/
518 : /* CreateCopy() */
519 : /************************************************************************/
520 :
521 38 : 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 38 : int nBands = poSrcDS->GetRasterCount();
531 38 : if (nBands == 0)
532 : {
533 : CPLError( CE_Failure, CPLE_NotSupported,
534 2 : "ZMap driver does not support source dataset with zero band.\n");
535 2 : return NULL;
536 : }
537 :
538 36 : if (nBands != 1)
539 : {
540 : CPLError( (bStrict) ? CE_Failure : CE_Warning, CPLE_NotSupported,
541 8 : "ZMap driver only uses the first band of the dataset.\n");
542 8 : if (bStrict)
543 8 : return NULL;
544 : }
545 :
546 28 : if( pfnProgress && !pfnProgress( 0.0, NULL, pProgressData ) )
547 0 : return NULL;
548 :
549 : /* -------------------------------------------------------------------- */
550 : /* Get source dataset info */
551 : /* -------------------------------------------------------------------- */
552 :
553 28 : int nXSize = poSrcDS->GetRasterXSize();
554 28 : int nYSize = poSrcDS->GetRasterYSize();
555 28 : if (nXSize == 1 || nYSize == 1)
556 : {
557 0 : return NULL;
558 : }
559 :
560 : double adfGeoTransform[6];
561 28 : poSrcDS->GetGeoTransform(adfGeoTransform);
562 28 : 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 28 : VSILFILE* fp = VSIFOpenL(pszFilename, "wb");
574 28 : if (fp == NULL)
575 : {
576 : CPLError( CE_Failure, CPLE_AppDefined,
577 4 : "Cannot create %s", pszFilename );
578 4 : return NULL;
579 : }
580 :
581 24 : int nFieldSize = 20;
582 24 : int nValuesPerLine = 4;
583 24 : int nDecimalCount = 7;
584 :
585 24 : int bHasNoDataValue = FALSE;
586 : double dfNoDataValue =
587 24 : poSrcDS->GetRasterBand(1)->GetNoDataValue(&bHasNoDataValue);
588 24 : if (!bHasNoDataValue)
589 24 : dfNoDataValue = 1.e30;
590 :
591 24 : VSIFPrintfL(fp, "!\n");
592 24 : VSIFPrintfL(fp, "! Created by GDAL.\n");
593 24 : VSIFPrintfL(fp, "!\n");
594 24 : VSIFPrintfL(fp, "@GRID FILE, GRID, %d\n", nValuesPerLine);
595 :
596 24 : WriteRightJustified(fp, nFieldSize, 10);
597 24 : VSIFPrintfL(fp, ",");
598 24 : WriteRightJustified(fp, dfNoDataValue, 10);
599 24 : VSIFPrintfL(fp, ",");
600 24 : WriteRightJustified(fp, "", 10);
601 24 : VSIFPrintfL(fp, ",");
602 24 : WriteRightJustified(fp, nDecimalCount, 10);
603 24 : VSIFPrintfL(fp, ",");
604 24 : WriteRightJustified(fp, 1, 10);
605 24 : VSIFPrintfL(fp, "\n");
606 :
607 24 : WriteRightJustified(fp, nYSize, 10);
608 24 : VSIFPrintfL(fp, ",");
609 24 : WriteRightJustified(fp, nXSize, 10);
610 24 : VSIFPrintfL(fp, ",");
611 :
612 24 : 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 24 : WriteRightJustified(fp, adfGeoTransform[0], 14, 7);
627 24 : VSIFPrintfL(fp, ",");
628 24 : WriteRightJustified(fp, adfGeoTransform[0] + adfGeoTransform[1] * nXSize, 14, 7);
629 24 : VSIFPrintfL(fp, ",");
630 24 : WriteRightJustified(fp, adfGeoTransform[3] + adfGeoTransform[5] * nYSize, 14, 7);
631 24 : VSIFPrintfL(fp, ",");
632 24 : WriteRightJustified(fp, adfGeoTransform[3], 14, 7);
633 : }
634 :
635 24 : VSIFPrintfL(fp, "\n");
636 :
637 24 : VSIFPrintfL(fp, "0.0, 0.0, 0.0\n");
638 24 : VSIFPrintfL(fp, "@\n");
639 :
640 : /* -------------------------------------------------------------------- */
641 : /* Copy imagery */
642 : /* -------------------------------------------------------------------- */
643 24 : double* padfLineBuffer = (double*) CPLMalloc(nYSize * sizeof(double));
644 : int i, j;
645 24 : CPLErr eErr = CE_None;
646 284 : 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 260 : GDT_Float64, 0, 0);
652 260 : if (eErr != CE_None)
653 0 : break;
654 260 : int bEOLPrinted = FALSE;
655 3260 : for(j=0;j<nYSize;j++)
656 : {
657 3000 : WriteRightJustified(fp, padfLineBuffer[j], nFieldSize, nDecimalCount);
658 3000 : if (((j + 1) % nValuesPerLine) == 0)
659 : {
660 640 : bEOLPrinted = TRUE;
661 640 : VSIFPrintfL(fp, "\n");
662 : }
663 : else
664 2360 : bEOLPrinted = FALSE;
665 : }
666 260 : if (!bEOLPrinted)
667 220 : VSIFPrintfL(fp, "\n");
668 :
669 260 : if (!pfnProgress( (j+1) * 1.0 / nYSize, NULL, pProgressData))
670 0 : break;
671 : }
672 24 : CPLFree(padfLineBuffer);
673 24 : VSIFCloseL(fp);
674 :
675 24 : if (eErr != CE_None)
676 0 : return NULL;
677 :
678 24 : return (GDALDataset*) GDALOpen(pszFilename, GA_ReadOnly);
679 : }
680 :
681 : /************************************************************************/
682 : /* GetGeoTransform() */
683 : /************************************************************************/
684 :
685 2 : CPLErr ZMapDataset::GetGeoTransform( double * padfTransform )
686 :
687 : {
688 2 : memcpy(padfTransform, adfGeoTransform, 6 * sizeof(double));
689 :
690 2 : return( CE_None );
691 : }
692 :
693 : /************************************************************************/
694 : /* GDALRegister_ZMap() */
695 : /************************************************************************/
696 :
697 1135 : void GDALRegister_ZMap()
698 :
699 : {
700 : GDALDriver *poDriver;
701 :
702 1135 : if( GDALGetDriverByName( "ZMap" ) == NULL )
703 : {
704 1093 : poDriver = new GDALDriver();
705 :
706 1093 : poDriver->SetDescription( "ZMap" );
707 : poDriver->SetMetadataItem( GDAL_DMD_LONGNAME,
708 1093 : "ZMap Plus Grid" );
709 : poDriver->SetMetadataItem( GDAL_DMD_HELPTOPIC,
710 1093 : "frmt_various.html#ZMap" );
711 1093 : poDriver->SetMetadataItem( GDAL_DMD_EXTENSION, "dat" );
712 :
713 1093 : poDriver->SetMetadataItem( GDAL_DCAP_VIRTUALIO, "YES" );
714 :
715 1093 : poDriver->pfnOpen = ZMapDataset::Open;
716 1093 : poDriver->pfnIdentify = ZMapDataset::Identify;
717 1093 : poDriver->pfnCreateCopy = ZMapDataset::CreateCopy;
718 :
719 1093 : GetGDALDriverManager()->RegisterDriver( poDriver );
720 : }
721 1135 : }
722 :
|
