Line data    Source code

       1             : /******************************************************************************
       2             :  *
       3             :  * Project:  GRC/GRD Reader
       4             :  * Purpose:  Northwood Format basic implementation
       5             :  * Author:   Perry Casson
       6             :  *
       7             :  ******************************************************************************
       8             :  * Copyright (c) 2007, Waypoint Information Technology
       9             :  * Copyright (c) 2009-2011, Even Rouault <even dot rouault at spatialys.com>
      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             : 
      32             : #include "northwood.h"
      33             : 
      34             : #include <algorithm>
      35             : #include <limits>
      36             : #include <string>
      37             : 
      38          14 : int nwt_ParseHeader(NWT_GRID *pGrd, const unsigned char *nwtHeader)
      39             : {
      40             :     /* double dfTmp; */
      41             : 
      42          14 :     if (nwtHeader[4] == '1')
      43          12 :         pGrd->cFormat = 0x00;  // grd - surface type
      44           2 :     else if (nwtHeader[4] == '8')
      45           2 :         pGrd->cFormat = 0x80;  //  grc classified type
      46             : 
      47          14 :     pGrd->stClassDict = nullptr;
      48             : 
      49          14 :     memcpy(&pGrd->fVersion, &nwtHeader[5], sizeof(pGrd->fVersion));
      50          14 :     CPL_LSBPTR32(&pGrd->fVersion);
      51             : 
      52             :     unsigned short usTmp;
      53          14 :     memcpy(&usTmp, &nwtHeader[9], 2);
      54          14 :     CPL_LSBPTR16(&usTmp);
      55          14 :     pGrd->nXSide = static_cast<unsigned int>(usTmp);
      56          14 :     if (pGrd->nXSide == 0)
      57             :     {
      58           0 :         memcpy(&pGrd->nXSide, &nwtHeader[128], sizeof(pGrd->nXSide));
      59           0 :         CPL_LSBPTR32(&pGrd->nXSide);
      60             :     }
      61          14 :     if (pGrd->nXSide <= 1)
      62           0 :         return FALSE;
      63             : 
      64          14 :     memcpy(&usTmp, &nwtHeader[11], 2);
      65          14 :     CPL_LSBPTR16(&usTmp);
      66          14 :     pGrd->nYSide = static_cast<unsigned int>(usTmp);
      67          14 :     if (pGrd->nYSide == 0)
      68             :     {
      69           0 :         memcpy(&pGrd->nYSide, &nwtHeader[132], sizeof(pGrd->nYSide));
      70           0 :         CPL_LSBPTR32(&pGrd->nYSide);
      71             :     }
      72             : 
      73          14 :     memcpy(&pGrd->dfMinX, &nwtHeader[13], sizeof(pGrd->dfMinX));
      74          14 :     CPL_LSBPTR64(&pGrd->dfMinX);
      75          14 :     memcpy(&pGrd->dfMaxX, &nwtHeader[21], sizeof(pGrd->dfMaxX));
      76          14 :     CPL_LSBPTR64(&pGrd->dfMaxX);
      77          14 :     memcpy(&pGrd->dfMinY, &nwtHeader[29], sizeof(pGrd->dfMinY));
      78          14 :     CPL_LSBPTR64(&pGrd->dfMinY);
      79          14 :     memcpy(&pGrd->dfMaxY, &nwtHeader[37], sizeof(pGrd->dfMaxY));
      80          14 :     CPL_LSBPTR64(&pGrd->dfMaxY);
      81             : 
      82          14 :     pGrd->dfStepSize = (pGrd->dfMaxX - pGrd->dfMinX) / (pGrd->nXSide - 1);
      83             :     /* dfTmp = (pGrd->dfMaxY - pGrd->dfMinY) / (pGrd->nYSide - 1); */
      84             : 
      85          14 :     memcpy(&pGrd->fZMin, &nwtHeader[45], sizeof(pGrd->fZMin));
      86          14 :     CPL_LSBPTR32(&pGrd->fZMin);
      87          14 :     memcpy(&pGrd->fZMax, &nwtHeader[49], sizeof(pGrd->fZMax));
      88          14 :     CPL_LSBPTR32(&pGrd->fZMax);
      89          14 :     memcpy(&pGrd->fZMinScale, &nwtHeader[53], sizeof(pGrd->fZMinScale));
      90          14 :     CPL_LSBPTR32(&pGrd->fZMinScale);
      91          14 :     memcpy(&pGrd->fZMaxScale, &nwtHeader[57], sizeof(pGrd->fZMaxScale));
      92          14 :     CPL_LSBPTR32(&pGrd->fZMaxScale);
      93             : 
      94          14 :     memcpy(&pGrd->cDescription, &nwtHeader[61], sizeof(pGrd->cDescription));
      95          14 :     memcpy(&pGrd->cZUnits, &nwtHeader[93], sizeof(pGrd->cZUnits));
      96             : 
      97             :     int i;
      98          14 :     memcpy(&i, &nwtHeader[136], 4);
      99          14 :     CPL_LSBPTR32(&i);
     100             : 
     101          14 :     if (i == 1129336130)
     102             :     {  // BMPC
     103           0 :         if (nwtHeader[140] & 0x01)
     104             :         {
     105           0 :             pGrd->cHillShadeBrightness = nwtHeader[144];
     106           0 :             pGrd->cHillShadeContrast = nwtHeader[145];
     107             :         }
     108             :     }
     109             : 
     110          14 :     memcpy(&pGrd->cMICoordSys, &nwtHeader[256], sizeof(pGrd->cMICoordSys));
     111          14 :     pGrd->cMICoordSys[sizeof(pGrd->cMICoordSys) - 1] = '\0';
     112             : 
     113          14 :     pGrd->iZUnits = nwtHeader[512];
     114             : 
     115          14 :     if (nwtHeader[513] & 0x80)
     116           0 :         pGrd->bShowGradient = true;
     117             : 
     118          14 :     if (nwtHeader[513] & 0x40)
     119           0 :         pGrd->bShowHillShade = true;
     120             : 
     121          14 :     if (nwtHeader[513] & 0x20)
     122           0 :         pGrd->bHillShadeExists = true;
     123             : 
     124          14 :     memcpy(&pGrd->iNumColorInflections, &nwtHeader[516], 2);
     125          14 :     CPL_LSBPTR16(&pGrd->iNumColorInflections);
     126             : 
     127          14 :     if (pGrd->iNumColorInflections > 32)
     128             :     {
     129           0 :         CPLError(CE_Failure, CPLE_AppDefined, "Corrupt header");
     130           0 :         pGrd->iNumColorInflections = 0;
     131           0 :         return FALSE;
     132             :     }
     133             : 
     134          68 :     for (i = 0; i < pGrd->iNumColorInflections; i++)
     135             :     {
     136          54 :         memcpy(&pGrd->stInflection[i].zVal, &nwtHeader[518 + (7 * i)], 4);
     137          54 :         CPL_LSBPTR32(&pGrd->stInflection[i].zVal);
     138          54 :         memcpy(&pGrd->stInflection[i].r, &nwtHeader[522 + (7 * i)], 1);
     139          54 :         memcpy(&pGrd->stInflection[i].g, &nwtHeader[523 + (7 * i)], 1);
     140          54 :         memcpy(&pGrd->stInflection[i].b, &nwtHeader[524 + (7 * i)], 1);
     141             :     }
     142             : 
     143          14 :     memcpy(&pGrd->fHillShadeAzimuth, &nwtHeader[966],
     144             :            sizeof(pGrd->fHillShadeAzimuth));
     145          14 :     CPL_LSBPTR32(&pGrd->fHillShadeAzimuth);
     146          14 :     memcpy(&pGrd->fHillShadeAngle, &nwtHeader[970],
     147             :            sizeof(pGrd->fHillShadeAngle));
     148          14 :     CPL_LSBPTR32(&pGrd->fHillShadeAngle);
     149             : 
     150          14 :     pGrd->cFormat += nwtHeader[1023];  // the msb for grd/grc was already set
     151             : 
     152             :     // there are more types than this - need to build other types for testing
     153          14 :     if (pGrd->cFormat & 0x80)
     154             :     {
     155           2 :         if (nwtHeader[1023] == 0)
     156           0 :             pGrd->nBitsPerPixel = 16;
     157             :         else
     158           2 :             pGrd->nBitsPerPixel = nwtHeader[1023] * 4;
     159             :     }
     160             :     else
     161          12 :         pGrd->nBitsPerPixel = nwtHeader[1023] * 8;
     162             : 
     163          14 :     if (pGrd->cFormat & 0x80)  // if is GRC load the Dictionary
     164             :     {
     165           2 :         vsi_l_offset nPixels =
     166           2 :             static_cast<vsi_l_offset>(pGrd->nXSide) * pGrd->nYSide;
     167           2 :         unsigned int nBytesPerPixel = pGrd->nBitsPerPixel / 8;
     168           4 :         if (nPixels > 0 &&
     169           2 :             (nBytesPerPixel >
     170           2 :                  std::numeric_limits<vsi_l_offset>::max() / nPixels ||
     171           2 :              nPixels * nBytesPerPixel >
     172           2 :                  std::numeric_limits<vsi_l_offset>::max() - 1024))
     173             :         {
     174           0 :             CPLError(CE_Failure, CPLE_FileIO,
     175             :                      "Invalid file dimension / bits per pixel");
     176           0 :             return FALSE;
     177             :         }
     178           2 :         VSIFSeekL(pGrd->fp, 1024 + nPixels * nBytesPerPixel, SEEK_SET);
     179             : 
     180           2 :         if (!VSIFReadL(&usTmp, 2, 1, pGrd->fp))
     181             :         {
     182           0 :             CPLError(CE_Failure, CPLE_FileIO, "Read failure, file short?");
     183           0 :             return FALSE;
     184             :         }
     185           2 :         CPL_LSBPTR16(&usTmp);
     186           2 :         pGrd->stClassDict = reinterpret_cast<NWT_CLASSIFIED_DICT *>(
     187           2 :             calloc(sizeof(NWT_CLASSIFIED_DICT), 1));
     188             : 
     189           2 :         pGrd->stClassDict->nNumClassifiedItems = usTmp;
     190             : 
     191           2 :         pGrd->stClassDict->stClassifiedItem =
     192             :             reinterpret_cast<NWT_CLASSIFIED_ITEM **>(
     193           2 :                 calloc(sizeof(NWT_CLASSIFIED_ITEM *),
     194           2 :                        pGrd->stClassDict->nNumClassifiedItems + 1));
     195             : 
     196             :         // load the dictionary
     197           8 :         for (unsigned int iItem = 0;
     198           8 :              iItem < pGrd->stClassDict->nNumClassifiedItems; iItem++)
     199             :         {
     200             :             NWT_CLASSIFIED_ITEM *psItem =
     201           6 :                 pGrd->stClassDict->stClassifiedItem[iItem] =
     202             :                     reinterpret_cast<NWT_CLASSIFIED_ITEM *>(
     203           6 :                         calloc(sizeof(NWT_CLASSIFIED_ITEM), 1));
     204             : 
     205             :             unsigned char cTmp[256];
     206           6 :             if (!VSIFReadL(&cTmp, 9, 1, pGrd->fp))
     207             :             {
     208           0 :                 CPLError(CE_Failure, CPLE_FileIO, "Read failure, file short?");
     209           0 :                 return FALSE;
     210             :             }
     211           6 :             memcpy(&psItem->usPixVal, &cTmp[0], 2);
     212           6 :             CPL_LSBPTR16(&psItem->usPixVal);
     213           6 :             memcpy(&psItem->res1, &cTmp[2], 1);
     214           6 :             memcpy(&psItem->r, &cTmp[3], 1);
     215           6 :             memcpy(&psItem->g, &cTmp[4], 1);
     216           6 :             memcpy(&psItem->b, &cTmp[5], 1);
     217           6 :             memcpy(&psItem->res2, &cTmp[6], 1);
     218           6 :             memcpy(&psItem->usLen, &cTmp[7], 2);
     219           6 :             CPL_LSBPTR16(&psItem->usLen);
     220             : 
     221           6 :             if (psItem->usLen > sizeof(psItem->szClassName) - 1)
     222             :             {
     223           0 :                 CPLError(CE_Failure, CPLE_AppDefined,
     224             :                          "Unexpected long class name, %d characters long - "
     225             :                          "unable to read file.",
     226           0 :                          psItem->usLen);
     227           0 :                 return FALSE;
     228             :             }
     229             : 
     230             :             // 0-len class names are possible
     231           6 :             psItem->szClassName[0] = '\0';
     232          12 :             if (psItem->usLen > 0 &&
     233           6 :                 !VSIFReadL(&psItem->szClassName, psItem->usLen, 1, pGrd->fp))
     234           0 :                 return FALSE;
     235             :         }
     236             :     }
     237             : 
     238          14 :     return TRUE;
     239             : }
     240             : 
     241             : // Create a color gradient ranging from ZMin to Zmax using the color
     242             : // inflections defined in grid
     243          12 : int nwt_LoadColors(NWT_RGB *pMap, int mapSize, NWT_GRID *pGrd)
     244             : {
     245             :     int i;
     246             :     NWT_RGB sColor;
     247          12 :     int nWarkerMark = 0;
     248             : 
     249          12 :     createIP(0, 255, 255, 255, pMap, &nWarkerMark);
     250          12 :     if (pGrd->iNumColorInflections == 0)
     251           0 :         return 0;
     252             : 
     253             :     // If Zmin is less than the 1st inflection use the 1st inflections color to
     254             :     // the start of the ramp
     255          12 :     if (pGrd->fZMin <= pGrd->stInflection[0].zVal)
     256             :     {
     257          12 :         createIP(1, pGrd->stInflection[0].r, pGrd->stInflection[0].g,
     258          12 :                  pGrd->stInflection[0].b, pMap, &nWarkerMark);
     259             :     }
     260             :     // find what inflections zmin is between
     261          12 :     for (i = 1; i < pGrd->iNumColorInflections; i++)
     262             :     {
     263          12 :         if (pGrd->fZMin < pGrd->stInflection[i].zVal)
     264             :         {
     265             :             // then we must be between i and i-1
     266          12 :             linearColor(&sColor, &pGrd->stInflection[i - 1],
     267             :                         &pGrd->stInflection[i], pGrd->fZMin);
     268          12 :             createIP(1, sColor.r, sColor.g, sColor.b, pMap, &nWarkerMark);
     269          12 :             break;
     270             :         }
     271             :     }
     272             :     // the interesting case of zmin beig higher than the max inflection value
     273          12 :     if (i >= pGrd->iNumColorInflections)
     274             :     {
     275           0 :         createIP(1, pGrd->stInflection[pGrd->iNumColorInflections - 1].r,
     276           0 :                  pGrd->stInflection[pGrd->iNumColorInflections - 1].g,
     277           0 :                  pGrd->stInflection[pGrd->iNumColorInflections - 1].b, pMap,
     278             :                  &nWarkerMark);
     279           0 :         createIP(mapSize - 1,
     280           0 :                  pGrd->stInflection[pGrd->iNumColorInflections - 1].r,
     281           0 :                  pGrd->stInflection[pGrd->iNumColorInflections - 1].g,
     282           0 :                  pGrd->stInflection[pGrd->iNumColorInflections - 1].b, pMap,
     283             :                  &nWarkerMark);
     284             :     }
     285             :     else
     286             :     {
     287          12 :         int index = 0;
     288          54 :         for (; i < pGrd->iNumColorInflections; i++)
     289             :         {
     290          42 :             if (pGrd->fZMax < pGrd->stInflection[i].zVal)
     291             :             {
     292             :                 // then we must be between i and i-1
     293           0 :                 linearColor(&sColor, &pGrd->stInflection[i - 1],
     294             :                             &pGrd->stInflection[i], pGrd->fZMax);
     295           0 :                 index = mapSize - 1;
     296           0 :                 createIP(index, sColor.r, sColor.g, sColor.b, pMap,
     297             :                          &nWarkerMark);
     298           0 :                 break;
     299             :             }
     300             :             // save the inflections between zmin and zmax
     301          42 :             index =
     302          42 :                 static_cast<int>(((pGrd->stInflection[i].zVal - pGrd->fZMin) /
     303          42 :                                   (pGrd->fZMax - pGrd->fZMin)) *
     304             :                                  mapSize);
     305             : 
     306          42 :             if (index >= mapSize)
     307          12 :                 index = mapSize - 1;
     308          42 :             createIP(index, pGrd->stInflection[i].r, pGrd->stInflection[i].g,
     309          42 :                      pGrd->stInflection[i].b, pMap, &nWarkerMark);
     310             :         }
     311          12 :         if (index < mapSize - 1)
     312           0 :             createIP(mapSize - 1,
     313           0 :                      pGrd->stInflection[pGrd->iNumColorInflections - 1].r,
     314           0 :                      pGrd->stInflection[pGrd->iNumColorInflections - 1].g,
     315           0 :                      pGrd->stInflection[pGrd->iNumColorInflections - 1].b, pMap,
     316             :                      &nWarkerMark);
     317             :     }
     318          12 :     return 0;
     319             : }
     320             : 
     321             : // solve for a color between pIPLow and pIPHigh
     322          12 : void linearColor(NWT_RGB *pRGB, NWT_INFLECTION *pIPLow, NWT_INFLECTION *pIPHigh,
     323             :                  float fMid)
     324             : {
     325          12 :     if (fMid < pIPLow->zVal)
     326             :     {
     327           0 :         pRGB->r = pIPLow->r;
     328           0 :         pRGB->g = pIPLow->g;
     329           0 :         pRGB->b = pIPLow->b;
     330             :     }
     331          12 :     else if (fMid > pIPHigh->zVal)
     332             :     {
     333           0 :         pRGB->r = pIPHigh->r;
     334           0 :         pRGB->g = pIPHigh->g;
     335           0 :         pRGB->b = pIPHigh->b;
     336             :     }
     337             :     else
     338             :     {
     339          12 :         float scale = (fMid - pIPLow->zVal) / (pIPHigh->zVal - pIPLow->zVal);
     340          12 :         pRGB->r = static_cast<unsigned char>(scale * (pIPHigh->r - pIPLow->r) +
     341          12 :                                              pIPLow->r + 0.5);
     342          12 :         pRGB->g = static_cast<unsigned char>(scale * (pIPHigh->g - pIPLow->g) +
     343          12 :                                              pIPLow->g + 0.5);
     344          12 :         pRGB->b = static_cast<unsigned char>(scale * (pIPHigh->b - pIPLow->b) +
     345          12 :                                              pIPLow->b + 0.5);
     346             :     }
     347          12 : }
     348             : 
     349             : // insert IP's into the map filling as we go
     350          78 : void createIP(int index, unsigned char r, unsigned char g, unsigned char b,
     351             :               NWT_RGB *map, int *pnWarkerMark)
     352             : {
     353             :     int i;
     354             : 
     355          78 :     if (index == 0)
     356             :     {
     357          12 :         map[0].r = r;
     358          12 :         map[0].g = g;
     359          12 :         map[0].b = b;
     360          12 :         *pnWarkerMark = 0;
     361          12 :         return;
     362             :     }
     363             : 
     364          66 :     if (index <= *pnWarkerMark)
     365          12 :         return;
     366             : 
     367          54 :     int wm = *pnWarkerMark;
     368             : 
     369          54 :     float rslope =
     370          54 :         static_cast<float>(r - map[wm].r) / static_cast<float>(index - wm);
     371          54 :     float gslope =
     372          54 :         static_cast<float>(g - map[wm].g) / static_cast<float>(index - wm);
     373          54 :     float bslope =
     374          54 :         static_cast<float>(b - map[wm].b) / static_cast<float>(index - wm);
     375       49140 :     for (i = wm + 1; i < index; i++)
     376             :     {
     377       49086 :         map[i].r =
     378       49086 :             static_cast<unsigned char>(map[wm].r + ((i - wm) * rslope) + 0.5);
     379       49086 :         map[i].g =
     380       49086 :             static_cast<unsigned char>(map[wm].g + ((i - wm) * gslope) + 0.5);
     381       49086 :         map[i].b =
     382       49086 :             static_cast<unsigned char>(map[wm].b + ((i - wm) * bslope) + 0.5);
     383             :     }
     384          54 :     map[index].r = r;
     385          54 :     map[index].g = g;
     386          54 :     map[index].b = b;
     387          54 :     *pnWarkerMark = index;
     388          54 :     return;
     389             : }
     390             : 
     391           0 : void nwt_HillShade(unsigned char *r, unsigned char *g, unsigned char *b,
     392             :                    unsigned char *h)
     393             : {
     394             :     HLS hls;
     395             :     NWT_RGB rgb;
     396           0 :     rgb.r = *r;
     397           0 :     rgb.g = *g;
     398           0 :     rgb.b = *b;
     399           0 :     hls = RGBtoHLS(rgb);
     400           0 :     hls.l = static_cast<short>(hls.l + (*h) * HLSMAX / 256);
     401           0 :     rgb = HLStoRGB(hls);
     402             : 
     403           0 :     *r = rgb.r;
     404           0 :     *g = rgb.g;
     405           0 :     *b = rgb.b;
     406           0 :     return;
     407             : }
     408             : 
     409           0 : NWT_GRID *nwtOpenGrid(char *filename)
     410             : {
     411             :     unsigned char nwtHeader[1024];
     412           0 :     VSILFILE *fp = VSIFOpenL(filename, "rb");
     413             : 
     414           0 :     if (fp == nullptr)
     415             :     {
     416           0 :         CPLError(CE_Failure, CPLE_OpenFailed, "Can't open %s", filename);
     417           0 :         return nullptr;
     418             :     }
     419             : 
     420           0 :     if (!VSIFReadL(nwtHeader, 1024, 1, fp))
     421           0 :         return nullptr;
     422             : 
     423           0 :     if (nwtHeader[0] != 'H' || nwtHeader[1] != 'G' || nwtHeader[2] != 'P' ||
     424           0 :         nwtHeader[3] != 'C')
     425           0 :         return nullptr;
     426             : 
     427           0 :     NWT_GRID *pGrd = reinterpret_cast<NWT_GRID *>(calloc(sizeof(NWT_GRID), 1));
     428             : 
     429           0 :     if (nwtHeader[4] == '1')
     430           0 :         pGrd->cFormat = 0x00;  // grd - surface type
     431           0 :     else if (nwtHeader[4] == '8')
     432           0 :         pGrd->cFormat = 0x80;  //  grc classified type
     433             :     else
     434             :     {
     435           0 :         CPLError(CE_Failure, CPLE_NotSupported,
     436           0 :                  "Unhandled Northwood format type = %0xd", nwtHeader[4]);
     437           0 :         if (pGrd)
     438           0 :             free(pGrd);
     439           0 :         return nullptr;
     440             :     }
     441             : 
     442           0 :     strncpy(pGrd->szFileName, filename, sizeof(pGrd->szFileName));
     443           0 :     pGrd->szFileName[sizeof(pGrd->szFileName) - 1] = '\0';
     444           0 :     pGrd->fp = fp;
     445           0 :     nwt_ParseHeader(pGrd, nwtHeader);
     446             : 
     447           0 :     return pGrd;
     448             : }
     449             : 
     450             : // close the file and free the mem
     451          28 : void nwtCloseGrid(NWT_GRID *pGrd)
     452             : {
     453          30 :     if ((pGrd->cFormat & 0x80) &&
     454          28 :         pGrd->stClassDict)  // if is GRC - free the Dictionary
     455             :     {
     456           8 :         for (unsigned int i = 0; i < pGrd->stClassDict->nNumClassifiedItems;
     457             :              i++)
     458             :         {
     459           6 :             free(pGrd->stClassDict->stClassifiedItem[i]);
     460             :         }
     461           2 :         free(pGrd->stClassDict->stClassifiedItem);
     462           2 :         free(pGrd->stClassDict);
     463             :     }
     464          28 :     if (pGrd->fp)
     465           0 :         VSIFCloseL(pGrd->fp);
     466          28 :     free(pGrd);
     467          28 :     return;
     468             : }
     469             : 
     470           0 : void nwtPrintGridHeader(NWT_GRID *pGrd)
     471             : {
     472           0 :     if (pGrd->cFormat & 0x80)
     473             :     {
     474           0 :         printf("\n%s\n\nGrid type is Classified ", pGrd->szFileName); /*ok*/
     475           0 :         if (pGrd->cFormat == 0x81)
     476           0 :             printf("4 bit (Less than 16 Classes)"); /*ok*/
     477           0 :         else if (pGrd->cFormat == 0x82)
     478           0 :             printf("8 bit (Less than 256 Classes)"); /*ok*/
     479           0 :         else if (pGrd->cFormat == 0x84)
     480           0 :             printf("16 bit (Less than 65536 Classes)"); /*ok*/
     481             :         else
     482             :         {
     483           0 :             printf("GRC - Unhandled Format or Type %d", pGrd->cFormat); /*ok*/
     484           0 :             return;
     485             :         }
     486             :     }
     487             :     else
     488             :     {
     489           0 :         printf("\n%s\n\nGrid type is Numeric ", pGrd->szFileName); /*ok*/
     490           0 :         if (pGrd->cFormat == 0x00)
     491           0 :             printf("16 bit (Standard Precision)"); /*ok*/
     492           0 :         else if (pGrd->cFormat == 0x01)
     493           0 :             printf("32 bit (High Precision)"); /*ok*/
     494             :         else
     495             :         {
     496           0 :             printf("GRD - Unhandled Format or Type %d", pGrd->cFormat); /*ok*/
     497           0 :             return;
     498             :         }
     499             :     }
     500           0 :     printf("\nDim (x,y) = (%u,%u)", pGrd->nXSide, pGrd->nYSide);     /*ok*/
     501           0 :     printf("\nStep Size = %f", pGrd->dfStepSize);                    /*ok*/
     502           0 :     printf("\nBounds = (%f,%f) (%f,%f)", pGrd->dfMinX, pGrd->dfMinY, /*ok*/
     503             :            pGrd->dfMaxX, pGrd->dfMaxY);
     504           0 :     printf("\nCoordinate System = %s", pGrd->cMICoordSys); /*ok*/
     505             : 
     506           0 :     if (!(pGrd->cFormat & 0x80))  // print the numeric specific stuff
     507             :     {
     508           0 :         printf("\nMin Z = %f Max Z = %f Z Units = %d \"%s\"", /*ok*/
     509           0 :                pGrd->fZMin, pGrd->fZMax, pGrd->iZUnits, pGrd->cZUnits);
     510             : 
     511           0 :         printf("\n\nDisplay Mode ="); /*ok*/
     512           0 :         if (pGrd->bShowGradient)
     513           0 :             printf(" Color Gradient"); /*ok*/
     514             : 
     515           0 :         if (pGrd->bShowGradient && pGrd->bShowHillShade)
     516           0 :             printf(" and"); /*ok*/
     517             : 
     518           0 :         if (pGrd->bShowHillShade)
     519           0 :             printf(" Hill Shading"); /*ok*/
     520             : 
     521           0 :         for (int i = 0; i < pGrd->iNumColorInflections; i++)
     522             :         {
     523           0 :             printf("\nColor Inflection %d - %f (%d,%d,%d)", i + 1, /*ok*/
     524           0 :                    pGrd->stInflection[i].zVal, pGrd->stInflection[i].r,
     525           0 :                    pGrd->stInflection[i].g, pGrd->stInflection[i].b);
     526             :         }
     527             : 
     528           0 :         if (pGrd->bHillShadeExists)
     529             :         {
     530           0 :             printf("\n\nHill Shade Azumith = %.1f Inclination = %.1f " /*ok*/
     531             :                    "Brightness = %d Contrast = %d",
     532           0 :                    pGrd->fHillShadeAzimuth, pGrd->fHillShadeAngle,
     533           0 :                    pGrd->cHillShadeBrightness, pGrd->cHillShadeContrast);
     534             :         }
     535             :         else
     536           0 :             printf("\n\nNo Hill Shade Data"); /*ok*/
     537             :     }
     538             :     else  // print the classified specific stuff
     539             :     {
     540           0 :         printf("\nNumber of Classes defined = %u", /*ok*/
     541           0 :                pGrd->stClassDict->nNumClassifiedItems);
     542           0 :         for (int i = 0;
     543           0 :              i < static_cast<int>(pGrd->stClassDict->nNumClassifiedItems); i++)
     544             :         {
     545           0 :             printf("\n%s - (%d,%d,%d)  Raw = %d  %d %d", /*ok*/
     546           0 :                    pGrd->stClassDict->stClassifiedItem[i]->szClassName,
     547           0 :                    pGrd->stClassDict->stClassifiedItem[i]->r,
     548           0 :                    pGrd->stClassDict->stClassifiedItem[i]->g,
     549           0 :                    pGrd->stClassDict->stClassifiedItem[i]->b,
     550           0 :                    pGrd->stClassDict->stClassifiedItem[i]->usPixVal,
     551           0 :                    pGrd->stClassDict->stClassifiedItem[i]->res1,
     552           0 :                    pGrd->stClassDict->stClassifiedItem[i]->res2);
     553             :         }
     554             :     }
     555             : }
     556             : 
     557           0 : HLS RGBtoHLS(NWT_RGB rgb)
     558             : {
     559             :     /* get R, G, and B out of DWORD */
     560           0 :     short R = rgb.r;
     561           0 :     short G = rgb.g;
     562           0 :     short B = rgb.b;
     563             : 
     564             :     /* calculate lightness */
     565             :     unsigned char cMax =
     566           0 :         static_cast<unsigned char>(std::max(std::max(R, G), B));
     567             :     unsigned char cMin =
     568           0 :         static_cast<unsigned char>(std::min(std::min(R, G), B));
     569             :     HLS hls;
     570           0 :     hls.l = (((cMax + cMin) * HLSMAX) + RGBMAX) / (2 * RGBMAX);
     571             : 
     572             :     short Rdelta, Gdelta, Bdelta; /* intermediate value: % of spread from max */
     573           0 :     if (cMax == cMin)
     574             :     {                      /* r=g=b --> achromatic case */
     575           0 :         hls.s = 0;         /* saturation */
     576           0 :         hls.h = UNDEFINED; /* hue */
     577             :     }
     578             :     else
     579             :     { /* chromatic case */
     580             :         /* saturation */
     581           0 :         if (hls.l <= (HLSMAX / 2))
     582           0 :             hls.s = (((cMax - cMin) * HLSMAX) + ((cMax + cMin) / 2)) /
     583           0 :                     (cMax + cMin);
     584             :         else
     585           0 :             hls.s =
     586           0 :                 (((cMax - cMin) * HLSMAX) + ((2 * RGBMAX - cMax - cMin) / 2)) /
     587           0 :                 (2 * RGBMAX - cMax - cMin);
     588             : 
     589             :         /* hue */
     590           0 :         Rdelta =
     591           0 :             (((cMax - R) * (HLSMAX / 6)) + ((cMax - cMin) / 2)) / (cMax - cMin);
     592           0 :         Gdelta =
     593           0 :             (((cMax - G) * (HLSMAX / 6)) + ((cMax - cMin) / 2)) / (cMax - cMin);
     594           0 :         Bdelta =
     595           0 :             (((cMax - B) * (HLSMAX / 6)) + ((cMax - cMin) / 2)) / (cMax - cMin);
     596             : 
     597           0 :         if (R == cMax)
     598           0 :             hls.h = Bdelta - Gdelta;
     599           0 :         else if (G == cMax)
     600           0 :             hls.h = (HLSMAX / 3) + Rdelta - Bdelta;
     601             :         else /* B == cMax */
     602           0 :             hls.h = ((2 * HLSMAX) / 3) + Gdelta - Rdelta;
     603             : 
     604           0 :         if (hls.h < 0)
     605           0 :             hls.h += HLSMAX;
     606           0 :         if (hls.h > HLSMAX)
     607           0 :             hls.h -= HLSMAX;
     608             :     }
     609           0 :     return hls;
     610             : }
     611             : 
     612             : /* utility routine for HLStoRGB */
     613           0 : static short HueToRGB(short n1, short n2, short hue)
     614             : {
     615             :     /* range check: note values passed add/subtract thirds of range */
     616           0 :     if (hue < 0)
     617           0 :         hue += HLSMAX;
     618             : 
     619           0 :     if (hue > HLSMAX)
     620           0 :         hue -= HLSMAX;
     621             : 
     622             :     /* return r,g, or b value from this tridrant */
     623           0 :     if (hue < (HLSMAX / 6))
     624           0 :         return n1 + (((n2 - n1) * hue + (HLSMAX / 12)) / (HLSMAX / 6));
     625           0 :     if (hue < (HLSMAX / 2))
     626           0 :         return n2;
     627           0 :     if (hue < ((HLSMAX * 2) / 3))
     628           0 :         return n1 + (((n2 - n1) * (((HLSMAX * 2) / 3) - hue) + (HLSMAX / 12)) /
     629           0 :                      (HLSMAX / 6));
     630             :     else
     631           0 :         return n1;
     632             : }
     633             : 
     634           0 : NWT_RGB HLStoRGB(HLS hls)
     635             : {
     636             :     NWT_RGB rgb;
     637             : 
     638           0 :     if (hls.s == 0)
     639             :     { /* achromatic case */
     640           0 :         rgb.r = static_cast<unsigned char>((hls.l * RGBMAX) / HLSMAX);
     641           0 :         rgb.g = rgb.r;
     642           0 :         rgb.b = rgb.r;
     643           0 :         if (hls.h != UNDEFINED)
     644             :         {
     645             :             /* ERROR */
     646             :         }
     647             :     }
     648             :     else
     649             :     { /* chromatic case */
     650             :         /* set up magic numbers */
     651             :         short Magic1, Magic2; /* calculated magic numbers (really!) */
     652           0 :         if (hls.l <= (HLSMAX / 2))
     653           0 :             Magic2 = (hls.l * (HLSMAX + hls.s) + (HLSMAX / 2)) / HLSMAX;
     654             :         else
     655           0 :             Magic2 = hls.l + hls.s - ((hls.l * hls.s) + (HLSMAX / 2)) / HLSMAX;
     656           0 :         Magic1 = 2 * hls.l - Magic2;
     657             : 
     658             :         /* get RGB, change units from HLSMAX to RGBMAX */
     659           0 :         rgb.r = static_cast<unsigned char>(
     660           0 :             (HueToRGB(Magic1, Magic2, hls.h + (HLSMAX / 3)) * RGBMAX +
     661           0 :              (HLSMAX / 2)) /
     662             :             HLSMAX);
     663           0 :         rgb.g = static_cast<unsigned char>(
     664           0 :             (HueToRGB(Magic1, Magic2, hls.h) * RGBMAX + (HLSMAX / 2)) / HLSMAX);
     665           0 :         rgb.b = static_cast<unsigned char>(
     666           0 :             (HueToRGB(Magic1, Magic2, hls.h - (HLSMAX / 3)) * RGBMAX +
     667           0 :              (HLSMAX / 2)) /
     668             :             HLSMAX);
     669             :     }
     670             : 
     671           0 :     return rgb;
     672             : }