1 : /******************************************************************************
2 : * $Id: hfaband.cpp 24206 2012-04-07 15:43:21Z rouault $
3 : *
4 : * Project: Erdas Imagine (.img) Translator
5 : * Purpose: Implementation of the HFABand, for accessing one Eimg_Layer.
6 : * Author: Frank Warmerdam, warmerdam@pobox.com
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 1999, Intergraph Corporation
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 "hfa_p.h"
31 : #include "cpl_conv.h"
32 :
33 : /* include the compression code */
34 :
35 : CPL_CVSID("$Id: hfaband.cpp 24206 2012-04-07 15:43:21Z rouault $");
36 :
37 : /************************************************************************/
38 : /* HFABand() */
39 : /************************************************************************/
40 :
41 690 : HFABand::HFABand( HFAInfo_t * psInfoIn, HFAEntry * poNodeIn )
42 :
43 : {
44 690 : psInfo = psInfoIn;
45 690 : poNode = poNodeIn;
46 :
47 690 : bOverviewsPending = TRUE;
48 :
49 690 : nBlockXSize = poNodeIn->GetIntField( "blockWidth" );
50 690 : nBlockYSize = poNodeIn->GetIntField( "blockHeight" );
51 690 : nDataType = poNodeIn->GetIntField( "pixelType" );
52 :
53 690 : nWidth = poNodeIn->GetIntField( "width" );
54 690 : nHeight = poNodeIn->GetIntField( "height" );
55 :
56 690 : panBlockStart = NULL;
57 690 : panBlockSize = NULL;
58 690 : panBlockFlag = NULL;
59 :
60 690 : nPCTColors = -1;
61 690 : apadfPCT[0] = apadfPCT[1] = apadfPCT[2] = apadfPCT[3] = NULL;
62 690 : padfPCTBins = NULL;
63 :
64 690 : nOverviews = 0;
65 690 : papoOverviews = NULL;
66 :
67 690 : fpExternal = NULL;
68 :
69 690 : bNoDataSet = FALSE;
70 690 : dfNoData = 0.0;
71 :
72 690 : if (nWidth <= 0 || nHeight <= 0 || nBlockXSize <= 0 || nBlockYSize <= 0)
73 : {
74 0 : nWidth = nHeight = 0;
75 : CPLError(CE_Failure, CPLE_AppDefined,
76 0 : "HFABand::HFABand : (nWidth <= 0 || nHeight <= 0 || nBlockXSize <= 0 || nBlockYSize <= 0)");
77 0 : return;
78 : }
79 690 : if (HFAGetDataTypeBits(nDataType) == 0)
80 : {
81 0 : nWidth = nHeight = 0;
82 : CPLError(CE_Failure, CPLE_AppDefined,
83 0 : "HFABand::HFABand : nDataType=%d unhandled", nDataType);
84 0 : return;
85 : }
86 :
87 : /* FIXME? : risk of overflow in additions and multiplication */
88 690 : nBlocksPerRow = (nWidth + nBlockXSize - 1) / nBlockXSize;
89 690 : nBlocksPerColumn = (nHeight + nBlockYSize - 1) / nBlockYSize;
90 690 : nBlocks = nBlocksPerRow * nBlocksPerColumn;
91 :
92 : /* -------------------------------------------------------------------- */
93 : /* Check for nodata. This is really an RDO (ESRI Raster Data */
94 : /* Objects?), not used by Imagine itself. */
95 : /* -------------------------------------------------------------------- */
96 690 : HFAEntry *poNDNode = poNode->GetNamedChild("Eimg_NonInitializedValue");
97 :
98 690 : if( poNDNode != NULL )
99 : {
100 18 : bNoDataSet = TRUE;
101 18 : dfNoData = poNDNode->GetDoubleField( "valueBD" );
102 : }
103 :
104 0 : }
105 :
106 : /************************************************************************/
107 : /* ~HFABand() */
108 : /************************************************************************/
109 :
110 690 : HFABand::~HFABand()
111 :
112 : {
113 748 : for( int iOverview = 0; iOverview < nOverviews; iOverview++ )
114 58 : delete papoOverviews[iOverview];
115 :
116 690 : if( nOverviews > 0 )
117 36 : CPLFree( papoOverviews );
118 :
119 690 : if ( panBlockStart )
120 141 : CPLFree( panBlockStart );
121 690 : if ( panBlockSize )
122 141 : CPLFree( panBlockSize );
123 690 : if ( panBlockFlag )
124 157 : CPLFree( panBlockFlag );
125 :
126 690 : CPLFree( apadfPCT[0] );
127 690 : CPLFree( apadfPCT[1] );
128 690 : CPLFree( apadfPCT[2] );
129 690 : CPLFree( apadfPCT[3] );
130 690 : CPLFree( padfPCTBins );
131 :
132 690 : if( fpExternal != NULL )
133 16 : VSIFCloseL( fpExternal );
134 690 : }
135 :
136 : /************************************************************************/
137 : /* LoadOverviews() */
138 : /************************************************************************/
139 :
140 168 : CPLErr HFABand::LoadOverviews()
141 :
142 : {
143 168 : if( !bOverviewsPending )
144 56 : return CE_None;
145 :
146 112 : bOverviewsPending = FALSE;
147 :
148 : /* -------------------------------------------------------------------- */
149 : /* Does this band have overviews? Try to find them. */
150 : /* -------------------------------------------------------------------- */
151 112 : HFAEntry *poRRDNames = poNode->GetNamedChild( "RRDNamesList" );
152 :
153 112 : if( poRRDNames != NULL )
154 : {
155 83 : for( int iName = 0; TRUE; iName++ )
156 : {
157 : char szField[128], *pszPath, *pszFilename, *pszEnd;
158 : const char *pszName;
159 : CPLErr eErr;
160 : HFAEntry *poOvEntry;
161 : int i;
162 : HFAInfo_t *psHFA;
163 :
164 83 : sprintf( szField, "nameList[%d].string", iName );
165 :
166 83 : pszName = poRRDNames->GetStringField( szField, &eErr );
167 83 : if( pszName == NULL || eErr != CE_None )
168 33 : break;
169 :
170 50 : pszFilename = CPLStrdup(pszName);
171 50 : pszEnd = strstr(pszFilename,"(:");
172 50 : if( pszEnd == NULL )
173 : {
174 0 : CPLFree( pszFilename );
175 0 : continue;
176 : }
177 :
178 50 : pszName = pszEnd + 2;
179 50 : pszEnd[0] = '\0';
180 :
181 : char *pszJustFilename;
182 :
183 50 : pszJustFilename = CPLStrdup(CPLGetFilename(pszFilename));
184 50 : psHFA = HFAGetDependent( psInfo, pszJustFilename );
185 50 : CPLFree( pszJustFilename );
186 :
187 : // Try finding the dependent file as this file with the
188 : // extension .rrd. This is intended to address problems
189 : // with users changing the names of their files.
190 50 : if( psHFA == NULL )
191 : {
192 : char *pszBasename =
193 2 : CPLStrdup(CPLGetBasename(psInfo->pszFilename));
194 :
195 : pszJustFilename =
196 2 : CPLStrdup(CPLFormFilename(NULL, pszBasename, "rrd"));
197 : CPLDebug( "HFA", "Failed to find overview file with expected name,\ntry %s instead.",
198 2 : pszJustFilename );
199 2 : psHFA = HFAGetDependent( psInfo, pszJustFilename );
200 2 : CPLFree( pszJustFilename );
201 2 : CPLFree( pszBasename );
202 : }
203 :
204 50 : if( psHFA == NULL )
205 : {
206 2 : CPLFree( pszFilename );
207 2 : continue;
208 : }
209 :
210 48 : pszPath = pszEnd + 2;
211 48 : if( pszPath[strlen(pszPath)-1] == ')' )
212 48 : pszPath[strlen(pszPath)-1] = '\0';
213 :
214 723 : for( i=0; pszPath[i] != '\0'; i++ )
215 : {
216 675 : if( pszPath[i] == ':' )
217 48 : pszPath[i] = '.';
218 : }
219 :
220 48 : poOvEntry = psHFA->poRoot->GetNamedChild( pszPath );
221 48 : CPLFree( pszFilename );
222 :
223 48 : if( poOvEntry == NULL )
224 0 : continue;
225 :
226 : /*
227 : * We have an overview node. Instanatiate a HFABand from it,
228 : * and add to the list.
229 : */
230 : papoOverviews = (HFABand **)
231 48 : CPLRealloc(papoOverviews, sizeof(void*) * ++nOverviews );
232 48 : papoOverviews[nOverviews-1] = new HFABand( psHFA, poOvEntry );
233 48 : if (papoOverviews[nOverviews-1]->nWidth == 0)
234 : {
235 0 : nWidth = nHeight = 0;
236 0 : delete papoOverviews[nOverviews-1];
237 0 : papoOverviews[nOverviews-1] = NULL;
238 0 : return CE_None;
239 : }
240 : }
241 : }
242 :
243 : /* -------------------------------------------------------------------- */
244 : /* If there are no overviews mentioned in this file, probe for */
245 : /* an .rrd file anyways. */
246 : /* -------------------------------------------------------------------- */
247 112 : HFAEntry *poBandProxyNode = poNode;
248 112 : HFAInfo_t *psOvHFA = psInfo;
249 :
250 112 : if( nOverviews == 0
251 : && EQUAL(CPLGetExtension(psInfo->pszFilename),"aux") )
252 : {
253 3 : CPLString osRRDFilename = CPLResetExtension( psInfo->pszFilename,"rrd");
254 : CPLString osFullRRD = CPLFormFilename( psInfo->pszPath, osRRDFilename,
255 3 : NULL );
256 : VSIStatBufL sStatBuf;
257 :
258 3 : if( VSIStatL( osFullRRD, &sStatBuf ) == 0 )
259 : {
260 0 : psOvHFA = HFAGetDependent( psInfo, osRRDFilename );
261 0 : if( psOvHFA )
262 : poBandProxyNode =
263 0 : psOvHFA->poRoot->GetNamedChild( poNode->GetName() );
264 : else
265 0 : psOvHFA = psInfo;
266 3 : }
267 : }
268 :
269 : /* -------------------------------------------------------------------- */
270 : /* If there are no named overviews, try looking for unnamed */
271 : /* overviews within the same layer, as occurs in floodplain.img */
272 : /* for instance, or in the not-referenced rrd mentioned in #3463. */
273 : /* -------------------------------------------------------------------- */
274 112 : if( nOverviews == 0 && poBandProxyNode != NULL )
275 : {
276 : HFAEntry *poChild;
277 :
278 371 : for( poChild = poBandProxyNode->GetChild();
279 : poChild != NULL;
280 : poChild = poChild->GetNext() )
281 : {
282 290 : if( EQUAL(poChild->GetType(),"Eimg_Layer_SubSample") )
283 : {
284 : papoOverviews = (HFABand **)
285 0 : CPLRealloc(papoOverviews, sizeof(void*) * ++nOverviews );
286 0 : papoOverviews[nOverviews-1] = new HFABand( psOvHFA, poChild );
287 0 : if (papoOverviews[nOverviews-1]->nWidth == 0)
288 : {
289 0 : nWidth = nHeight = 0;
290 0 : delete papoOverviews[nOverviews-1];
291 0 : papoOverviews[nOverviews-1] = NULL;
292 0 : return CE_None;
293 : }
294 : }
295 : }
296 :
297 : int i1, i2;
298 :
299 : // bubble sort into biggest to smallest order.
300 81 : for( i1 = 0; i1 < nOverviews; i1++ )
301 : {
302 0 : for( i2 = 0; i2 < nOverviews-1; i2++ )
303 : {
304 0 : if( papoOverviews[i2]->nWidth <
305 0 : papoOverviews[i2+1]->nWidth )
306 : {
307 0 : HFABand *poTemp = papoOverviews[i2+1];
308 0 : papoOverviews[i2+1] = papoOverviews[i2];
309 0 : papoOverviews[i2] = poTemp;
310 : }
311 : }
312 : }
313 : }
314 112 : return CE_None;
315 : }
316 :
317 : /************************************************************************/
318 : /* LoadBlockInfo() */
319 : /************************************************************************/
320 :
321 505 : CPLErr HFABand::LoadBlockInfo()
322 :
323 : {
324 : int iBlock;
325 : HFAEntry *poDMS;
326 :
327 505 : if( panBlockFlag != NULL )
328 348 : return( CE_None );
329 :
330 157 : poDMS = poNode->GetNamedChild( "RasterDMS" );
331 157 : if( poDMS == NULL )
332 : {
333 16 : if( poNode->GetNamedChild( "ExternalRasterDMS" ) != NULL )
334 16 : return LoadExternalBlockInfo();
335 :
336 : CPLError( CE_Failure, CPLE_AppDefined,
337 0 : "Can't find RasterDMS field in Eimg_Layer with block list.\n");
338 :
339 0 : return CE_Failure;
340 : }
341 :
342 141 : panBlockStart = (vsi_l_offset *)VSIMalloc2(sizeof(vsi_l_offset), nBlocks);
343 141 : panBlockSize = (int *) VSIMalloc2(sizeof(int), nBlocks);
344 141 : panBlockFlag = (int *) VSIMalloc2(sizeof(int), nBlocks);
345 :
346 141 : if (panBlockStart == NULL || panBlockSize == NULL || panBlockFlag == NULL)
347 : {
348 : CPLError( CE_Failure, CPLE_OutOfMemory,
349 0 : "HFABand::LoadBlockInfo : Out of memory\n");
350 :
351 0 : CPLFree(panBlockStart);
352 0 : CPLFree(panBlockSize);
353 0 : CPLFree(panBlockFlag);
354 0 : panBlockStart = NULL;
355 0 : panBlockSize = NULL;
356 0 : panBlockFlag = NULL;
357 0 : return CE_Failure;
358 : }
359 :
360 508 : for( iBlock = 0; iBlock < nBlocks; iBlock++ )
361 : {
362 : char szVarName[64];
363 : int nLogvalid, nCompressType;
364 :
365 367 : sprintf( szVarName, "blockinfo[%d].offset", iBlock );
366 367 : panBlockStart[iBlock] = (GUInt32)poDMS->GetIntField( szVarName );
367 :
368 367 : sprintf( szVarName, "blockinfo[%d].size", iBlock );
369 367 : panBlockSize[iBlock] = poDMS->GetIntField( szVarName );
370 :
371 367 : sprintf( szVarName, "blockinfo[%d].logvalid", iBlock );
372 367 : nLogvalid = poDMS->GetIntField( szVarName );
373 :
374 367 : sprintf( szVarName, "blockinfo[%d].compressionType", iBlock );
375 367 : nCompressType = poDMS->GetIntField( szVarName );
376 :
377 367 : panBlockFlag[iBlock] = 0;
378 367 : if( nLogvalid )
379 200 : panBlockFlag[iBlock] |= BFLG_VALID;
380 367 : if( nCompressType != 0 )
381 114 : panBlockFlag[iBlock] |= BFLG_COMPRESSED;
382 : }
383 :
384 141 : return( CE_None );
385 : }
386 :
387 : /************************************************************************/
388 : /* LoadExternalBlockInfo() */
389 : /************************************************************************/
390 :
391 16 : CPLErr HFABand::LoadExternalBlockInfo()
392 :
393 : {
394 : int iBlock;
395 : HFAEntry *poDMS;
396 :
397 16 : if( panBlockFlag != NULL )
398 0 : return( CE_None );
399 :
400 : /* -------------------------------------------------------------------- */
401 : /* Get the info structure. */
402 : /* -------------------------------------------------------------------- */
403 16 : poDMS = poNode->GetNamedChild( "ExternalRasterDMS" );
404 16 : CPLAssert( poDMS != NULL );
405 :
406 16 : nLayerStackCount = poDMS->GetIntField( "layerStackCount" );
407 16 : nLayerStackIndex = poDMS->GetIntField( "layerStackIndex" );
408 :
409 : /* -------------------------------------------------------------------- */
410 : /* Open raw data file. */
411 : /* -------------------------------------------------------------------- */
412 16 : const char *pszFullFilename = HFAGetIGEFilename( psInfo );
413 16 : if (pszFullFilename == NULL)
414 : {
415 : CPLError( CE_Failure, CPLE_OpenFailed,
416 0 : "Cannot find external data file name" );
417 0 : return CE_Failure;
418 : }
419 :
420 16 : if( psInfo->eAccess == HFA_ReadOnly )
421 9 : fpExternal = VSIFOpenL( pszFullFilename, "rb" );
422 : else
423 7 : fpExternal = VSIFOpenL( pszFullFilename, "r+b" );
424 16 : if( fpExternal == NULL )
425 : {
426 : CPLError( CE_Failure, CPLE_OpenFailed,
427 : "Unable to open external data file:\n%s\n",
428 0 : pszFullFilename );
429 0 : return CE_Failure;
430 : }
431 :
432 : /* -------------------------------------------------------------------- */
433 : /* Verify header. */
434 : /* -------------------------------------------------------------------- */
435 : char szHeader[49];
436 :
437 16 : VSIFReadL( szHeader, 49, 1, fpExternal );
438 :
439 16 : if( strncmp( szHeader, "ERDAS_IMG_EXTERNAL_RASTER", 26 ) != 0 )
440 : {
441 : CPLError( CE_Failure, CPLE_AppDefined,
442 : "Raw data file %s appears to be corrupt.\n",
443 0 : pszFullFilename );
444 0 : return CE_Failure;
445 : }
446 :
447 : /* -------------------------------------------------------------------- */
448 : /* Allocate blockmap. */
449 : /* -------------------------------------------------------------------- */
450 16 : panBlockFlag = (int *) VSIMalloc2(sizeof(int), nBlocks);
451 16 : if (panBlockFlag == NULL)
452 : {
453 : CPLError( CE_Failure, CPLE_OutOfMemory,
454 0 : "HFABand::LoadExternalBlockInfo : Out of memory\n");
455 0 : return CE_Failure;
456 : }
457 :
458 : /* -------------------------------------------------------------------- */
459 : /* Load the validity bitmap. */
460 : /* -------------------------------------------------------------------- */
461 : unsigned char *pabyBlockMap;
462 : int nBytesPerRow;
463 :
464 16 : nBytesPerRow = (nBlocksPerRow + 7) / 8;
465 : pabyBlockMap = (unsigned char *)
466 16 : VSIMalloc(nBytesPerRow*nBlocksPerColumn+20);
467 16 : if (pabyBlockMap == NULL)
468 : {
469 : CPLError( CE_Failure, CPLE_OutOfMemory,
470 0 : "HFABand::LoadExternalBlockInfo : Out of memory\n");
471 0 : return CE_Failure;
472 : }
473 :
474 : VSIFSeekL( fpExternal,
475 : poDMS->GetBigIntField( "layerStackValidFlagsOffset" ),
476 16 : SEEK_SET );
477 :
478 16 : if( VSIFReadL( pabyBlockMap, nBytesPerRow * nBlocksPerColumn + 20, 1,
479 : fpExternal ) != 1 )
480 : {
481 : CPLError( CE_Failure, CPLE_FileIO,
482 0 : "Failed to read block validity map." );
483 0 : return CE_Failure;
484 : }
485 :
486 : /* -------------------------------------------------------------------- */
487 : /* Establish block information. Block position is computed */
488 : /* from data base address. Blocks are never compressed. */
489 : /* Validity is determined from the validity bitmap. */
490 : /* -------------------------------------------------------------------- */
491 16 : nBlockStart = poDMS->GetBigIntField( "layerStackDataOffset" );
492 16 : nBlockSize = (nBlockXSize*nBlockYSize*HFAGetDataTypeBits(nDataType)+7) / 8;
493 :
494 32 : for( iBlock = 0; iBlock < nBlocks; iBlock++ )
495 : {
496 : int nRow, nColumn, nBit;
497 :
498 16 : nColumn = iBlock % nBlocksPerRow;
499 16 : nRow = iBlock / nBlocksPerRow;
500 16 : nBit = nRow * nBytesPerRow * 8 + nColumn + 20 * 8;
501 :
502 16 : if( (pabyBlockMap[nBit>>3] >> (nBit&7)) & 0x1 )
503 16 : panBlockFlag[iBlock] = BFLG_VALID;
504 : else
505 0 : panBlockFlag[iBlock] = 0;
506 : }
507 :
508 16 : CPLFree( pabyBlockMap );
509 :
510 16 : return( CE_None );
511 : }
512 :
513 : /************************************************************************/
514 : /* UncompressBlock() */
515 : /* */
516 : /* Uncompress ESRI Grid compression format block. */
517 : /************************************************************************/
518 :
519 : #define CHECK_ENOUGH_BYTES(n) \
520 : if (nSrcBytes < (n)) goto not_enough_bytes;
521 :
522 108 : static CPLErr UncompressBlock( GByte *pabyCData, int nSrcBytes,
523 : GByte *pabyDest, int nMaxPixels,
524 : int nDataType )
525 :
526 : {
527 : GUInt32 nDataMin;
528 108 : int nNumBits, nPixelsOutput=0;
529 : GInt32 nNumRuns, nDataOffset;
530 : GByte *pabyCounter, *pabyValues;
531 : int nValueBitOffset;
532 : int nCounterOffset;
533 :
534 108 : CHECK_ENOUGH_BYTES(13);
535 :
536 108 : memcpy( &nDataMin, pabyCData, 4 );
537 108 : nDataMin = CPL_LSBWORD32( nDataMin );
538 :
539 108 : memcpy( &nNumRuns, pabyCData+4, 4 );
540 108 : nNumRuns = CPL_LSBWORD32( nNumRuns );
541 :
542 108 : memcpy( &nDataOffset, pabyCData+8, 4 );
543 108 : nDataOffset = CPL_LSBWORD32( nDataOffset );
544 :
545 108 : nNumBits = pabyCData[12];
546 :
547 : /* ==================================================================== */
548 : /* If this is not run length encoded, but just reduced */
549 : /* precision, handle it now. */
550 : /* ==================================================================== */
551 108 : if( nNumRuns == -1 )
552 : {
553 20 : pabyValues = pabyCData + 13;
554 20 : nValueBitOffset = 0;
555 :
556 20 : if (nNumBits > INT_MAX / nMaxPixels ||
557 : nNumBits * nMaxPixels > INT_MAX - 7 ||
558 : (nNumBits * nMaxPixels + 7)/8 > INT_MAX - 13)
559 : {
560 0 : CPLError(CE_Failure, CPLE_AppDefined, "Integer overflow : nNumBits * nMaxPixels + 7");
561 0 : return CE_Failure;
562 : }
563 20 : CHECK_ENOUGH_BYTES(13 + (nNumBits * nMaxPixels + 7)/8);
564 :
565 : /* -------------------------------------------------------------------- */
566 : /* Loop over block pixels. */
567 : /* -------------------------------------------------------------------- */
568 81940 : for( nPixelsOutput = 0; nPixelsOutput < nMaxPixels; nPixelsOutput++ )
569 : {
570 : int nDataValue, nRawValue;
571 :
572 : /* -------------------------------------------------------------------- */
573 : /* Extract the data value in a way that depends on the number */
574 : /* of bits in it. */
575 : /* -------------------------------------------------------------------- */
576 81920 : if( nNumBits == 0 )
577 : {
578 0 : nRawValue = 0;
579 : }
580 81920 : else if( nNumBits == 1 )
581 : {
582 : nRawValue =
583 0 : (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0x1;
584 0 : nValueBitOffset++;
585 : }
586 81920 : else if( nNumBits == 2 )
587 : {
588 : nRawValue =
589 0 : (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0x3;
590 0 : nValueBitOffset += 2;
591 : }
592 81920 : else if( nNumBits == 4 )
593 : {
594 : nRawValue =
595 0 : (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0xf;
596 0 : nValueBitOffset += 4;
597 : }
598 81920 : else if( nNumBits == 8 )
599 : {
600 16384 : nRawValue = *pabyValues;
601 16384 : pabyValues++;
602 : }
603 65536 : else if( nNumBits == 16 )
604 : {
605 65536 : nRawValue = 256 * *(pabyValues++);
606 65536 : nRawValue += *(pabyValues++);
607 : }
608 0 : else if( nNumBits == 32 )
609 : {
610 0 : nRawValue = 256 * 256 * 256 * *(pabyValues++);
611 0 : nRawValue += 256 * 256 * *(pabyValues++);
612 0 : nRawValue += 256 * *(pabyValues++);
613 0 : nRawValue += *(pabyValues++);
614 : }
615 : else
616 : {
617 : CPLError(CE_Failure, CPLE_NotSupported,
618 0 : "Unsupported nNumBits value : %d", nNumBits);
619 0 : return CE_Failure;
620 : }
621 :
622 : /* -------------------------------------------------------------------- */
623 : /* Offset by the minimum value. */
624 : /* -------------------------------------------------------------------- */
625 81920 : nDataValue = nRawValue + nDataMin;
626 :
627 : /* -------------------------------------------------------------------- */
628 : /* Now apply to the output buffer in a type specific way. */
629 : /* -------------------------------------------------------------------- */
630 81920 : if( nDataType == EPT_u8 )
631 : {
632 0 : ((GByte *) pabyDest)[nPixelsOutput] = (GByte) nDataValue;
633 : }
634 81920 : else if( nDataType == EPT_u1 )
635 : {
636 0 : if( nDataValue == 1 )
637 0 : pabyDest[nPixelsOutput>>3] |= (1 << (nPixelsOutput & 0x7));
638 : else
639 0 : pabyDest[nPixelsOutput>>3] &= ~(1<<(nPixelsOutput & 0x7));
640 : }
641 81920 : else if( nDataType == EPT_u2 )
642 : {
643 0 : if( (nPixelsOutput & 0x3) == 0 )
644 0 : pabyDest[nPixelsOutput>>2] = (GByte) nDataValue;
645 0 : else if( (nPixelsOutput & 0x3) == 1 )
646 0 : pabyDest[nPixelsOutput>>2] |= (GByte) (nDataValue<<2);
647 0 : else if( (nPixelsOutput & 0x3) == 2 )
648 0 : pabyDest[nPixelsOutput>>2] |= (GByte) (nDataValue<<4);
649 : else
650 0 : pabyDest[nPixelsOutput>>2] |= (GByte) (nDataValue<<6);
651 : }
652 81920 : else if( nDataType == EPT_u4 )
653 : {
654 0 : if( (nPixelsOutput & 0x1) == 0 )
655 0 : pabyDest[nPixelsOutput>>1] = (GByte) nDataValue;
656 : else
657 0 : pabyDest[nPixelsOutput>>1] |= (GByte) (nDataValue<<4);
658 : }
659 81920 : else if( nDataType == EPT_s8 )
660 : {
661 0 : ((GByte *) pabyDest)[nPixelsOutput] = (GByte) nDataValue;
662 : }
663 81920 : else if( nDataType == EPT_u16 )
664 : {
665 0 : ((GUInt16 *) pabyDest)[nPixelsOutput] = (GUInt16) nDataValue;
666 : }
667 81920 : else if( nDataType == EPT_s16 )
668 : {
669 0 : ((GInt16 *) pabyDest)[nPixelsOutput] = (GInt16) nDataValue;
670 : }
671 81920 : else if( nDataType == EPT_s32 )
672 : {
673 16384 : ((GInt32 *) pabyDest)[nPixelsOutput] = nDataValue;
674 : }
675 65536 : else if( nDataType == EPT_u32 )
676 : {
677 0 : ((GUInt32 *) pabyDest)[nPixelsOutput] = nDataValue;
678 : }
679 65536 : else if( nDataType == EPT_f32 )
680 : {
681 : /* -------------------------------------------------------------------- */
682 : /* Note, floating point values are handled as if they were signed */
683 : /* 32-bit integers (bug #1000). */
684 : /* -------------------------------------------------------------------- */
685 65536 : ((float *) pabyDest)[nPixelsOutput] = *((float*)( &nDataValue ));
686 : }
687 : else
688 : {
689 0 : CPLAssert( FALSE );
690 : }
691 : }
692 :
693 20 : return CE_None;
694 : }
695 :
696 : /* ==================================================================== */
697 : /* Establish data pointers for runs. */
698 : /* ==================================================================== */
699 88 : if (nNumRuns < 0 || nDataOffset < 0)
700 : {
701 : CPLError(CE_Failure, CPLE_AppDefined, "nNumRuns=%d, nDataOffset=%d",
702 0 : nNumRuns, nDataOffset);
703 0 : return CE_Failure;
704 : }
705 :
706 88 : if (nNumBits > INT_MAX / nNumRuns ||
707 : nNumBits * nNumRuns > INT_MAX - 7 ||
708 : (nNumBits * nNumRuns + 7)/8 > INT_MAX - nDataOffset)
709 : {
710 : CPLError(CE_Failure, CPLE_AppDefined,
711 0 : "Integer overflow : nDataOffset + (nNumBits * nNumRuns + 7)/8");
712 0 : return CE_Failure;
713 : }
714 88 : CHECK_ENOUGH_BYTES(nDataOffset + (nNumBits * nNumRuns + 7)/8);
715 :
716 88 : pabyCounter = pabyCData + 13;
717 88 : nCounterOffset = 13;
718 88 : pabyValues = pabyCData + nDataOffset;
719 88 : nValueBitOffset = 0;
720 :
721 : /* -------------------------------------------------------------------- */
722 : /* Loop over runs. */
723 : /* -------------------------------------------------------------------- */
724 : int iRun;
725 :
726 62911 : for( iRun = 0; iRun < nNumRuns; iRun++ )
727 : {
728 62823 : int nRepeatCount = 0;
729 : int nDataValue;
730 :
731 : /* -------------------------------------------------------------------- */
732 : /* Get the repeat count. This can be stored as one, two, three */
733 : /* or four bytes depending on the low order two bits of the */
734 : /* first byte. */
735 : /* -------------------------------------------------------------------- */
736 62823 : CHECK_ENOUGH_BYTES(nCounterOffset+1);
737 62823 : if( ((*pabyCounter) & 0xc0) == 0x00 )
738 : {
739 62754 : nRepeatCount = (*(pabyCounter++)) & 0x3f;
740 62754 : nCounterOffset ++;
741 : }
742 69 : else if( ((*pabyCounter) & 0xc0) == 0x40 )
743 : {
744 69 : CHECK_ENOUGH_BYTES(nCounterOffset + 2);
745 69 : nRepeatCount = (*(pabyCounter++)) & 0x3f;
746 69 : nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
747 69 : nCounterOffset += 2;
748 : }
749 0 : else if( ((*pabyCounter) & 0xc0) == 0x80 )
750 : {
751 0 : CHECK_ENOUGH_BYTES(nCounterOffset + 3);
752 0 : nRepeatCount = (*(pabyCounter++)) & 0x3f;
753 0 : nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
754 0 : nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
755 0 : nCounterOffset += 3;
756 : }
757 0 : else if( ((*pabyCounter) & 0xc0) == 0xc0 )
758 : {
759 0 : CHECK_ENOUGH_BYTES(nCounterOffset + 4);
760 0 : nRepeatCount = (*(pabyCounter++)) & 0x3f;
761 0 : nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
762 0 : nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
763 0 : nRepeatCount = nRepeatCount * 256 + (*(pabyCounter++));
764 0 : nCounterOffset += 4;
765 : }
766 :
767 : /* -------------------------------------------------------------------- */
768 : /* Extract the data value in a way that depends on the number */
769 : /* of bits in it. */
770 : /* -------------------------------------------------------------------- */
771 62823 : if( nNumBits == 0 )
772 : {
773 0 : nDataValue = 0;
774 : }
775 62823 : else if( nNumBits == 1 )
776 : {
777 : nDataValue =
778 0 : (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0x1;
779 0 : nValueBitOffset++;
780 : }
781 62823 : else if( nNumBits == 2 )
782 : {
783 : nDataValue =
784 0 : (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0x3;
785 0 : nValueBitOffset += 2;
786 : }
787 62823 : else if( nNumBits == 4 )
788 : {
789 : nDataValue =
790 0 : (pabyValues[nValueBitOffset>>3] >> (nValueBitOffset&7)) & 0xf;
791 0 : nValueBitOffset += 4;
792 : }
793 62823 : else if( nNumBits == 8 )
794 : {
795 3795 : nDataValue = *pabyValues;
796 3795 : pabyValues++;
797 : }
798 59028 : else if( nNumBits == 16 )
799 : {
800 9864 : nDataValue = 256 * *(pabyValues++);
801 9864 : nDataValue += *(pabyValues++);
802 : }
803 49164 : else if( nNumBits == 32 )
804 : {
805 49164 : nDataValue = 256 * 256 * 256 * *(pabyValues++);
806 49164 : nDataValue += 256 * 256 * *(pabyValues++);
807 49164 : nDataValue += 256 * *(pabyValues++);
808 49164 : nDataValue += *(pabyValues++);
809 : }
810 : else
811 : {
812 : CPLError( CE_Failure, CPLE_NotSupported,
813 0 : "nNumBits = %d", nNumBits );
814 0 : return CE_Failure;
815 : }
816 :
817 : /* -------------------------------------------------------------------- */
818 : /* Offset by the minimum value. */
819 : /* -------------------------------------------------------------------- */
820 62823 : nDataValue += nDataMin;
821 :
822 : /* -------------------------------------------------------------------- */
823 : /* Now apply to the output buffer in a type specific way. */
824 : /* -------------------------------------------------------------------- */
825 62823 : if( nPixelsOutput + nRepeatCount > nMaxPixels )
826 : {
827 0 : CPLDebug("HFA", "Repeat count too big : %d", nRepeatCount);
828 0 : nRepeatCount = nMaxPixels - nPixelsOutput;
829 : }
830 :
831 62823 : if( nDataType == EPT_u8 )
832 : {
833 : int i;
834 :
835 5597 : for( i = 0; i < nRepeatCount; i++ )
836 : {
837 : //CPLAssert( nDataValue < 256 );
838 4864 : ((GByte *) pabyDest)[nPixelsOutput++] = (GByte)nDataValue;
839 : }
840 : }
841 62090 : else if( nDataType == EPT_u16 )
842 : {
843 : int i;
844 :
845 4458 : for( i = 0; i < nRepeatCount; i++ )
846 : {
847 4096 : ((GUInt16 *) pabyDest)[nPixelsOutput++] = (GUInt16)nDataValue;
848 : }
849 : }
850 61728 : else if( nDataType == EPT_s8 )
851 : {
852 : int i;
853 :
854 0 : for( i = 0; i < nRepeatCount; i++ )
855 : {
856 : //CPLAssert( nDataValue < 256 );
857 0 : ((GByte *) pabyDest)[nPixelsOutput++] = (GByte)nDataValue;
858 : }
859 : }
860 61728 : else if( nDataType == EPT_s16 )
861 : {
862 : int i;
863 :
864 0 : for( i = 0; i < nRepeatCount; i++ )
865 : {
866 0 : ((GInt16 *) pabyDest)[nPixelsOutput++] = (GInt16)nDataValue;
867 : }
868 : }
869 61728 : else if( nDataType == EPT_u32 )
870 : {
871 : int i;
872 :
873 0 : for( i = 0; i < nRepeatCount; i++ )
874 : {
875 0 : ((GUInt32 *) pabyDest)[nPixelsOutput++] = (GUInt32)nDataValue;
876 : }
877 : }
878 61728 : else if( nDataType == EPT_s32 )
879 : {
880 : int i;
881 :
882 58292 : for( i = 0; i < nRepeatCount; i++ )
883 : {
884 49152 : ((GInt32 *) pabyDest)[nPixelsOutput++] = (GInt32)nDataValue;
885 : }
886 : }
887 52588 : else if( nDataType == EPT_f32 )
888 : {
889 : int i;
890 : float fDataValue;
891 :
892 49164 : memcpy( &fDataValue, &nDataValue, 4);
893 245772 : for( i = 0; i < nRepeatCount; i++ )
894 : {
895 196608 : ((float *) pabyDest)[nPixelsOutput++] = fDataValue;
896 : }
897 : }
898 3424 : else if( nDataType == EPT_u1 )
899 : {
900 : int i;
901 :
902 : //CPLAssert( nDataValue == 0 || nDataValue == 1 );
903 :
904 2538 : if( nDataValue == 1 )
905 : {
906 69233 : for( i = 0; i < nRepeatCount; i++ )
907 : {
908 67962 : pabyDest[nPixelsOutput>>3] |= (1 << (nPixelsOutput & 0x7));
909 67962 : nPixelsOutput++;
910 : }
911 : }
912 : else
913 : {
914 23417 : for( i = 0; i < nRepeatCount; i++ )
915 : {
916 22150 : pabyDest[nPixelsOutput>>3] &= ~(1<<(nPixelsOutput & 0x7));
917 22150 : nPixelsOutput++;
918 : }
919 : }
920 : }
921 886 : else if( nDataType == EPT_u2 )
922 : {
923 : int i;
924 :
925 : //CPLAssert( nDataValue >= 0 && nDataValue < 4 );
926 :
927 13174 : for( i = 0; i < nRepeatCount; i++ )
928 : {
929 12288 : if( (nPixelsOutput & 0x3) == 0 )
930 3072 : pabyDest[nPixelsOutput>>2] = (GByte) nDataValue;
931 9216 : else if( (nPixelsOutput & 0x3) == 1 )
932 3072 : pabyDest[nPixelsOutput>>2] |= (GByte) (nDataValue<<2);
933 6144 : else if( (nPixelsOutput & 0x3) == 2 )
934 3072 : pabyDest[nPixelsOutput>>2] |= (GByte) (nDataValue<<4);
935 : else
936 3072 : pabyDest[nPixelsOutput>>2] |= (GByte) (nDataValue<<6);
937 12288 : nPixelsOutput++;
938 : }
939 : }
940 0 : else if( nDataType == EPT_u4 )
941 : {
942 : int i;
943 :
944 : //CPLAssert( nDataValue >= 0 && nDataValue < 16 );
945 :
946 0 : for( i = 0; i < nRepeatCount; i++ )
947 : {
948 0 : if( (nPixelsOutput & 0x1) == 0 )
949 0 : pabyDest[nPixelsOutput>>1] = (GByte) nDataValue;
950 : else
951 0 : pabyDest[nPixelsOutput>>1] |= (GByte) (nDataValue<<4);
952 :
953 0 : nPixelsOutput++;
954 : }
955 : }
956 : else
957 : {
958 : CPLError( CE_Failure, CPLE_AppDefined,
959 0 : "Attempt to uncompress an unsupported pixel data type.");
960 0 : return CE_Failure;
961 : }
962 : }
963 :
964 88 : return CE_None;
965 :
966 : not_enough_bytes:
967 :
968 0 : CPLError(CE_Failure, CPLE_AppDefined, "Not enough bytes in compressed block");
969 0 : return CE_Failure;
970 : }
971 :
972 : /************************************************************************/
973 : /* NullBlock() */
974 : /* */
975 : /* Set the block buffer to zero or the nodata value as */
976 : /* appropriate. */
977 : /************************************************************************/
978 :
979 159 : void HFABand::NullBlock( void *pData )
980 :
981 : {
982 159 : int nChunkSize = MAX(1,HFAGetDataTypeBits(nDataType)/8);
983 159 : int nWords = nBlockXSize * nBlockYSize;
984 :
985 159 : if( !bNoDataSet )
986 : {
987 : #ifdef ESRI_BUILD
988 : // We want special defaulting for 1 bit data in ArcGIS.
989 : if ( nDataType >= EPT_u2 )
990 : memset( pData, 0, nChunkSize*nWords );
991 : else
992 : memset( pData, 255, nChunkSize*nWords );
993 : #else
994 91 : memset( pData, 0, nChunkSize*nWords );
995 : #endif
996 : }
997 : else
998 : {
999 : double adfND[2];
1000 : int i;
1001 :
1002 68 : switch( nDataType )
1003 : {
1004 : case EPT_u1:
1005 : {
1006 0 : nWords = (nWords + 7)/8;
1007 0 : if( dfNoData != 0.0 )
1008 0 : ((unsigned char *) &adfND)[0] = 0xff;
1009 : else
1010 0 : ((unsigned char *) &adfND)[0] = 0x00;
1011 : }
1012 0 : break;
1013 :
1014 : case EPT_u2:
1015 : {
1016 0 : nWords = (nWords + 3)/4;
1017 0 : if( dfNoData == 0.0 )
1018 0 : ((unsigned char *) &adfND)[0] = 0x00;
1019 0 : else if( dfNoData == 1.0 )
1020 0 : ((unsigned char *) &adfND)[0] = 0x55;
1021 0 : else if( dfNoData == 2.0 )
1022 0 : ((unsigned char *) &adfND)[0] = 0xaa;
1023 : else
1024 0 : ((unsigned char *) &adfND)[0] = 0xff;
1025 : }
1026 0 : break;
1027 :
1028 : case EPT_u4:
1029 : {
1030 : unsigned char byVal =
1031 0 : (unsigned char) MAX(0,MIN(15,(int)dfNoData));
1032 :
1033 0 : nWords = (nWords + 1)/2;
1034 :
1035 0 : ((unsigned char *) &adfND)[0] = byVal + (byVal << 4);
1036 : }
1037 0 : break;
1038 :
1039 : case EPT_u8:
1040 65 : ((unsigned char *) &adfND)[0] =
1041 65 : (unsigned char) MAX(0,MIN(255,(int)dfNoData));
1042 65 : break;
1043 :
1044 : case EPT_s8:
1045 0 : ((signed char *) &adfND)[0] =
1046 0 : (signed char) MAX(-128,MIN(127,(int)dfNoData));
1047 0 : break;
1048 :
1049 : case EPT_u16:
1050 0 : ((GUInt16 *) &adfND)[0] = (GUInt16) dfNoData;
1051 0 : break;
1052 :
1053 : case EPT_s16:
1054 0 : ((GInt16 *) &adfND)[0] = (GInt16) dfNoData;
1055 0 : break;
1056 :
1057 : case EPT_u32:
1058 0 : ((GUInt32 *) &adfND)[0] = (GUInt32) dfNoData;
1059 0 : break;
1060 :
1061 : case EPT_s32:
1062 2 : ((GInt32 *) &adfND)[0] = (GInt32) dfNoData;
1063 2 : break;
1064 :
1065 : case EPT_f32:
1066 1 : ((float *) &adfND)[0] = (float) dfNoData;
1067 1 : break;
1068 :
1069 : case EPT_f64:
1070 0 : ((double *) &adfND)[0] = dfNoData;
1071 0 : break;
1072 :
1073 : case EPT_c64:
1074 0 : ((float *) &adfND)[0] = (float) dfNoData;
1075 0 : ((float *) &adfND)[1] = 0;
1076 0 : break;
1077 :
1078 : case EPT_c128:
1079 0 : ((double *) &adfND)[0] = dfNoData;
1080 0 : ((double *) &adfND)[1] = 0;
1081 : break;
1082 : }
1083 :
1084 278596 : for( i = 0; i < nWords; i++ )
1085 : memcpy( ((GByte *) pData) + nChunkSize * i,
1086 278528 : &adfND[0], nChunkSize );
1087 : }
1088 :
1089 159 : }
1090 :
1091 : /************************************************************************/
1092 : /* GetRasterBlock() */
1093 : /************************************************************************/
1094 :
1095 388 : CPLErr HFABand::GetRasterBlock( int nXBlock, int nYBlock, void * pData, int nDataSize )
1096 :
1097 : {
1098 : int iBlock;
1099 : VSILFILE *fpData;
1100 :
1101 388 : if( LoadBlockInfo() != CE_None )
1102 0 : return CE_Failure;
1103 :
1104 388 : iBlock = nXBlock + nYBlock * nBlocksPerRow;
1105 :
1106 : /* -------------------------------------------------------------------- */
1107 : /* If the block isn't valid, we just return all zeros, and an */
1108 : /* indication of success. */
1109 : /* -------------------------------------------------------------------- */
1110 388 : if( (panBlockFlag[iBlock] & BFLG_VALID) == 0 )
1111 : {
1112 159 : NullBlock( pData );
1113 159 : return( CE_None );
1114 : }
1115 :
1116 : /* -------------------------------------------------------------------- */
1117 : /* Otherwise we really read the data. */
1118 : /* -------------------------------------------------------------------- */
1119 : vsi_l_offset nBlockOffset;
1120 :
1121 : // Calculate block offset in case we have spill file. Use predefined
1122 : // block map otherwise.
1123 229 : if ( fpExternal )
1124 : {
1125 16 : fpData = fpExternal;
1126 : nBlockOffset = nBlockStart + nBlockSize * iBlock * nLayerStackCount
1127 16 : + nLayerStackIndex * nBlockSize;
1128 : }
1129 : else
1130 : {
1131 213 : fpData = psInfo->fp;
1132 213 : nBlockOffset = panBlockStart[iBlock];
1133 213 : nBlockSize = panBlockSize[iBlock];
1134 : }
1135 :
1136 229 : if( VSIFSeekL( fpData, nBlockOffset, SEEK_SET ) != 0 )
1137 : {
1138 : // XXX: We will not report error here, because file just may be
1139 : // in update state and data for this block will be available later
1140 0 : if ( psInfo->eAccess == HFA_Update )
1141 : {
1142 : memset( pData, 0,
1143 0 : HFAGetDataTypeBits(nDataType)*nBlockXSize*nBlockYSize/8 );
1144 0 : return CE_None;
1145 : }
1146 : else
1147 : {
1148 : CPLError( CE_Failure, CPLE_FileIO,
1149 : "Seek to %x:%08x on %p failed\n%s",
1150 : (int) (nBlockOffset >> 32),
1151 : (int) (nBlockOffset & 0xffffffff),
1152 0 : fpData, VSIStrerror(errno) );
1153 0 : return CE_Failure;
1154 : }
1155 : }
1156 :
1157 : /* -------------------------------------------------------------------- */
1158 : /* If the block is compressed, read into an intermediate buffer */
1159 : /* and convert. */
1160 : /* -------------------------------------------------------------------- */
1161 229 : if( panBlockFlag[iBlock] & BFLG_COMPRESSED )
1162 : {
1163 : GByte *pabyCData;
1164 : CPLErr eErr;
1165 :
1166 108 : pabyCData = (GByte *) VSIMalloc( (size_t) nBlockSize );
1167 108 : if (pabyCData == NULL)
1168 : {
1169 : CPLError( CE_Failure, CPLE_OutOfMemory,
1170 0 : "HFABand::GetRasterBlock : Out of memory\n");
1171 0 : return CE_Failure;
1172 : }
1173 :
1174 108 : if( VSIFReadL( pabyCData, (size_t) nBlockSize, 1, fpData ) != 1 )
1175 : {
1176 0 : CPLFree( pabyCData );
1177 :
1178 : // XXX: Suppose that file in update state
1179 0 : if ( psInfo->eAccess == HFA_Update )
1180 : {
1181 : memset( pData, 0,
1182 0 : HFAGetDataTypeBits(nDataType)*nBlockXSize*nBlockYSize/8 );
1183 0 : return CE_None;
1184 : }
1185 : else
1186 : {
1187 : CPLError( CE_Failure, CPLE_FileIO,
1188 : "Read of %d bytes at %x:%08x on %p failed.\n%s",
1189 : (int) nBlockSize,
1190 : (int) (nBlockOffset >> 32),
1191 : (int) (nBlockOffset & 0xffffffff),
1192 0 : fpData, VSIStrerror(errno) );
1193 0 : return CE_Failure;
1194 : }
1195 : }
1196 :
1197 : eErr = UncompressBlock( pabyCData, (int) nBlockSize,
1198 : (GByte *) pData, nBlockXSize*nBlockYSize,
1199 108 : nDataType );
1200 :
1201 108 : CPLFree( pabyCData );
1202 :
1203 108 : return eErr;
1204 : }
1205 :
1206 : /* -------------------------------------------------------------------- */
1207 : /* Read uncompressed data directly into the return buffer. */
1208 : /* -------------------------------------------------------------------- */
1209 121 : if ( nDataSize != -1 && (nBlockSize > INT_MAX ||
1210 : (int)nBlockSize > nDataSize) )
1211 : {
1212 : CPLError( CE_Failure, CPLE_AppDefined,
1213 0 : "Invalid block size : %d", (int)nBlockSize);
1214 0 : return CE_Failure;
1215 : }
1216 :
1217 121 : if( VSIFReadL( pData, (size_t) nBlockSize, 1, fpData ) != 1 )
1218 : {
1219 : memset( pData, 0,
1220 0 : HFAGetDataTypeBits(nDataType)*nBlockXSize*nBlockYSize/8 );
1221 :
1222 0 : if( fpData != fpExternal )
1223 : CPLDebug( "HFABand",
1224 : "Read of %x:%08x bytes at %d on %p failed.\n%s",
1225 : (int) nBlockSize,
1226 : (int) (nBlockOffset >> 32),
1227 : (int) (nBlockOffset & 0xffffffff),
1228 0 : fpData, VSIStrerror(errno) );
1229 :
1230 0 : return CE_None;
1231 : }
1232 :
1233 : /* -------------------------------------------------------------------- */
1234 : /* Byte swap to local byte order if required. It appears that */
1235 : /* raster data is always stored in Intel byte order in Imagine */
1236 : /* files. */
1237 : /* -------------------------------------------------------------------- */
1238 :
1239 : #ifdef CPL_MSB
1240 : if( HFAGetDataTypeBits(nDataType) == 16 )
1241 : {
1242 : for( int ii = 0; ii < nBlockXSize*nBlockYSize; ii++ )
1243 : CPL_SWAP16PTR( ((unsigned char *) pData) + ii*2 );
1244 : }
1245 : else if( HFAGetDataTypeBits(nDataType) == 32 )
1246 : {
1247 : for( int ii = 0; ii < nBlockXSize*nBlockYSize; ii++ )
1248 : CPL_SWAP32PTR( ((unsigned char *) pData) + ii*4 );
1249 : }
1250 : else if( nDataType == EPT_f64 )
1251 : {
1252 : for( int ii = 0; ii < nBlockXSize*nBlockYSize; ii++ )
1253 : CPL_SWAP64PTR( ((unsigned char *) pData) + ii*8 );
1254 : }
1255 : else if( nDataType == EPT_c64 )
1256 : {
1257 : for( int ii = 0; ii < nBlockXSize*nBlockYSize*2; ii++ )
1258 : CPL_SWAP32PTR( ((unsigned char *) pData) + ii*4 );
1259 : }
1260 : else if( nDataType == EPT_c128 )
1261 : {
1262 : for( int ii = 0; ii < nBlockXSize*nBlockYSize*2; ii++ )
1263 : CPL_SWAP64PTR( ((unsigned char *) pData) + ii*8 );
1264 : }
1265 : #endif /* def CPL_MSB */
1266 :
1267 121 : return( CE_None );
1268 : }
1269 :
1270 : /************************************************************************/
1271 : /* ReAllocBlock() */
1272 : /************************************************************************/
1273 :
1274 5 : void HFABand::ReAllocBlock( int iBlock, int nSize )
1275 : {
1276 : /* For compressed files - need to realloc the space for the block */
1277 :
1278 : // TODO: Should check to see if panBlockStart[iBlock] is not zero then do a HFAFreeSpace()
1279 : // but that doesn't exist yet.
1280 : // Instead as in interim measure it will reuse the existing block if
1281 : // the new data will fit in.
1282 5 : if( ( panBlockStart[iBlock] != 0 ) && ( nSize <= panBlockSize[iBlock] ) )
1283 : {
1284 0 : panBlockSize[iBlock] = nSize;
1285 : //fprintf( stderr, "Reusing block %d\n", iBlock );
1286 : }
1287 : else
1288 : {
1289 5 : panBlockStart[iBlock] = HFAAllocateSpace( psInfo, nSize );
1290 :
1291 5 : panBlockSize[iBlock] = nSize;
1292 :
1293 : // need to re - write this info to the RasterDMS node
1294 5 : HFAEntry *poDMS = poNode->GetNamedChild( "RasterDMS" );
1295 :
1296 : char szVarName[64];
1297 5 : sprintf( szVarName, "blockinfo[%d].offset", iBlock );
1298 5 : poDMS->SetIntField( szVarName, (int) panBlockStart[iBlock] );
1299 :
1300 5 : sprintf( szVarName, "blockinfo[%d].size", iBlock );
1301 5 : poDMS->SetIntField( szVarName, panBlockSize[iBlock] );
1302 : }
1303 :
1304 5 : }
1305 :
1306 :
1307 : /************************************************************************/
1308 : /* SetRasterBlock() */
1309 : /************************************************************************/
1310 :
1311 117 : CPLErr HFABand::SetRasterBlock( int nXBlock, int nYBlock, void * pData )
1312 :
1313 : {
1314 : int iBlock;
1315 : VSILFILE *fpData;
1316 :
1317 117 : if( psInfo->eAccess == HFA_ReadOnly )
1318 : {
1319 : CPLError( CE_Failure, CPLE_NoWriteAccess,
1320 0 : "Attempt to write block to read-only HFA file failed." );
1321 0 : return CE_Failure;
1322 : }
1323 :
1324 117 : if( LoadBlockInfo() != CE_None )
1325 0 : return CE_Failure;
1326 :
1327 117 : iBlock = nXBlock + nYBlock * nBlocksPerRow;
1328 :
1329 : /* -------------------------------------------------------------------- */
1330 : /* For now we don't support write invalid uncompressed blocks. */
1331 : /* To do so we will need logic to make space at the end of the */
1332 : /* file in the right size. */
1333 : /* -------------------------------------------------------------------- */
1334 316 : if( (panBlockFlag[iBlock] & BFLG_VALID) == 0
1335 102 : && !(panBlockFlag[iBlock] & BFLG_COMPRESSED)
1336 97 : && panBlockStart[iBlock] == 0 )
1337 : {
1338 : CPLError( CE_Failure, CPLE_AppDefined,
1339 : "Attempt to write to invalid tile with number %d "
1340 : "(X position %d, Y position %d). This\n operation currently "
1341 : "unsupported by HFABand::SetRasterBlock().\n",
1342 0 : iBlock, nXBlock, nYBlock );
1343 :
1344 0 : return CE_Failure;
1345 : }
1346 :
1347 : /* -------------------------------------------------------------------- */
1348 : /* Move to the location that the data sits. */
1349 : /* -------------------------------------------------------------------- */
1350 : vsi_l_offset nBlockOffset;
1351 :
1352 : // Calculate block offset in case we have spill file. Use predefined
1353 : // block map otherwise.
1354 117 : if ( fpExternal )
1355 : {
1356 7 : fpData = fpExternal;
1357 : nBlockOffset = nBlockStart + nBlockSize * iBlock * nLayerStackCount
1358 7 : + nLayerStackIndex * nBlockSize;
1359 : }
1360 : else
1361 : {
1362 110 : fpData = psInfo->fp;
1363 110 : nBlockOffset = panBlockStart[iBlock];
1364 110 : nBlockSize = panBlockSize[iBlock];
1365 : }
1366 :
1367 : /* ==================================================================== */
1368 : /* Compressed Tile Handling. */
1369 : /* ==================================================================== */
1370 117 : if( panBlockFlag[iBlock] & BFLG_COMPRESSED )
1371 : {
1372 : /* ------------------------------------------------------------ */
1373 : /* Write compressed data. */
1374 : /* ------------------------------------------------------------ */
1375 5 : int nInBlockSize = (nBlockXSize * nBlockYSize * HFAGetDataTypeBits(nDataType) + 7 ) / 8;
1376 :
1377 : /* create the compressor object */
1378 5 : HFACompress compress( pData, nInBlockSize, nDataType );
1379 5 : if( compress.getCounts() == NULL ||
1380 : compress.getValues() == NULL)
1381 : {
1382 0 : CPLError( CE_Failure, CPLE_OutOfMemory, "Out of memory");
1383 0 : return CE_Failure;
1384 : }
1385 :
1386 : /* compress the data */
1387 5 : if( compress.compressBlock() )
1388 : {
1389 : /* get the data out of the object */
1390 3 : GByte *pCounts = compress.getCounts();
1391 3 : GUInt32 nSizeCount = compress.getCountSize();
1392 3 : GByte *pValues = compress.getValues();
1393 3 : GUInt32 nSizeValues = compress.getValueSize();
1394 3 : GUInt32 nMin = compress.getMin();
1395 3 : GUInt32 nNumRuns = compress.getNumRuns();
1396 3 : GByte nNumBits = compress.getNumBits();
1397 :
1398 : /* Compensate for the header info */
1399 3 : GUInt32 nDataOffset = nSizeCount + 13;
1400 3 : int nTotalSize = nSizeCount + nSizeValues + 13;
1401 :
1402 : //fprintf( stderr, "sizecount = %d sizevalues = %d min = %d numruns = %d numbits = %d\n", nSizeCount, nSizeValues, nMin, nNumRuns, (int)nNumBits );
1403 :
1404 : // Allocate space for the compressed block and seek to it.
1405 3 : ReAllocBlock( iBlock, nTotalSize );
1406 :
1407 3 : nBlockOffset = panBlockStart[iBlock];
1408 3 : nBlockSize = panBlockSize[iBlock];
1409 :
1410 : // Seek to offset
1411 3 : if( VSIFSeekL( fpData, nBlockOffset, SEEK_SET ) != 0 )
1412 : {
1413 : CPLError( CE_Failure, CPLE_FileIO, "Seek to %x:%08x on %p failed\n%s",
1414 : (int) (nBlockOffset >> 32),
1415 : (int) (nBlockOffset & 0xffffffff),
1416 0 : fpData, VSIStrerror(errno) );
1417 0 : return CE_Failure;
1418 : }
1419 :
1420 : /* -------------------------------------------------------------------- */
1421 : /* Byte swap to local byte order if required. It appears that */
1422 : /* raster data is always stored in Intel byte order in Imagine */
1423 : /* files. */
1424 : /* -------------------------------------------------------------------- */
1425 :
1426 : #ifdef CPL_MSB
1427 :
1428 : CPL_SWAP32PTR( &nMin );
1429 : CPL_SWAP32PTR( &nNumRuns );
1430 : CPL_SWAP32PTR( &nDataOffset );
1431 :
1432 : #endif /* def CPL_MSB */
1433 :
1434 : /* Write out the Minimum value */
1435 3 : VSIFWriteL( &nMin, (size_t) sizeof( nMin ), 1, fpData );
1436 :
1437 : /* the number of runs */
1438 3 : VSIFWriteL( &nNumRuns, (size_t) sizeof( nNumRuns ), 1, fpData );
1439 :
1440 : /* The offset to the data */
1441 3 : VSIFWriteL( &nDataOffset, (size_t) sizeof( nDataOffset ), 1, fpData );
1442 :
1443 : /* The number of bits */
1444 3 : VSIFWriteL( &nNumBits, (size_t) sizeof( nNumBits ), 1, fpData );
1445 :
1446 : /* The counters - MSB stuff handled in HFACompress */
1447 3 : VSIFWriteL( pCounts, (size_t) sizeof( GByte ), nSizeCount, fpData );
1448 :
1449 : /* The values - MSB stuff handled in HFACompress */
1450 3 : VSIFWriteL( pValues, (size_t) sizeof( GByte ), nSizeValues, fpData );
1451 :
1452 : /* Compressed data is freed in the HFACompress destructor */
1453 : }
1454 : else
1455 : {
1456 : /* If we have actually made the block bigger - ie does not compress well */
1457 2 : panBlockFlag[iBlock] ^= BFLG_COMPRESSED;
1458 : // alloc more space for the uncompressed block
1459 2 : ReAllocBlock( iBlock, nInBlockSize );
1460 :
1461 2 : nBlockOffset = panBlockStart[iBlock];
1462 2 : nBlockSize = panBlockSize[iBlock];
1463 :
1464 : /* Need to change the RasterDMS entry */
1465 2 : HFAEntry *poDMS = poNode->GetNamedChild( "RasterDMS" );
1466 :
1467 : char szVarName[64];
1468 2 : sprintf( szVarName, "blockinfo[%d].compressionType", iBlock );
1469 2 : poDMS->SetIntField( szVarName, 0 );
1470 : }
1471 :
1472 : /* -------------------------------------------------------------------- */
1473 : /* If the block was previously invalid, mark it as valid now. */
1474 : /* -------------------------------------------------------------------- */
1475 5 : if( (panBlockFlag[iBlock] & BFLG_VALID) == 0 )
1476 : {
1477 : char szVarName[64];
1478 5 : HFAEntry *poDMS = poNode->GetNamedChild( "RasterDMS" );
1479 :
1480 5 : sprintf( szVarName, "blockinfo[%d].logvalid", iBlock );
1481 5 : poDMS->SetStringField( szVarName, "true" );
1482 :
1483 5 : panBlockFlag[iBlock] |= BFLG_VALID;
1484 0 : }
1485 : }
1486 :
1487 : /* ==================================================================== */
1488 : /* Uncompressed TILE handling. */
1489 : /* ==================================================================== */
1490 117 : if( ( panBlockFlag[iBlock] & BFLG_COMPRESSED ) == 0 )
1491 : {
1492 :
1493 114 : if( VSIFSeekL( fpData, nBlockOffset, SEEK_SET ) != 0 )
1494 : {
1495 : CPLError( CE_Failure, CPLE_FileIO, "Seek to %x:%08x on %p failed\n%s",
1496 : (int) (nBlockOffset >> 32),
1497 : (int) (nBlockOffset & 0xffffffff),
1498 0 : fpData, VSIStrerror(errno) );
1499 0 : return CE_Failure;
1500 : }
1501 :
1502 : /* -------------------------------------------------------------------- */
1503 : /* Byte swap to local byte order if required. It appears that */
1504 : /* raster data is always stored in Intel byte order in Imagine */
1505 : /* files. */
1506 : /* -------------------------------------------------------------------- */
1507 :
1508 : #ifdef CPL_MSB
1509 : if( HFAGetDataTypeBits(nDataType) == 16 )
1510 : {
1511 : for( int ii = 0; ii < nBlockXSize*nBlockYSize; ii++ )
1512 : CPL_SWAP16PTR( ((unsigned char *) pData) + ii*2 );
1513 : }
1514 : else if( HFAGetDataTypeBits(nDataType) == 32 )
1515 : {
1516 : for( int ii = 0; ii < nBlockXSize*nBlockYSize; ii++ )
1517 : CPL_SWAP32PTR( ((unsigned char *) pData) + ii*4 );
1518 : }
1519 : else if( nDataType == EPT_f64 )
1520 : {
1521 : for( int ii = 0; ii < nBlockXSize*nBlockYSize; ii++ )
1522 : CPL_SWAP64PTR( ((unsigned char *) pData) + ii*8 );
1523 : }
1524 : else if( nDataType == EPT_c64 )
1525 : {
1526 : for( int ii = 0; ii < nBlockXSize*nBlockYSize*2; ii++ )
1527 : CPL_SWAP32PTR( ((unsigned char *) pData) + ii*4 );
1528 : }
1529 : else if( nDataType == EPT_c128 )
1530 : {
1531 : for( int ii = 0; ii < nBlockXSize*nBlockYSize*2; ii++ )
1532 : CPL_SWAP64PTR( ((unsigned char *) pData) + ii*8 );
1533 : }
1534 : #endif /* def CPL_MSB */
1535 :
1536 : /* -------------------------------------------------------------------- */
1537 : /* Write uncompressed data. */
1538 : /* -------------------------------------------------------------------- */
1539 114 : if( VSIFWriteL( pData, (size_t) nBlockSize, 1, fpData ) != 1 )
1540 : {
1541 : CPLError( CE_Failure, CPLE_FileIO,
1542 : "Write of %d bytes at %x:%08x on %p failed.\n%s",
1543 : (int) nBlockSize,
1544 : (int) (nBlockOffset >> 32),
1545 : (int) (nBlockOffset & 0xffffffff),
1546 0 : fpData, VSIStrerror(errno) );
1547 0 : return CE_Failure;
1548 : }
1549 :
1550 : /* -------------------------------------------------------------------- */
1551 : /* If the block was previously invalid, mark it as valid now. */
1552 : /* -------------------------------------------------------------------- */
1553 114 : if( (panBlockFlag[iBlock] & BFLG_VALID) == 0 )
1554 : {
1555 : char szVarName[64];
1556 97 : HFAEntry *poDMS = poNode->GetNamedChild( "RasterDMS" );
1557 :
1558 97 : sprintf( szVarName, "blockinfo[%d].logvalid", iBlock );
1559 97 : poDMS->SetStringField( szVarName, "true" );
1560 :
1561 97 : panBlockFlag[iBlock] |= BFLG_VALID;
1562 : }
1563 : }
1564 : /* -------------------------------------------------------------------- */
1565 : /* Swap back, since we don't really have permission to change */
1566 : /* the callers buffer. */
1567 : /* -------------------------------------------------------------------- */
1568 :
1569 : #ifdef CPL_MSB
1570 : if( HFAGetDataTypeBits(nDataType) == 16 )
1571 : {
1572 : for( int ii = 0; ii < nBlockXSize*nBlockYSize; ii++ )
1573 : CPL_SWAP16PTR( ((unsigned char *) pData) + ii*2 );
1574 : }
1575 : else if( HFAGetDataTypeBits(nDataType) == 32 )
1576 : {
1577 : for( int ii = 0; ii < nBlockXSize*nBlockYSize; ii++ )
1578 : CPL_SWAP32PTR( ((unsigned char *) pData) + ii*4 );
1579 : }
1580 : else if( nDataType == EPT_f64 )
1581 : {
1582 : for( int ii = 0; ii < nBlockXSize*nBlockYSize; ii++ )
1583 : CPL_SWAP64PTR( ((unsigned char *) pData) + ii*8 );
1584 : }
1585 : else if( nDataType == EPT_c64 )
1586 : {
1587 : for( int ii = 0; ii < nBlockXSize*nBlockYSize*2; ii++ )
1588 : CPL_SWAP32PTR( ((unsigned char *) pData) + ii*4 );
1589 : }
1590 : else if( nDataType == EPT_c128 )
1591 : {
1592 : for( int ii = 0; ii < nBlockXSize*nBlockYSize*2; ii++ )
1593 : CPL_SWAP64PTR( ((unsigned char *) pData) + ii*8 );
1594 : }
1595 : #endif /* def CPL_MSB */
1596 :
1597 117 : return( CE_None );
1598 : }
1599 :
1600 : /************************************************************************/
1601 : /* GetBandName() */
1602 : /* */
1603 : /* Return the Layer Name */
1604 : /************************************************************************/
1605 :
1606 179 : const char * HFABand::GetBandName()
1607 : {
1608 179 : if( strlen(poNode->GetName()) > 0 )
1609 179 : return( poNode->GetName() );
1610 : else
1611 : {
1612 : int iBand;
1613 0 : for( iBand = 0; iBand < psInfo->nBands; iBand++ )
1614 : {
1615 0 : if( psInfo->papoBand[iBand] == this )
1616 : {
1617 0 : osOverName.Printf( "Layer_%d", iBand+1 );
1618 0 : return osOverName;
1619 : }
1620 : }
1621 :
1622 0 : osOverName.Printf( "Layer_%x", poNode->GetFilePos() );
1623 0 : return osOverName;
1624 : }
1625 : }
1626 :
1627 : /************************************************************************/
1628 : /* SetBandName() */
1629 : /* */
1630 : /* Set the Layer Name */
1631 : /************************************************************************/
1632 :
1633 5 : void HFABand::SetBandName(const char *pszName)
1634 : {
1635 5 : if( psInfo->eAccess == HFA_Update )
1636 : {
1637 5 : poNode->SetName(pszName);
1638 : }
1639 5 : }
1640 :
1641 : /************************************************************************/
1642 : /* SetNoDataValue() */
1643 : /* */
1644 : /* Set the band no-data value */
1645 : /************************************************************************/
1646 :
1647 8 : CPLErr HFABand::SetNoDataValue( double dfValue )
1648 : {
1649 8 : CPLErr eErr = CE_Failure;
1650 :
1651 8 : if ( psInfo->eAccess == HFA_Update )
1652 : {
1653 8 : HFAEntry *poNDNode = poNode->GetNamedChild( "Eimg_NonInitializedValue" );
1654 :
1655 8 : if ( poNDNode == NULL )
1656 : {
1657 : poNDNode = new HFAEntry( psInfo,
1658 : "Eimg_NonInitializedValue",
1659 : "Eimg_NonInitializedValue",
1660 8 : poNode );
1661 : }
1662 :
1663 8 : poNDNode->MakeData( 8 + 12 + 8 );
1664 8 : poNDNode->SetPosition();
1665 :
1666 8 : poNDNode->SetIntField( "valueBD[-3]", EPT_f64 );
1667 8 : poNDNode->SetIntField( "valueBD[-2]", 1 );
1668 8 : poNDNode->SetIntField( "valueBD[-1]", 1 );
1669 8 : if ( poNDNode->SetDoubleField( "valueBD[0]", dfValue) != CE_Failure )
1670 : {
1671 8 : bNoDataSet = TRUE;
1672 8 : dfNoData = dfValue;
1673 8 : eErr = CE_None;
1674 : }
1675 : }
1676 :
1677 8 : return eErr;
1678 : }
1679 :
1680 : /************************************************************************/
1681 : /* HFAReadBFUniqueBins() */
1682 : /* */
1683 : /* Attempt to read the bins used for a PCT or RAT from a */
1684 : /* BinFunction node. On failure just return NULL. */
1685 : /************************************************************************/
1686 :
1687 24 : double *HFAReadBFUniqueBins( HFAEntry *poBinFunc, int nPCTColors )
1688 :
1689 : {
1690 : /* -------------------------------------------------------------------- */
1691 : /* First confirm this is a "BFUnique" bin function. We don't */
1692 : /* know what to do with any other types. */
1693 : /* -------------------------------------------------------------------- */
1694 : const char *pszBinFunctionType =
1695 24 : poBinFunc->GetStringField( "binFunction.type.string" );
1696 :
1697 24 : if( pszBinFunctionType == NULL
1698 : || !EQUAL(pszBinFunctionType,"BFUnique") )
1699 0 : return NULL;
1700 :
1701 : /* -------------------------------------------------------------------- */
1702 : /* Process dictionary. */
1703 : /* -------------------------------------------------------------------- */
1704 : const char *pszDict =
1705 24 : poBinFunc->GetStringField( "binFunction.MIFDictionary.string" );
1706 24 : if( pszDict == NULL )
1707 0 : poBinFunc->GetStringField( "binFunction.MIFDictionary" );
1708 :
1709 24 : HFADictionary oMiniDict( pszDict );
1710 :
1711 24 : HFAType *poBFUnique = oMiniDict.FindType( "BFUnique" );
1712 24 : if( poBFUnique == NULL )
1713 0 : return NULL;
1714 :
1715 : /* -------------------------------------------------------------------- */
1716 : /* Field the MIFObject raw data pointer. */
1717 : /* -------------------------------------------------------------------- */
1718 : const GByte *pabyMIFObject = (const GByte *)
1719 24 : poBinFunc->GetStringField("binFunction.MIFObject");
1720 :
1721 24 : if( pabyMIFObject == NULL )
1722 0 : return NULL;
1723 :
1724 : /* -------------------------------------------------------------------- */
1725 : /* Confirm that this is a 64bit floating point basearray. */
1726 : /* -------------------------------------------------------------------- */
1727 24 : if( pabyMIFObject[20] != 0x0a || pabyMIFObject[21] != 0x00 )
1728 : {
1729 0 : CPLDebug( "HFA", "HFAReadPCTBins(): The basedata does not appear to be EGDA_TYPE_F64." );
1730 0 : return NULL;
1731 : }
1732 :
1733 : /* -------------------------------------------------------------------- */
1734 : /* Decode bins. */
1735 : /* -------------------------------------------------------------------- */
1736 24 : double *padfBins = (double *) CPLCalloc(sizeof(double),nPCTColors);
1737 : int i;
1738 :
1739 24 : memcpy( padfBins, pabyMIFObject + 24, sizeof(double) * nPCTColors );
1740 :
1741 24 : for( i = 0; i < nPCTColors; i++ )
1742 : {
1743 : HFAStandard( 8, padfBins + i );
1744 : // CPLDebug( "HFA", "Bin[%d] = %g", i, padfBins[i] );
1745 : }
1746 :
1747 24 : return padfBins;
1748 : }
1749 :
1750 : /************************************************************************/
1751 : /* GetPCT() */
1752 : /* */
1753 : /* Return PCT information, if any exists. */
1754 : /************************************************************************/
1755 :
1756 466 : CPLErr HFABand::GetPCT( int * pnColors,
1757 : double **ppadfRed,
1758 : double **ppadfGreen,
1759 : double **ppadfBlue,
1760 : double **ppadfAlpha,
1761 : double **ppadfBins )
1762 :
1763 : {
1764 466 : *pnColors = 0;
1765 466 : *ppadfRed = NULL;
1766 466 : *ppadfGreen = NULL;
1767 466 : *ppadfBlue = NULL;
1768 466 : *ppadfAlpha = NULL;
1769 466 : *ppadfBins = NULL;
1770 :
1771 : /* -------------------------------------------------------------------- */
1772 : /* If we haven't already tried to load the colors, do so now. */
1773 : /* -------------------------------------------------------------------- */
1774 466 : if( nPCTColors == -1 )
1775 : {
1776 : HFAEntry *poColumnEntry;
1777 : int i, iColumn;
1778 :
1779 466 : nPCTColors = 0;
1780 :
1781 466 : poColumnEntry = poNode->GetNamedChild("Descriptor_Table.Red");
1782 466 : if( poColumnEntry == NULL )
1783 459 : return( CE_Failure );
1784 :
1785 : /* FIXME? : we could also check that nPCTColors is not too big */
1786 7 : nPCTColors = poColumnEntry->GetIntField( "numRows" );
1787 35 : for( iColumn = 0; iColumn < 4; iColumn++ )
1788 : {
1789 28 : apadfPCT[iColumn] = (double *)VSIMalloc2(sizeof(double),nPCTColors);
1790 28 : if (apadfPCT[iColumn] == NULL)
1791 : {
1792 0 : CPLError(CE_Failure, CPLE_OutOfMemory, "Color palette will be ignored");
1793 0 : return CE_Failure;
1794 : }
1795 :
1796 28 : if( iColumn == 0 )
1797 7 : poColumnEntry = poNode->GetNamedChild("Descriptor_Table.Red");
1798 21 : else if( iColumn == 1 )
1799 7 : poColumnEntry= poNode->GetNamedChild("Descriptor_Table.Green");
1800 14 : else if( iColumn == 2 )
1801 7 : poColumnEntry = poNode->GetNamedChild("Descriptor_Table.Blue");
1802 7 : else if( iColumn == 3 ) {
1803 7 : poColumnEntry = poNode->GetNamedChild("Descriptor_Table.Opacity");
1804 : }
1805 :
1806 28 : if( poColumnEntry == NULL )
1807 : {
1808 0 : double *pdCol = apadfPCT[iColumn];
1809 0 : for( i = 0; i < nPCTColors; i++ )
1810 0 : pdCol[i] = 1.0;
1811 : }
1812 : else
1813 : {
1814 28 : if (VSIFSeekL( psInfo->fp, poColumnEntry->GetIntField("columnDataPtr"),
1815 : SEEK_SET ) < 0)
1816 : {
1817 : CPLError( CE_Failure, CPLE_FileIO,
1818 0 : "VSIFSeekL() failed in HFABand::GetPCT()." );
1819 0 : return CE_Failure;
1820 : }
1821 28 : if (VSIFReadL( apadfPCT[iColumn], sizeof(double), nPCTColors,
1822 : psInfo->fp) != (size_t)nPCTColors)
1823 : {
1824 : CPLError( CE_Failure, CPLE_FileIO,
1825 0 : "VSIFReadL() failed in HFABand::GetPCT()." );
1826 0 : return CE_Failure;
1827 : }
1828 :
1829 28 : for( i = 0; i < nPCTColors; i++ )
1830 : HFAStandard( 8, apadfPCT[iColumn] + i );
1831 : }
1832 : }
1833 :
1834 : /* -------------------------------------------------------------------- */
1835 : /* Do we have a custom binning function? If so, try reading it. */
1836 : /* -------------------------------------------------------------------- */
1837 : HFAEntry *poBinFunc =
1838 7 : poNode->GetNamedChild("Descriptor_Table.#Bin_Function840#");
1839 :
1840 7 : if( poBinFunc != NULL )
1841 : {
1842 4 : padfPCTBins = HFAReadBFUniqueBins( poBinFunc, nPCTColors );
1843 : }
1844 : }
1845 :
1846 : /* -------------------------------------------------------------------- */
1847 : /* Return the values. */
1848 : /* -------------------------------------------------------------------- */
1849 7 : if( nPCTColors == 0 )
1850 0 : return( CE_Failure );
1851 :
1852 7 : *pnColors = nPCTColors;
1853 7 : *ppadfRed = apadfPCT[0];
1854 7 : *ppadfGreen = apadfPCT[1];
1855 7 : *ppadfBlue = apadfPCT[2];
1856 7 : *ppadfAlpha = apadfPCT[3];
1857 7 : *ppadfBins = padfPCTBins;
1858 :
1859 7 : return( CE_None );
1860 : }
1861 :
1862 : /************************************************************************/
1863 : /* SetPCT() */
1864 : /* */
1865 : /* Set the PCT information for this band. */
1866 : /************************************************************************/
1867 :
1868 3 : CPLErr HFABand::SetPCT( int nColors,
1869 : double *padfRed,
1870 : double *padfGreen,
1871 : double *padfBlue ,
1872 : double *padfAlpha)
1873 :
1874 : {
1875 : static const char *apszColNames[4] = {"Red", "Green", "Blue", "Opacity"};
1876 : HFAEntry *poEdsc_Table;
1877 : int iColumn;
1878 :
1879 : /* -------------------------------------------------------------------- */
1880 : /* Do we need to try and clear any existing color table? */
1881 : /* -------------------------------------------------------------------- */
1882 3 : if( nColors == 0 )
1883 : {
1884 2 : poEdsc_Table = poNode->GetNamedChild( "Descriptor_Table" );
1885 2 : if( poEdsc_Table == NULL )
1886 0 : return CE_None;
1887 :
1888 10 : for( iColumn = 0; iColumn < 4; iColumn++ )
1889 : {
1890 : HFAEntry *poEdsc_Column;
1891 :
1892 8 : poEdsc_Column = poEdsc_Table->GetNamedChild(apszColNames[iColumn]);
1893 8 : if( poEdsc_Column )
1894 8 : poEdsc_Column->RemoveAndDestroy();
1895 : }
1896 :
1897 2 : return CE_None;
1898 : }
1899 :
1900 : /* -------------------------------------------------------------------- */
1901 : /* Create the Descriptor table. */
1902 : /* -------------------------------------------------------------------- */
1903 1 : poEdsc_Table = poNode->GetNamedChild( "Descriptor_Table" );
1904 1 : if( poEdsc_Table == NULL
1905 : || !EQUAL(poEdsc_Table->GetType(),"Edsc_Table") )
1906 : poEdsc_Table = new HFAEntry( psInfo, "Descriptor_Table",
1907 1 : "Edsc_Table", poNode );
1908 :
1909 1 : poEdsc_Table->SetIntField( "numrows", nColors );
1910 :
1911 : /* -------------------------------------------------------------------- */
1912 : /* Create the Binning function node. I am not sure that we */
1913 : /* really need this though. */
1914 : /* -------------------------------------------------------------------- */
1915 : HFAEntry *poEdsc_BinFunction;
1916 :
1917 1 : poEdsc_BinFunction = poEdsc_Table->GetNamedChild( "#Bin_Function#" );
1918 1 : if( poEdsc_BinFunction == NULL
1919 : || !EQUAL(poEdsc_BinFunction->GetType(),"Edsc_BinFunction") )
1920 : poEdsc_BinFunction = new HFAEntry( psInfo, "#Bin_Function#",
1921 : "Edsc_BinFunction",
1922 1 : poEdsc_Table );
1923 :
1924 : // Because of the BaseData we have to hardcode the size.
1925 1 : poEdsc_BinFunction->MakeData( 30 );
1926 :
1927 1 : poEdsc_BinFunction->SetIntField( "numBins", nColors );
1928 1 : poEdsc_BinFunction->SetStringField( "binFunction", "direct" );
1929 1 : poEdsc_BinFunction->SetDoubleField( "minLimit", 0.0 );
1930 1 : poEdsc_BinFunction->SetDoubleField( "maxLimit", nColors - 1.0 );
1931 :
1932 : /* -------------------------------------------------------------------- */
1933 : /* Process each color component */
1934 : /* -------------------------------------------------------------------- */
1935 5 : for( iColumn = 0; iColumn < 4; iColumn++ )
1936 : {
1937 : HFAEntry *poEdsc_Column;
1938 4 : double *padfValues=NULL;
1939 4 : const char *pszName = apszColNames[iColumn];
1940 :
1941 4 : if( iColumn == 0 )
1942 1 : padfValues = padfRed;
1943 3 : else if( iColumn == 1 )
1944 1 : padfValues = padfGreen;
1945 2 : else if( iColumn == 2 )
1946 1 : padfValues = padfBlue;
1947 1 : else if( iColumn == 3 )
1948 1 : padfValues = padfAlpha;
1949 :
1950 : /* -------------------------------------------------------------------- */
1951 : /* Create the Edsc_Column. */
1952 : /* -------------------------------------------------------------------- */
1953 4 : poEdsc_Column = poEdsc_Table->GetNamedChild( pszName );
1954 4 : if( poEdsc_Column == NULL
1955 : || !EQUAL(poEdsc_Column->GetType(),"Edsc_Column") )
1956 : poEdsc_Column = new HFAEntry( psInfo, pszName, "Edsc_Column",
1957 4 : poEdsc_Table );
1958 :
1959 4 : poEdsc_Column->SetIntField( "numRows", nColors );
1960 4 : poEdsc_Column->SetStringField( "dataType", "real" );
1961 4 : poEdsc_Column->SetIntField( "maxNumChars", 0 );
1962 :
1963 : /* -------------------------------------------------------------------- */
1964 : /* Write the data out. */
1965 : /* -------------------------------------------------------------------- */
1966 4 : int nOffset = HFAAllocateSpace( psInfo, 8*nColors);
1967 : double *padfFileData;
1968 :
1969 4 : poEdsc_Column->SetIntField( "columnDataPtr", nOffset );
1970 :
1971 4 : padfFileData = (double *) CPLMalloc(nColors*sizeof(double));
1972 1028 : for( int iColor = 0; iColor < nColors; iColor++ )
1973 : {
1974 1024 : padfFileData[iColor] = padfValues[iColor];
1975 : HFAStandard( 8, padfFileData + iColor );
1976 : }
1977 4 : VSIFSeekL( psInfo->fp, nOffset, SEEK_SET );
1978 4 : VSIFWriteL( padfFileData, 8, nColors, psInfo->fp );
1979 4 : CPLFree( padfFileData );
1980 : }
1981 :
1982 : /* -------------------------------------------------------------------- */
1983 : /* Update the layer type to be thematic. */
1984 : /* -------------------------------------------------------------------- */
1985 1 : poNode->SetStringField( "layerType", "thematic" );
1986 :
1987 1 : return( CE_None );
1988 : }
1989 :
1990 : /************************************************************************/
1991 : /* CreateOverview() */
1992 : /************************************************************************/
1993 :
1994 10 : int HFABand::CreateOverview( int nOverviewLevel, const char *pszResampling )
1995 :
1996 : {
1997 :
1998 10 : CPLString osLayerName;
1999 : int nOXSize, nOYSize;
2000 :
2001 10 : nOXSize = (psInfo->nXSize + nOverviewLevel - 1) / nOverviewLevel;
2002 10 : nOYSize = (psInfo->nYSize + nOverviewLevel - 1) / nOverviewLevel;
2003 :
2004 : /* -------------------------------------------------------------------- */
2005 : /* Do we want to use a dependent file (.rrd) for the overviews? */
2006 : /* Or just create them directly in this file? */
2007 : /* -------------------------------------------------------------------- */
2008 10 : HFAInfo_t *psRRDInfo = psInfo;
2009 10 : HFAEntry *poParent = poNode;
2010 :
2011 10 : if( CSLTestBoolean( CPLGetConfigOption( "HFA_USE_RRD", "NO" ) ) )
2012 : {
2013 2 : psRRDInfo = HFACreateDependent( psInfo );
2014 :
2015 2 : poParent = psRRDInfo->poRoot->GetNamedChild( GetBandName() );
2016 :
2017 : // Need to create layer object.
2018 2 : if( poParent == NULL )
2019 : {
2020 : poParent =
2021 : new HFAEntry( psRRDInfo, GetBandName(),
2022 1 : "Eimg_Layer", psRRDInfo->poRoot );
2023 : }
2024 : }
2025 :
2026 : /* -------------------------------------------------------------------- */
2027 : /* What pixel type should we use for the overview. Usually */
2028 : /* this is the same as the base layer, but when */
2029 : /* AVERAGE_BIT2GRAYSCALE is in effect we force it to u8 from u1. */
2030 : /* -------------------------------------------------------------------- */
2031 10 : int nOverviewDataType = nDataType;
2032 :
2033 10 : if( EQUALN(pszResampling,"AVERAGE_BIT2GR",14) )
2034 1 : nOverviewDataType = EPT_u8;
2035 :
2036 : /* -------------------------------------------------------------------- */
2037 : /* Eventually we need to decide on the whether to use the spill */
2038 : /* file, primarily on the basis of whether the new overview */
2039 : /* will drive our .img file size near 4GB. For now, just base */
2040 : /* it on the config options. */
2041 : /* -------------------------------------------------------------------- */
2042 : int bCreateLargeRaster = CSLTestBoolean(
2043 10 : CPLGetConfigOption("USE_SPILL","NO") );
2044 10 : GIntBig nValidFlagsOffset = 0, nDataOffset = 0;
2045 :
2046 10 : if( (psRRDInfo->nEndOfFile
2047 : + (nOXSize * (double) nOYSize)
2048 : * (HFAGetDataTypeBits(nOverviewDataType) / 8)) > 2000000000.0 )
2049 0 : bCreateLargeRaster = TRUE;
2050 :
2051 10 : if( bCreateLargeRaster )
2052 : {
2053 0 : if( !HFACreateSpillStack( psRRDInfo, nOXSize, nOYSize, 1,
2054 : 64, nOverviewDataType,
2055 : &nValidFlagsOffset, &nDataOffset ) )
2056 : {
2057 0 : return -1;
2058 : }
2059 : }
2060 :
2061 : /* -------------------------------------------------------------------- */
2062 : /* Create the layer. */
2063 : /* -------------------------------------------------------------------- */
2064 10 : osLayerName.Printf( "_ss_%d_", nOverviewLevel );
2065 :
2066 10 : if( !HFACreateLayer( psRRDInfo, poParent, osLayerName,
2067 : TRUE, 64, FALSE, bCreateLargeRaster, FALSE,
2068 : nOXSize, nOYSize, nOverviewDataType, NULL,
2069 : nValidFlagsOffset, nDataOffset, 1, 0 ) )
2070 0 : return -1;
2071 :
2072 10 : HFAEntry *poOverLayer = poParent->GetNamedChild( osLayerName );
2073 10 : if( poOverLayer == NULL )
2074 0 : return -1;
2075 :
2076 : /* -------------------------------------------------------------------- */
2077 : /* Create RRDNamesList list if it does not yet exist. */
2078 : /* -------------------------------------------------------------------- */
2079 10 : HFAEntry *poRRDNamesList = poNode->GetNamedChild("RRDNamesList");
2080 10 : if( poRRDNamesList == NULL )
2081 : {
2082 : poRRDNamesList = new HFAEntry( psInfo, "RRDNamesList",
2083 : "Eimg_RRDNamesList",
2084 4 : poNode );
2085 4 : poRRDNamesList->MakeData( 23+16+8+ 3000 /* hack for growth room*/ );
2086 :
2087 : /* we need to hardcode file offset into the data, so locate it now */
2088 4 : poRRDNamesList->SetPosition();
2089 :
2090 : poRRDNamesList->SetStringField( "algorithm.string",
2091 4 : "IMAGINE 2X2 Resampling" );
2092 : }
2093 :
2094 : /* -------------------------------------------------------------------- */
2095 : /* Add new overview layer to RRDNamesList. */
2096 : /* -------------------------------------------------------------------- */
2097 10 : int iNextName = poRRDNamesList->GetFieldCount( "nameList" );
2098 : char szName[50];
2099 10 : CPLString osNodeName;
2100 :
2101 10 : sprintf( szName, "nameList[%d].string", iNextName );
2102 :
2103 : osLayerName.Printf( "%s(:%s:_ss_%d_)",
2104 : psRRDInfo->pszFilename, GetBandName(),
2105 10 : nOverviewLevel );
2106 :
2107 : // TODO: Need to add to end of array (thats pretty hard).
2108 10 : if( poRRDNamesList->SetStringField( szName, osLayerName ) != CE_None )
2109 : {
2110 1 : poRRDNamesList->MakeData( poRRDNamesList->GetDataSize() + 3000 );
2111 1 : if( poRRDNamesList->SetStringField( szName, osLayerName ) != CE_None )
2112 0 : return -1;
2113 : }
2114 :
2115 : /* -------------------------------------------------------------------- */
2116 : /* Add to the list of overviews for this band. */
2117 : /* -------------------------------------------------------------------- */
2118 : papoOverviews = (HFABand **)
2119 10 : CPLRealloc(papoOverviews, sizeof(void*) * ++nOverviews );
2120 10 : papoOverviews[nOverviews-1] = new HFABand( psRRDInfo, poOverLayer );
2121 :
2122 : /* -------------------------------------------------------------------- */
2123 : /* If there is a nodata value, copy it to the overview band. */
2124 : /* -------------------------------------------------------------------- */
2125 10 : if( bNoDataSet )
2126 1 : papoOverviews[nOverviews-1]->SetNoDataValue( dfNoData );
2127 :
2128 10 : return nOverviews-1;
2129 : }
|
