1 : /******************************************************************************
2 : * $Id: grddataset.cpp 23597 2011-12-18 23:29:48Z rouault $
3 : *
4 : * Project: GRD Reader
5 : * Purpose: GDAL driver for Northwood Grid Format
6 : * Author: Perry Casson
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 2006, Waypoint Information Technology
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_pam.h"
31 : #include "northwood.h"
32 :
33 : #ifdef OGR_ENABLED
34 : #ifdef MSVC
35 : #include "..\..\ogr\ogrsf_frmts\mitab\mitab.h"
36 : #else
37 : #include "../../ogr/ogrsf_frmts/mitab/mitab.h"
38 : #endif
39 : #endif
40 :
41 :
42 : CPL_C_START void GDALRegister_NWT_GRD( void );
43 : CPL_C_END
44 : /************************************************************************/
45 : /* ==================================================================== */
46 : /* NWT_GRDDataset */
47 : /* ==================================================================== */
48 : /************************************************************************/
49 : class NWT_GRDRasterBand;
50 :
51 : class NWT_GRDDataset:public GDALPamDataset
52 : {
53 : friend class NWT_GRDRasterBand;
54 :
55 : VSILFILE *fp;
56 : GByte abyHeader[1024];
57 : NWT_GRID *pGrd;
58 : NWT_RGB ColorMap[4096];
59 : char *pszProjection;
60 :
61 : public:
62 : NWT_GRDDataset();
63 : ~NWT_GRDDataset();
64 :
65 : static GDALDataset *Open( GDALOpenInfo * );
66 : static int Identify( GDALOpenInfo * );
67 :
68 : CPLErr GetGeoTransform( double *padfTransform );
69 : const char *GetProjectionRef();
70 : };
71 :
72 : /************************************************************************/
73 : /* ==================================================================== */
74 : /* NWT_GRDRasterBand */
75 : /* ==================================================================== */
76 : /************************************************************************/
77 :
78 : class NWT_GRDRasterBand:public GDALPamRasterBand
79 24 : {
80 : friend class NWT_GRDDataset;
81 :
82 : int bHaveOffsetScale;
83 : double dfOffset;
84 : double dfScale;
85 :
86 : public:
87 :
88 : NWT_GRDRasterBand( NWT_GRDDataset *, int );
89 :
90 : virtual CPLErr IReadBlock( int, int, void * );
91 : virtual double GetNoDataValue( int *pbSuccess );
92 :
93 : /* FIXME. I don't believe it is correct to advertize offset and */
94 : /* scale because IReadBlock() already apply them. */
95 : virtual double GetOffset( int *pbSuccess = NULL );
96 : virtual CPLErr SetOffset( double dfNewValue );
97 : virtual double GetScale( int *pbSuccess = NULL );
98 : virtual CPLErr SetScale( double dfNewValue );
99 :
100 : virtual GDALColorInterp GetColorInterpretation();
101 : };
102 :
103 :
104 : /************************************************************************/
105 : /* NWT_GRDRasterBand() */
106 : /************************************************************************/
107 24 : NWT_GRDRasterBand::NWT_GRDRasterBand( NWT_GRDDataset * poDS, int nBand )
108 : {
109 24 : this->poDS = poDS;
110 24 : this->nBand = nBand;
111 :
112 24 : if( nBand == 4 )
113 : {
114 6 : bHaveOffsetScale = TRUE;
115 6 : dfOffset = poDS->pGrd->fZMin;
116 :
117 6 : if( poDS->pGrd->cFormat == 0x01 )
118 : {
119 0 : eDataType = GDT_Float32;
120 0 : dfScale =( poDS->pGrd->fZMax - poDS->pGrd->fZMin ) / 4294967294.0;
121 : }
122 : else
123 : {
124 6 : eDataType = GDT_Float32;
125 6 : dfScale =( poDS->pGrd->fZMax - poDS->pGrd->fZMin ) / 65534.0;
126 : }
127 : }
128 : else
129 : {
130 18 : bHaveOffsetScale = FALSE;
131 18 : dfOffset = 0;
132 18 : dfScale = 1.0;
133 18 : eDataType = GDT_Byte;
134 : }
135 24 : nBlockXSize = poDS->GetRasterXSize();
136 24 : nBlockYSize = 1;
137 24 : }
138 :
139 0 : double NWT_GRDRasterBand::GetNoDataValue( int *pbSuccess )
140 : {
141 0 : if (nBand == 4)
142 : {
143 0 : if( pbSuccess != NULL )
144 0 : *pbSuccess = TRUE;
145 :
146 0 : return (float)-1.e37;
147 : }
148 : else
149 : {
150 0 : if( pbSuccess != NULL )
151 0 : *pbSuccess = FALSE;
152 :
153 0 : return 0;
154 : }
155 : }
156 :
157 0 : GDALColorInterp NWT_GRDRasterBand::GetColorInterpretation()
158 : {
159 : //return GCI_RGB;
160 0 : if( nBand == 4 )
161 0 : return GCI_Undefined;
162 0 : else if( nBand == 1 )
163 0 : return GCI_RedBand;
164 0 : else if( nBand == 2 )
165 0 : return GCI_GreenBand;
166 0 : else if( nBand == 3 )
167 0 : return GCI_BlueBand;
168 : else
169 0 : return GCI_Undefined;
170 : }
171 :
172 : /************************************************************************/
173 : /* IReadBlock() */
174 : /************************************************************************/
175 282 : CPLErr NWT_GRDRasterBand::IReadBlock( int nBlockXOff, int nBlockYOff, void *pImage )
176 : {
177 282 : NWT_GRDDataset *poGDS = (NWT_GRDDataset *) poDS;
178 : char *pszRecord;
179 282 : int nRecordSize = nBlockXSize * 2;
180 : int i;
181 : unsigned short raw1;
182 :
183 282 : VSIFSeekL( poGDS->fp, 1024 + nRecordSize * nBlockYOff, SEEK_SET );
184 :
185 282 : pszRecord = (char *) CPLMalloc( nRecordSize );
186 282 : VSIFReadL( pszRecord, 1, nRecordSize, poGDS->fp );
187 :
188 282 : if( nBand == 4 ) //Z values
189 : {
190 0 : for( i = 0; i < nBlockXSize; i++ )
191 : {
192 0 : memcpy( (void *) &raw1, (void *)(pszRecord + 2 * i), 2 );
193 : CPL_LSBPTR16(&raw1);
194 0 : if( raw1 == 0 )
195 : {
196 0 : ((float *)pImage)[i] = (float)-1.e37; // null value
197 : }
198 : else
199 : {
200 0 : ((float *)pImage)[i] = (float) (dfOffset + ((raw1 - 1) * dfScale));
201 : }
202 : }
203 : }
204 282 : else if( nBand == 1 ) // red values
205 : {
206 6016 : for( i = 0; i < nBlockXSize; i++ )
207 : {
208 5922 : memcpy( (void *) &raw1, (void *)(pszRecord + 2 * i), 2 );
209 : CPL_LSBPTR16(&raw1);
210 5922 : ((char *)pImage)[i] = poGDS->ColorMap[raw1 / 16].r;
211 : }
212 : }
213 188 : else if( nBand == 2 ) // green
214 : {
215 6016 : for( i = 0; i < nBlockXSize; i++ )
216 : {
217 5922 : memcpy( (void *) &raw1, (void *)(pszRecord + 2 * i), 2 );
218 : CPL_LSBPTR16(&raw1);
219 5922 : ((char *) pImage)[i] = poGDS->ColorMap[raw1 / 16].g;
220 : }
221 : }
222 94 : else if( nBand == 3 ) // blue
223 : {
224 6016 : for( i = 0; i < nBlockXSize; i++ )
225 : {
226 5922 : memcpy( (void *) &raw1, (void *) (pszRecord + 2 * i), 2 );
227 : CPL_LSBPTR16(&raw1);
228 5922 : ((char *) pImage)[i] = poGDS->ColorMap[raw1 / 16].b;
229 : }
230 : }
231 : else
232 : {
233 : CPLError( CE_Failure, CPLE_IllegalArg,
234 : "No band number %d",
235 0 : nBand );
236 0 : if( pszRecord != NULL )
237 0 : CPLFree( pszRecord );
238 0 : return CE_Failure;
239 : }
240 282 : if( pszRecord != NULL )
241 : {
242 282 : CPLFree( pszRecord );
243 : }
244 282 : return CE_None;
245 : }
246 :
247 : /************************************************************************/
248 : /* GetOffset() */
249 : /************************************************************************/
250 0 : double NWT_GRDRasterBand::GetOffset( int *pbSuccess )
251 : {
252 0 : if( pbSuccess )
253 0 : *pbSuccess = bHaveOffsetScale;
254 0 : return dfOffset;
255 : }
256 :
257 : /************************************************************************/
258 : /* SetOffset() */
259 : /************************************************************************/
260 0 : CPLErr NWT_GRDRasterBand::SetOffset( double dfNewValue )
261 : {
262 : //poGDS->bMetadataChanged = TRUE;
263 :
264 0 : bHaveOffsetScale = TRUE;
265 0 : dfOffset = dfNewValue;
266 0 : return CE_None;
267 : }
268 :
269 : /************************************************************************/
270 : /* GetScale() */
271 : /************************************************************************/
272 0 : double NWT_GRDRasterBand::GetScale( int *pbSuccess )
273 : {
274 0 : if( pbSuccess )
275 0 : *pbSuccess = bHaveOffsetScale;
276 0 : return dfScale;
277 : }
278 :
279 : /************************************************************************/
280 : /* SetScale() */
281 : /************************************************************************/
282 0 : CPLErr NWT_GRDRasterBand::SetScale( double dfNewValue )
283 : {
284 0 : bHaveOffsetScale = TRUE;
285 0 : dfScale = dfNewValue;
286 0 : return CE_None;
287 : }
288 :
289 : /************************************************************************/
290 : /* ==================================================================== */
291 : /* NWT_GRDDataset */
292 : /* ==================================================================== */
293 : /************************************************************************/
294 6 : NWT_GRDDataset::NWT_GRDDataset()
295 : {
296 6 : pszProjection = NULL;
297 : //poCT = NULL;
298 6 : }
299 :
300 :
301 : /************************************************************************/
302 : /* ~NWT_GRDDataset() */
303 : /************************************************************************/
304 :
305 6 : NWT_GRDDataset::~NWT_GRDDataset()
306 : {
307 6 : FlushCache();
308 6 : pGrd->fp = NULL; // this prevents nwtCloseGrid from closing the fp
309 6 : nwtCloseGrid( pGrd );
310 :
311 6 : if( fp != NULL )
312 6 : VSIFCloseL( fp );
313 :
314 6 : if( pszProjection != NULL )
315 : {
316 0 : CPLFree( pszProjection );
317 : }
318 : /*if( poCT != NULL )
319 : delete poCT;*/
320 6 : }
321 :
322 : /************************************************************************/
323 : /* GetGeoTransform() */
324 : /************************************************************************/
325 0 : CPLErr NWT_GRDDataset::GetGeoTransform( double *padfTransform )
326 : {
327 0 : padfTransform[0] = pGrd->dfMinX - ( pGrd->dfStepSize * 0.5 );
328 0 : padfTransform[3] = pGrd->dfMaxY + ( pGrd->dfStepSize * 0.5 );
329 0 : padfTransform[1] = pGrd->dfStepSize;
330 0 : padfTransform[2] = 0.0;
331 :
332 0 : padfTransform[4] = 0.0;
333 0 : padfTransform[5] = -1 * pGrd->dfStepSize;
334 :
335 0 : return CE_None;
336 : }
337 :
338 : /************************************************************************/
339 : /* GetProjectionRef() */
340 : /************************************************************************/
341 0 : const char *NWT_GRDDataset::GetProjectionRef()
342 : {
343 : #ifdef OGR_ENABLED
344 0 : if (pszProjection == NULL)
345 : {
346 : OGRSpatialReference *poSpatialRef;
347 0 : poSpatialRef = MITABCoordSys2SpatialRef( pGrd->cMICoordSys );
348 0 : if (poSpatialRef)
349 : {
350 0 : poSpatialRef->exportToWkt( &pszProjection );
351 0 : poSpatialRef->Release();
352 : }
353 : }
354 : #endif
355 0 : return( (const char *)pszProjection );
356 : }
357 :
358 : /************************************************************************/
359 : /* Identify() */
360 : /************************************************************************/
361 :
362 22240 : int NWT_GRDDataset::Identify( GDALOpenInfo * poOpenInfo )
363 : {
364 : /* -------------------------------------------------------------------- */
365 : /* Look for the header */
366 : /* -------------------------------------------------------------------- */
367 22240 : if( poOpenInfo->nHeaderBytes < 50 )
368 21244 : return FALSE;
369 :
370 1070 : if( poOpenInfo->pabyHeader[0] != 'H' ||
371 50 : poOpenInfo->pabyHeader[1] != 'G' ||
372 8 : poOpenInfo->pabyHeader[2] != 'P' ||
373 8 : poOpenInfo->pabyHeader[3] != 'C' ||
374 8 : poOpenInfo->pabyHeader[4] != '1' )
375 990 : return FALSE;
376 :
377 6 : return TRUE;
378 : }
379 :
380 : /************************************************************************/
381 : /* Open() */
382 : /************************************************************************/
383 3332 : GDALDataset *NWT_GRDDataset::Open( GDALOpenInfo * poOpenInfo )
384 : {
385 3332 : if( !Identify(poOpenInfo) )
386 3326 : return NULL;
387 :
388 : /* -------------------------------------------------------------------- */
389 : /* Create a corresponding GDALDataset. */
390 : /* -------------------------------------------------------------------- */
391 : NWT_GRDDataset *poDS;
392 :
393 6 : poDS = new NWT_GRDDataset();
394 :
395 6 : poDS->fp = VSIFOpenL(poOpenInfo->pszFilename, "rb");
396 6 : if (poDS->fp == NULL)
397 : {
398 0 : delete poDS;
399 0 : return NULL;
400 : }
401 :
402 : /* -------------------------------------------------------------------- */
403 : /* Read the header. */
404 : /* -------------------------------------------------------------------- */
405 6 : VSIFSeekL( poDS->fp, 0, SEEK_SET );
406 6 : VSIFReadL( poDS->abyHeader, 1, 1024, poDS->fp );
407 6 : poDS->pGrd = (NWT_GRID *) malloc(sizeof(NWT_GRID));
408 :
409 6 : poDS->pGrd->fp = poDS->fp;
410 :
411 6 : if (!nwt_ParseHeader( poDS->pGrd, (char *) poDS->abyHeader ) ||
412 : !GDALCheckDatasetDimensions(poDS->pGrd->nXSide, poDS->pGrd->nYSide) )
413 : {
414 0 : delete poDS;
415 0 : return NULL;
416 : }
417 :
418 6 : poDS->nRasterXSize = poDS->pGrd->nXSide;
419 6 : poDS->nRasterYSize = poDS->pGrd->nYSide;
420 :
421 : // create a colorTable
422 : // if( poDS->pGrd->iNumColorInflections > 0 )
423 : // poDS->CreateColorTable();
424 6 : nwt_LoadColors( poDS->ColorMap, 4096, poDS->pGrd );
425 : /* -------------------------------------------------------------------- */
426 : /* Create band information objects. */
427 : /* -------------------------------------------------------------------- */
428 6 : poDS->SetBand( 1, new NWT_GRDRasterBand( poDS, 1 ) ); //r
429 12 : poDS->SetBand( 2, new NWT_GRDRasterBand( poDS, 2 ) ); //g
430 12 : poDS->SetBand( 3, new NWT_GRDRasterBand( poDS, 3 ) ); //b
431 12 : poDS->SetBand( 4, new NWT_GRDRasterBand( poDS, 4 ) ); //z
432 :
433 : /* -------------------------------------------------------------------- */
434 : /* Initialize any PAM information. */
435 : /* -------------------------------------------------------------------- */
436 6 : poDS->SetDescription( poOpenInfo->pszFilename );
437 6 : poDS->TryLoadXML();
438 :
439 : /* -------------------------------------------------------------------- */
440 : /* Check for external overviews. */
441 : /* -------------------------------------------------------------------- */
442 6 : poDS->oOvManager.Initialize( poDS, poOpenInfo->pszFilename, poOpenInfo->papszSiblingFiles );
443 :
444 6 : return (poDS);
445 : }
446 :
447 :
448 : /************************************************************************/
449 : /* GDALRegister_GRD() */
450 : /************************************************************************/
451 1135 : void GDALRegister_NWT_GRD()
452 : {
453 : GDALDriver *poDriver;
454 :
455 1135 : if( GDALGetDriverByName( "NWT_GRD" ) == NULL )
456 : {
457 1093 : poDriver = new GDALDriver();
458 :
459 1093 : poDriver->SetDescription( "NWT_GRD" );
460 : poDriver->SetMetadataItem( GDAL_DMD_LONGNAME,
461 1093 : "Northwood Numeric Grid Format .grd/.tab" );
462 1093 : poDriver->SetMetadataItem( GDAL_DMD_HELPTOPIC, "frmt_various.html#grd");
463 1093 : poDriver->SetMetadataItem( GDAL_DMD_EXTENSION, "grd" );
464 1093 : poDriver->SetMetadataItem( GDAL_DCAP_VIRTUALIO, "YES" );
465 :
466 1093 : poDriver->pfnOpen = NWT_GRDDataset::Open;
467 1093 : poDriver->pfnIdentify = NWT_GRDDataset::Identify;
468 :
469 1093 : GetGDALDriverManager()->RegisterDriver( poDriver );
470 : }
471 1135 : }
|
