LCOV - code coverage report
Current view: directory - frmts/pds - isis2dataset.cpp (source / functions) Found Hit Coverage
Test: gdal_filtered.info Lines: 480 382 79.6 %
Date: 2011-12-18 Functions: 21 18 85.7 %

       1                 : /******************************************************************************
       2                 :  * $Id: isis2dataset.cpp 23371 2011-11-13 14:24:00Z rouault $
       3                 :  *
       4                 :  * Project:  ISIS Version 2 Driver
       5                 :  * Purpose:  Implementation of ISIS2Dataset
       6                 :  * Author:   Trent Hare (thare@usgs.gov),
       7                 :  *           Robert Soricone (rsoricone@usgs.gov)
       8                 :  *           Ludovic Mercier (ludovic.mercier@gmail.com) 
       9                 :  *           Frank Warmerdam (warmerdam@pobox.com)
      10                 :  *
      11                 :  * NOTE: Original code authored by Trent and Robert and placed in the public 
      12                 :  * domain as per US government policy.  I have (within my rights) appropriated 
      13                 :  * it and placed it under the following license.  This is not intended to 
      14                 :  * diminish Trent and Roberts contribution. 
      15                 :  ******************************************************************************
      16                 :  * Copyright (c) 2006, Frank Warmerdam <warmerdam@pobox.com>
      17                 :  * 
      18                 :  * Permission is hereby granted, free of charge, to any person obtaining a
      19                 :  * copy of this software and associated documentation files (the "Software"),
      20                 :  * to deal in the Software without restriction, including without limitation
      21                 :  * the rights to use, copy, modify, merge, publish, distribute, sublicense,
      22                 :  * and/or sell copies of the Software, and to permit persons to whom the
      23                 :  * Software is furnished to do so, subject to the following conditions:
      24                 :  *
      25                 :  * The above copyright notice and this permission notice shall be included
      26                 :  * in all copies or substantial portions of the Software.
      27                 :  *
      28                 :  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
      29                 :  * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
      30                 :  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
      31                 :  * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
      32                 :  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
      33                 :  * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
      34                 :  * DEALINGS IN THE SOFTWARE.
      35                 :  ****************************************************************************/
      36                 : 
      37                 : #define NULL1 0
      38                 : #define NULL2 -32768
      39                 : #define NULL3 -3.4028226550889044521e+38
      40                 : 
      41                 : #ifndef PI
      42                 : #  define PI 3.1415926535897932384626433832795
      43                 : #endif
      44                 : 
      45                 : #define RECORD_SIZE 512
      46                 : 
      47                 : #include "rawdataset.h"
      48                 : #include "ogr_spatialref.h"
      49                 : #include "cpl_string.h" 
      50                 : #include "nasakeywordhandler.h"
      51                 : 
      52                 : CPL_CVSID("$Id: isis2dataset.cpp 23371 2011-11-13 14:24:00Z rouault $");
      53                 : 
      54                 : CPL_C_START
      55                 : void  GDALRegister_ISIS2(void);
      56                 : CPL_C_END
      57                 : 
      58                 : /************************************************************************/
      59                 : /* ==================================================================== */
      60                 : /*      ISISDataset version2                  */
      61                 : /* ==================================================================== */
      62                 : /************************************************************************/
      63                 : 
      64                 : class ISIS2Dataset : public RawDataset
      65                 : {
      66                 :     VSILFILE  *fpImage; // image data file.
      67                 :     CPLString    osExternalCube;
      68                 : 
      69                 :     NASAKeywordHandler  oKeywords;
      70                 :   
      71                 :     int         bGotTransform;
      72                 :     double      adfGeoTransform[6];
      73                 :   
      74                 :     CPLString   osProjection;
      75                 : 
      76                 :     int parse_label(const char *file, char *keyword, char *value);
      77                 :     int strstrip(char instr[], char outstr[], int position);
      78                 : 
      79                 :     CPLString   oTempResult;
      80                 : 
      81                 :     void        CleanString( CPLString &osInput );
      82                 : 
      83                 :     const char *GetKeyword( const char *pszPath, 
      84                 :                             const char *pszDefault = "");
      85                 :     const char *GetKeywordSub( const char *pszPath, 
      86                 :                                int iSubscript, 
      87                 :                                const char *pszDefault = "");
      88                 : 
      89                 : public:
      90                 :     ISIS2Dataset();
      91                 :     ~ISIS2Dataset();
      92                 : 
      93                 :     virtual CPLErr GetGeoTransform( double * padfTransform );
      94                 :     virtual const char *GetProjectionRef(void);
      95                 :   
      96                 :     virtual char **GetFileList();
      97                 : 
      98                 :     static int          Identify( GDALOpenInfo * );
      99                 :     static GDALDataset *Open( GDALOpenInfo * );
     100                 :     static GDALDataset *Create( const char * pszFilename,
     101                 :                                 int nXSize, int nYSize, int nBands,
     102                 :                                 GDALDataType eType, char ** papszParmList );
     103                 : 
     104                 :     // Write related.
     105                 :     static int WriteRaster(CPLString osFilename, bool includeLabel, GUIntBig iRecord, GUIntBig iLabelRecords, GDALDataType eType, const char * pszInterleaving);
     106                 : 
     107                 :     static int WriteLabel(CPLString osFilename, CPLString osRasterFile, CPLString sObjectTag, unsigned int nXSize, unsigned int nYSize, unsigned int nBands, GDALDataType eType,
     108                 :                 GUIntBig iRecords, const char * pszInterleaving, GUIntBig & iLabelRecords, bool bRelaunch=false);
     109                 :     static int WriteQUBE_Information(VSILFILE *fpLabel, unsigned int iLevel, unsigned int & nWritingBytes,
     110                 :                                      unsigned int nXSize, unsigned int nYSize, unsigned int nBands, GDALDataType eType, const char * pszInterleaving);
     111                 : 
     112                 :     static unsigned int WriteKeyword(VSILFILE *fpLabel, unsigned int iLevel, CPLString key, CPLString value);
     113                 :     static unsigned int WriteFormatting(VSILFILE *fpLabel, CPLString data);
     114                 :     static GUIntBig RecordSizeCalculation(unsigned int nXSize, unsigned int nYSize, unsigned int nBands, GDALDataType eType );
     115                 : };
     116                 : 
     117                 : /************************************************************************/
     118                 : /*                            ISIS2Dataset()                            */
     119                 : /************************************************************************/
     120                 : 
     121              44 : ISIS2Dataset::ISIS2Dataset()
     122                 : {
     123              44 :     fpImage = NULL;
     124              44 :     bGotTransform = FALSE;
     125              44 :     adfGeoTransform[0] = 0.0;
     126              44 :     adfGeoTransform[1] = 1.0;
     127              44 :     adfGeoTransform[2] = 0.0;
     128              44 :     adfGeoTransform[3] = 0.0;
     129              44 :     adfGeoTransform[4] = 0.0;
     130              44 :     adfGeoTransform[5] = 1.0;
     131              44 : }
     132                 : 
     133                 : /************************************************************************/
     134                 : /*                            ~ISIS2Dataset()                            */
     135                 : /************************************************************************/
     136                 : 
     137              44 : ISIS2Dataset::~ISIS2Dataset()
     138                 : 
     139                 : {
     140              44 :     FlushCache();
     141              44 :     if( fpImage != NULL )
     142              42 :         VSIFCloseL( fpImage );
     143              44 : }
     144                 : 
     145                 : /************************************************************************/
     146                 : /*                            GetFileList()                             */
     147                 : /************************************************************************/
     148                 : 
     149               2 : char **ISIS2Dataset::GetFileList()
     150                 : 
     151                 : {
     152               2 :     char **papszFileList = NULL;
     153                 : 
     154               2 :     papszFileList = GDALPamDataset::GetFileList();
     155                 : 
     156               2 :     if( strlen(osExternalCube) > 0 )
     157               1 :         papszFileList = CSLAddString( papszFileList, osExternalCube );
     158                 : 
     159               2 :     return papszFileList;
     160                 : }
     161                 : 
     162                 : /************************************************************************/
     163                 : /*                          GetProjectionRef()                          */
     164                 : /************************************************************************/
     165                 : 
     166               1 : const char *ISIS2Dataset::GetProjectionRef()
     167                 : 
     168                 : {
     169               1 :     if( strlen(osProjection) > 0 )
     170               1 :         return osProjection;
     171                 :     else
     172               0 :         return GDALPamDataset::GetProjectionRef();
     173                 : }
     174                 : 
     175                 : /************************************************************************/
     176                 : /*                          GetGeoTransform()                           */
     177                 : /************************************************************************/
     178                 : 
     179              14 : CPLErr ISIS2Dataset::GetGeoTransform( double * padfTransform )
     180                 : 
     181                 : {
     182              14 :     if( bGotTransform )
     183                 :     {
     184               1 :         memcpy( padfTransform, adfGeoTransform, sizeof(double) * 6 );
     185               1 :         return CE_None;
     186                 :     }
     187                 :     else
     188                 :     {
     189              13 :         return GDALPamDataset::GetGeoTransform( padfTransform );
     190                 :     }
     191                 : }
     192                 : 
     193                 : /************************************************************************/
     194                 : /*                              Identify()                              */
     195                 : /************************************************************************/
     196                 : 
     197                 : 
     198           11540 : int ISIS2Dataset::Identify( GDALOpenInfo * poOpenInfo )
     199                 : {
     200           11540 :     if( poOpenInfo->pabyHeader == NULL )
     201           10392 :         return FALSE;
     202                 : 
     203            1148 :     if( strstr((const char *)poOpenInfo->pabyHeader,"^QUBE") == NULL )
     204            1104 :         return FALSE;
     205                 :     else
     206              44 :         return TRUE;
     207                 : }
     208                 : 
     209                 : /************************************************************************/
     210                 : /*                                Open()                                */
     211                 : /************************************************************************/
     212                 : 
     213            2295 : GDALDataset *ISIS2Dataset::Open( GDALOpenInfo * poOpenInfo )
     214                 : {
     215                 : /* -------------------------------------------------------------------- */
     216                 : /*      Does this look like a CUBE or an IMAGE Primary Data Object?     */
     217                 : /* -------------------------------------------------------------------- */
     218            2295 :     if( !Identify( poOpenInfo ) )
     219            2251 :         return NULL;
     220                 : 
     221                 : /* -------------------------------------------------------------------- */
     222                 : /*      Open the file using the large file API.                         */
     223                 : /* -------------------------------------------------------------------- */
     224              44 :     VSILFILE *fpQube = VSIFOpenL( poOpenInfo->pszFilename, "rb" );
     225                 : 
     226              44 :     if( fpQube == NULL )
     227               0 :         return NULL;
     228                 : 
     229                 :     ISIS2Dataset  *poDS;
     230                 : 
     231              44 :     poDS = new ISIS2Dataset();
     232                 : 
     233              44 :     if( ! poDS->oKeywords.Ingest( fpQube, 0 ) )
     234                 :     {
     235               0 :         VSIFCloseL( fpQube );
     236               0 :         delete poDS;
     237               0 :         return NULL;
     238                 :     }
     239                 :     
     240              44 :     VSIFCloseL( fpQube );
     241                 : 
     242                 : /* -------------------------------------------------------------------- */
     243                 : /*  We assume the user is pointing to the label (ie. .lab) file.    */
     244                 : /* -------------------------------------------------------------------- */
     245                 :     // QUBE can be inline or detached and point to an image name
     246                 :     // ^QUBE = 76
     247                 :     // ^QUBE = ("ui31s015.img",6441<BYTES>) - has another label on the image
     248                 :     // ^QUBE = "ui31s015.img" - which implies no label or skip value
     249                 : 
     250              44 :     const char *pszQube = poDS->GetKeyword( "^QUBE" );
     251              44 :     GUIntBig nQube = 0;
     252              44 :     int bByteLocation = FALSE;
     253              44 :     CPLString osTargetFile = poOpenInfo->pszFilename;
     254                 : 
     255              44 :     if( pszQube[0] == '"' )
     256                 :     { 
     257               0 :         CPLString osTPath = CPLGetPath(poOpenInfo->pszFilename); 
     258               0 :         CPLString osFilename = pszQube;
     259               0 :         poDS->CleanString( osFilename ); 
     260               0 :         osTargetFile = CPLFormCIFilename( osTPath, osFilename, NULL ); 
     261               0 :         poDS->osExternalCube = osTargetFile;
     262                 :     }
     263              44 :     else if( pszQube[0] == '(' ) 
     264                 :     { 
     265               3 :         CPLString osTPath = CPLGetPath(poOpenInfo->pszFilename); 
     266               3 :         CPLString osFilename = poDS->GetKeywordSub("^QUBE",1,""); 
     267               3 :         poDS->CleanString( osFilename ); 
     268               3 :         osTargetFile = CPLFormCIFilename( osTPath, osFilename, NULL ); 
     269               3 :         poDS->osExternalCube = osTargetFile;
     270                 : 
     271               3 :         nQube = atoi(poDS->GetKeywordSub("^QUBE",2,"1"));
     272               3 :         if( strstr(poDS->GetKeywordSub("^QUBE",2,"1"),"<BYTES>") != NULL )
     273               0 :             bByteLocation = true;
     274                 :     } 
     275                 :     else 
     276                 :     { 
     277              41 :         nQube = atoi(pszQube);
     278              41 :         if( strstr(pszQube,"<BYTES>") != NULL )
     279               0 :             bByteLocation = true;
     280                 :     }
     281                 : 
     282                 : /* -------------------------------------------------------------------- */
     283                 : /*      Check if file an ISIS2 header file?  Read a few lines of text   */
     284                 : /*      searching for something starting with nrows or ncols.           */
     285                 : /* -------------------------------------------------------------------- */
     286              44 :     GDALDataType eDataType = GDT_Byte;
     287              44 :     OGRSpatialReference oSRS;
     288                 : 
     289                 :     //image parameters
     290              44 :     int nRows, nCols, nBands = 1;
     291              44 :     GUIntBig nSkipBytes = 0;
     292                 :     int itype;
     293                 :     int  s_ix, s_iy, s_iz; // check SUFFIX_ITEMS params.
     294                 :     int record_bytes;
     295              44 :     int bNoDataSet = FALSE;
     296              44 :     char chByteOrder = 'M';  //default to MSB
     297                 :  
     298                 :     //Georef parameters
     299              44 :     double dfULXMap=0.5;
     300              44 :     double dfULYMap = 0.5;
     301              44 :     double dfXDim = 1.0;
     302              44 :     double dfYDim = 1.0;
     303              44 :     double dfNoData = 0.0;
     304              44 :     double xulcenter = 0.0;
     305              44 :     double yulcenter = 0.0;
     306                 : 
     307                 :     //projection parameters
     308              44 :     int bProjectionSet = TRUE;
     309              44 :     double semi_major = 0.0;
     310              44 :     double semi_minor = 0.0;
     311              44 :     double iflattening = 0.0;
     312              44 :     double center_lat = 0.0;
     313              44 :     double center_lon = 0.0;
     314              44 :     double first_std_parallel = 0.0;
     315              44 :     double second_std_parallel = 0.0;
     316                 :     VSILFILE  *fp;
     317                 : 
     318                 :     /* -------------------------------------------------------------------- */
     319                 :     /*      Checks to see if this is valid ISIS2 cube                       */
     320                 :     /*      SUFFIX_ITEM tag in .cub file should be (0,0,0); no side-planes  */
     321                 :     /* -------------------------------------------------------------------- */
     322              44 :     s_ix = atoi(poDS->GetKeywordSub( "QUBE.SUFFIX_ITEMS", 1 ));
     323              44 :     s_iy = atoi(poDS->GetKeywordSub( "QUBE.SUFFIX_ITEMS", 2 ));
     324              44 :     s_iz = atoi(poDS->GetKeywordSub( "QUBE.SUFFIX_ITEMS", 3 ));
     325                 :      
     326              44 :     if( s_ix != 0 || s_iy != 0 || s_iz != 0 ) 
     327                 :     {
     328                 :         CPLError( CE_Failure, CPLE_OpenFailed, 
     329                 :                   "*** ISIS 2 cube file has invalid SUFFIX_ITEMS parameters:\n"
     330                 :                   "*** gdal isis2 driver requires (0, 0, 0), thus no sideplanes or backplanes\n"
     331               0 :                   "found: (%i, %i, %i)\n\n", s_ix, s_iy, s_iz );
     332               0 :         delete poDS;
     333               0 :         return NULL;
     334                 :     } 
     335                 : 
     336                 :     /**************** end SUFFIX_ITEM check ***********************/
     337                 :     
     338                 :     
     339                 :     /***********   Grab layout type (BSQ, BIP, BIL) ************/
     340                 :     //  AXIS_NAME = (SAMPLE,LINE,BAND)
     341                 :     /***********************************************************/
     342                 :     const char *value;
     343                 : 
     344              44 :     char szLayout[10] = "BSQ"; //default to band seq.
     345              44 :     value = poDS->GetKeyword( "QUBE.AXIS_NAME", "" );
     346              44 :     if (EQUAL(value,"(SAMPLE,LINE,BAND)") )
     347              44 :         strcpy(szLayout,"BSQ");
     348               0 :     else if (EQUAL(value,"(BAND,LINE,SAMPLE)") )
     349               0 :         strcpy(szLayout,"BIP");
     350               0 :     else if (EQUAL(value,"(SAMPLE,BAND,LINE)") || EQUAL(value,"") )
     351               0 :         strcpy(szLayout,"BSQ");
     352                 :     else {
     353                 :         CPLError( CE_Failure, CPLE_OpenFailed, 
     354               0 :                   "%s layout not supported. Abort\n\n", value);
     355               0 :         delete poDS;
     356               0 :         return NULL;
     357                 :     }
     358                 : 
     359                 :     /***********   Grab samples lines band ************/
     360              44 :     nCols = atoi(poDS->GetKeywordSub("QUBE.CORE_ITEMS",1));
     361              44 :     nRows = atoi(poDS->GetKeywordSub("QUBE.CORE_ITEMS",2));
     362              44 :     nBands = atoi(poDS->GetKeywordSub("QUBE.CORE_ITEMS",3));
     363                 :     
     364                 :     /***********   Grab Qube record bytes  **********/
     365              44 :     record_bytes = atoi(poDS->GetKeyword("RECORD_BYTES"));
     366                 : 
     367              44 :     if (nQube > 0 && bByteLocation )
     368               0 :         nSkipBytes = (nQube - 1);
     369              44 :     else if( nQube > 0 )
     370              44 :         nSkipBytes = (nQube - 1) * record_bytes;     
     371                 :     else
     372               0 :         nSkipBytes = 0;     
     373                 :      
     374                 :     /***********   Grab samples lines band ************/
     375              44 :     CPLString osCoreItemType = poDS->GetKeyword( "QUBE.CORE_ITEM_TYPE" );
     376              44 :     if( (EQUAL(osCoreItemType,"PC_INTEGER")) || 
     377                 :         (EQUAL(osCoreItemType,"PC_UNSIGNED_INTEGER")) || 
     378                 :         (EQUAL(osCoreItemType,"PC_REAL")) ) {
     379              43 :         chByteOrder = 'I';
     380                 :     }
     381                 :     
     382                 :     /********   Grab format type - isis2 only supports 8,16,32 *******/
     383              44 :     itype = atoi(poDS->GetKeyword("QUBE.CORE_ITEM_BYTES",""));
     384              44 :     switch(itype) {
     385                 :       case 1 :
     386              20 :         eDataType = GDT_Byte;
     387              20 :         dfNoData = NULL1;
     388              20 :         bNoDataSet = TRUE;
     389              20 :         break;
     390                 :       case 2 :
     391              10 :         if( strstr(osCoreItemType,"UNSIGNED") != NULL )
     392                 :         {
     393               5 :             dfNoData = 0;
     394               5 :             eDataType = GDT_UInt16;
     395                 :         }
     396                 :         else
     397                 :         {
     398               5 :             dfNoData = NULL2;
     399               5 :             eDataType = GDT_Int16;
     400                 :         }
     401              10 :         bNoDataSet = TRUE;
     402              10 :         break;
     403                 :       case 4 :
     404               9 :         eDataType = GDT_Float32;
     405               9 :         dfNoData = NULL3;
     406               9 :         bNoDataSet = TRUE;
     407               9 :         break;
     408                 :       case 8 :
     409               5 :         eDataType = GDT_Float64;
     410               5 :         dfNoData = NULL3;
     411               5 :         bNoDataSet = TRUE;
     412               5 :         break;
     413                 :       default :
     414                 :         CPLError( CE_Failure, CPLE_AppDefined,
     415                 :                   "Itype of %d is not supported in ISIS 2.",
     416               0 :                   itype); 
     417               0 :         delete poDS;
     418               0 :         return NULL;
     419                 :     }
     420                 : 
     421                 :     /***********   Grab Cellsize ************/
     422              44 :     value = poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.MAP_SCALE");
     423              44 :     if (strlen(value) > 0 ) {
     424               1 :         dfXDim = (float) atof(value) * 1000.0; /* convert from km to m */
     425               1 :         dfYDim = (float) atof(value) * 1000.0 * -1;
     426                 :     }
     427                 :     
     428                 :     /***********   Grab LINE_PROJECTION_OFFSET ************/
     429              44 :     value = poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.LINE_PROJECTION_OFFSET");
     430              44 :     if (strlen(value) > 0) {
     431               1 :         yulcenter = (float) atof(value);
     432               1 :         yulcenter = ((yulcenter) * dfYDim);
     433               1 :         dfULYMap = yulcenter - (dfYDim/2);
     434                 :     }
     435                 :      
     436                 :     /***********   Grab SAMPLE_PROJECTION_OFFSET ************/
     437              44 :     value = poDS->GetKeyword("QUBE.IMAGE_MAP_PROJECTION.SAMPLE_PROJECTION_OFFSET");
     438              44 :     if( strlen(value) > 0 ) {
     439               1 :         xulcenter = (float) atof(value);
     440               1 :         xulcenter = ((xulcenter) * dfXDim);
     441               1 :         dfULXMap = xulcenter - (dfXDim/2);
     442                 :     }
     443                 :      
     444                 :     /***********  Grab TARGET_NAME  ************/
     445                 :     /**** This is the planets name i.e. MARS ***/
     446              44 :     CPLString target_name = poDS->GetKeyword("QUBE.TARGET_NAME");
     447                 :      
     448                 :     /***********   Grab MAP_PROJECTION_TYPE ************/
     449                 :     CPLString map_proj_name = 
     450              44 :         poDS->GetKeyword( "QUBE.IMAGE_MAP_PROJECTION.MAP_PROJECTION_TYPE");
     451              44 :     poDS->CleanString( map_proj_name );
     452                 : 
     453                 :     /***********   Grab SEMI-MAJOR ************/
     454                 :     semi_major = 
     455              44 :         atof(poDS->GetKeyword( "QUBE.IMAGE_MAP_PROJECTION.A_AXIS_RADIUS")) * 1000.0;
     456                 : 
     457                 :     /***********   Grab semi-minor ************/
     458                 :     semi_minor = 
     459              44 :         atof(poDS->GetKeyword( "QUBE.IMAGE_MAP_PROJECTION.C_AXIS_RADIUS")) * 1000.0;
     460                 : 
     461                 :     /***********   Grab CENTER_LAT ************/
     462                 :     center_lat = 
     463              44 :         atof(poDS->GetKeyword( "QUBE.IMAGE_MAP_PROJECTION.CENTER_LATITUDE"));
     464                 : 
     465                 :     /***********   Grab CENTER_LON ************/
     466                 :     center_lon = 
     467              44 :         atof(poDS->GetKeyword( "QUBE.IMAGE_MAP_PROJECTION.CENTER_LONGITUDE"));
     468                 : 
     469                 :     /***********   Grab 1st std parallel ************/
     470                 :     first_std_parallel = 
     471              44 :         atof(poDS->GetKeyword( "QUBE.IMAGE_MAP_PROJECTION.FIRST_STANDARD_PARALLEL"));
     472                 : 
     473                 :     /***********   Grab 2nd std parallel ************/
     474                 :     second_std_parallel = 
     475              44 :         atof(poDS->GetKeyword( "QUBE.IMAGE_MAP_PROJECTION.SECOND_STANDARD_PARALLEL"));
     476                 :      
     477                 :     /*** grab  PROJECTION_LATITUDE_TYPE = "PLANETOCENTRIC" ****/
     478                 :     // Need to further study how ocentric/ographic will effect the gdal library.
     479                 :     // So far we will use this fact to define a sphere or ellipse for some projections
     480                 :     // Frank - may need to talk this over
     481              44 :     char bIsGeographic = TRUE;
     482              44 :     value = poDS->GetKeyword("CUBE.IMAGE_MAP_PROJECTION.PROJECTION_LATITUDE_TYPE");
     483              44 :     if (EQUAL( value, "\"PLANETOCENTRIC\"" ))
     484               0 :         bIsGeographic = FALSE; 
     485                 :      
     486              44 :     CPLDebug("ISIS2","using projection %s", map_proj_name.c_str() );
     487                 : 
     488                 :     //Set oSRS projection and parameters
     489              44 :     if ((EQUAL( map_proj_name, "EQUIRECTANGULAR_CYLINDRICAL" )) ||
     490                 :         (EQUAL( map_proj_name, "EQUIRECTANGULAR" )) ||
     491                 :         (EQUAL( map_proj_name, "SIMPLE_CYLINDRICAL" )) ) {
     492               1 :         oSRS.OGRSpatialReference::SetEquirectangular2 ( 0.0, center_lon, center_lat, 0, 0 );
     493              43 :     } else if (EQUAL( map_proj_name, "ORTHOGRAPHIC" )) { 
     494               0 :         oSRS.OGRSpatialReference::SetOrthographic ( center_lat, center_lon, 0, 0 );
     495              43 :     } else if ((EQUAL( map_proj_name, "SINUSOIDAL" )) ||
     496                 :                (EQUAL( map_proj_name, "SINUSOIDAL_EQUAL-AREA" ))) {
     497               0 :         oSRS.OGRSpatialReference::SetSinusoidal ( center_lon, 0, 0 );
     498              43 :     } else if (EQUAL( map_proj_name, "MERCATOR" )) {
     499               0 :         oSRS.OGRSpatialReference::SetMercator ( center_lat, center_lon, 1, 0, 0 );
     500              43 :     } else if (EQUAL( map_proj_name, "POLAR_STEREOGRAPHIC" )) {
     501               0 :         oSRS.OGRSpatialReference::SetPS ( center_lat, center_lon, 1, 0, 0 );
     502              43 :     } else if (EQUAL( map_proj_name, "TRANSVERSE_MERCATOR" )) {
     503               0 :         oSRS.OGRSpatialReference::SetTM ( center_lat, center_lon, 1, 0, 0 );
     504              43 :     } else if (EQUAL( map_proj_name, "LAMBERT_CONFORMAL_CONIC" )) {
     505               0 :         oSRS.OGRSpatialReference::SetLCC ( first_std_parallel, second_std_parallel, center_lat, center_lon, 0, 0 );
     506              43 :     } else if (EQUAL( map_proj_name, "") ) {
     507                 :         /* no projection */
     508              43 :         bProjectionSet = FALSE;
     509                 :     } else {
     510                 :         CPLDebug( "ISIS2",
     511                 :                   "Dataset projection %s is not supported. Continuing...",
     512               0 :                   map_proj_name.c_str() );
     513               0 :         bProjectionSet = FALSE;
     514                 :     }
     515                 : 
     516              44 :     if (bProjectionSet) {
     517                 :         //Create projection name, i.e. MERCATOR MARS and set as ProjCS keyword
     518               1 :         CPLString proj_target_name = map_proj_name + " " + target_name;
     519               1 :         oSRS.SetProjCS(proj_target_name); //set ProjCS keyword
     520                 :      
     521                 :         //The geographic/geocentric name will be the same basic name as the body name
     522                 :         //'GCS' = Geographic/Geocentric Coordinate System
     523               1 :         CPLString geog_name = "GCS_" + target_name;
     524                 :         
     525                 :         //The datum and sphere names will be the same basic name aas the planet
     526               1 :         CPLString datum_name = "D_" + target_name;
     527               1 :         CPLString sphere_name = target_name; // + "_IAU_IAG");  //Might not be IAU defined so don't add
     528                 :           
     529                 :         //calculate inverse flattening from major and minor axis: 1/f = a/(a-b)
     530               1 :         if ((semi_major - semi_minor) < 0.0000001) 
     531               1 :             iflattening = 0;
     532                 :         else
     533               0 :             iflattening = semi_major / (semi_major - semi_minor);
     534                 :      
     535                 :         //Set the body size but take into consideration which proj is being used to help w/ proj4 compatibility
     536                 :         //The use of a Sphere, polar radius or ellipse here is based on how ISIS does it internally
     537               1 :         if ( ( (EQUAL( map_proj_name, "STEREOGRAPHIC" ) && (fabs(center_lat) == 90)) ) || 
     538                 :              (EQUAL( map_proj_name, "POLAR_STEREOGRAPHIC" )))  
     539                 :         {
     540               0 :             if (bIsGeographic) { 
     541                 :                 //Geograpraphic, so set an ellipse
     542                 :                 oSRS.SetGeogCS( geog_name, datum_name, sphere_name,
     543                 :                                 semi_major, iflattening, 
     544               0 :                                 "Reference_Meridian", 0.0 );
     545                 :             } else {
     546                 :                 //Geocentric, so force a sphere using the semi-minor axis. I hope... 
     547               0 :                 sphere_name += "_polarRadius";
     548                 :                 oSRS.SetGeogCS( geog_name, datum_name, sphere_name,
     549                 :                                 semi_minor, 0.0, 
     550               0 :                                 "Reference_Meridian", 0.0 );
     551                 :             }
     552                 :         }
     553               1 :         else if ( (EQUAL( map_proj_name, "SIMPLE_CYLINDRICAL" )) || 
     554                 :                   (EQUAL( map_proj_name, "ORTHOGRAPHIC" )) ||
     555                 :                   (EQUAL( map_proj_name, "STEREOGRAPHIC" )) ||
     556                 :                   (EQUAL( map_proj_name, "SINUSOIDAL_EQUAL-AREA" )) ||
     557                 :                   (EQUAL( map_proj_name, "SINUSOIDAL" ))  ) {
     558                 :             //isis uses the sphereical equation for these projections so force a sphere
     559                 :             oSRS.SetGeogCS( geog_name, datum_name, sphere_name,
     560                 :                             semi_major, 0.0, 
     561               1 :                             "Reference_Meridian", 0.0 );
     562                 :         } 
     563               0 :         else if  ((EQUAL( map_proj_name, "EQUIRECTANGULAR_CYLINDRICAL" )) || 
     564                 :                   (EQUAL( map_proj_name, "EQUIRECTANGULAR" )) ) {
     565                 :             //Calculate localRadius using ISIS3 simple elliptical method 
     566                 :             //  not the more standard Radius of Curvature method
     567                 :             //PI = 4 * atan(1);
     568                 :             double radLat, localRadius;
     569               0 :             if (center_lon == 0) { //No need to calculate local radius
     570                 :                 oSRS.SetGeogCS( geog_name, datum_name, sphere_name,
     571                 :                                 semi_major, 0.0, 
     572               0 :                                 "Reference_Meridian", 0.0 );
     573                 :             } else {  
     574               0 :                 radLat = center_lat * PI / 180;  // in radians
     575                 :                 localRadius = semi_major * semi_minor / sqrt(pow(semi_minor*cos(radLat),2)
     576               0 :                                                              + pow(semi_major*sin(radLat),2) );
     577               0 :                 sphere_name += "_localRadius";
     578                 :                 oSRS.SetGeogCS( geog_name, datum_name, sphere_name,
     579                 :                                 localRadius, 0.0, 
     580               0 :                                 "Reference_Meridian", 0.0 );
     581               0 :                 CPLDebug( "ISIS2", "local radius: %f", localRadius);
     582                 :             }
     583                 :         } 
     584                 :         else { 
     585                 :             //All other projections: Mercator, Transverse Mercator, Lambert Conformal, etc.
     586                 :             //Geographic, so set an ellipse
     587               0 :             if (bIsGeographic) {
     588                 :                 oSRS.SetGeogCS( geog_name, datum_name, sphere_name,
     589                 :                                 semi_major, iflattening, 
     590               0 :                                 "Reference_Meridian", 0.0 );
     591                 :             } else { 
     592                 :                 //Geocentric, so force a sphere. I hope... 
     593                 :                 oSRS.SetGeogCS( geog_name, datum_name, sphere_name,
     594                 :                                 semi_major, 0.0, 
     595               0 :                                 "Reference_Meridian", 0.0 );
     596                 :             }
     597                 :         }
     598                 :         
     599                 : 
     600                 :         // translate back into a projection string.
     601               1 :         char *pszResult = NULL;
     602               1 :         oSRS.exportToWkt( &pszResult );
     603               1 :         poDS->osProjection = pszResult;
     604               1 :         CPLFree( pszResult );
     605                 :     }
     606                 : 
     607                 : /* END ISIS2 Label Read */
     608                 : /*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*/
     609                 :     
     610                 : /* -------------------------------------------------------------------- */
     611                 : /*      Did we get the required keywords?  If not we return with        */
     612                 : /*      this never having been considered to be a match. This isn't     */
     613                 : /*      an error!                                                       */
     614                 : /* -------------------------------------------------------------------- */
     615              44 :     if( nRows < 1 || nCols < 1 || nBands < 1 )
     616                 :     {
     617               2 :         delete poDS;
     618               2 :         return NULL;
     619                 :     }
     620                 : 
     621                 : /* -------------------------------------------------------------------- */
     622                 : /*      Capture some information from the file that is of interest.     */
     623                 : /* -------------------------------------------------------------------- */
     624              42 :     poDS->nRasterXSize = nCols;
     625              42 :     poDS->nRasterYSize = nRows;
     626                 : 
     627                 : /* -------------------------------------------------------------------- */
     628                 : /*      Open target binary file.                                        */
     629                 : /* -------------------------------------------------------------------- */
     630                 :     
     631              42 :     if( poOpenInfo->eAccess == GA_ReadOnly )
     632              18 :         poDS->fpImage = VSIFOpenL( osTargetFile, "rb" );
     633                 :     else
     634              24 :         poDS->fpImage = VSIFOpenL( osTargetFile, "r+b" );
     635                 : 
     636              42 :     if( poDS->fpImage == NULL )
     637                 :     {
     638                 :         CPLError( CE_Failure, CPLE_OpenFailed, 
     639                 :                   "Failed to open %s with write permission.\n%s", 
     640               0 :                   osTargetFile.c_str(), VSIStrerror( errno ) );
     641               0 :         delete poDS;
     642               0 :         return NULL;
     643                 :     }
     644                 : 
     645              42 :     poDS->eAccess = poOpenInfo->eAccess;
     646                 : 
     647                 : /* -------------------------------------------------------------------- */
     648                 : /*      Compute the line offset.                                        */
     649                 : /* -------------------------------------------------------------------- */
     650              42 :     int     nItemSize = GDALGetDataTypeSize(eDataType)/8;
     651                 :     int   nLineOffset, nPixelOffset, nBandOffset;
     652                 :     
     653              42 :     if( EQUAL(szLayout,"BIP") )
     654                 :     {
     655               0 :         nPixelOffset = nItemSize * nBands;
     656               0 :         nLineOffset = nPixelOffset * nCols;
     657               0 :         nBandOffset = nItemSize;
     658                 :     }
     659              42 :     else if( EQUAL(szLayout,"BSQ") )
     660                 :     {
     661              42 :         nPixelOffset = nItemSize;
     662              42 :         nLineOffset = nPixelOffset * nCols;
     663              42 :         nBandOffset = nLineOffset * nRows;
     664                 :     }
     665                 :     else /* assume BIL */
     666                 :     {
     667               0 :         nPixelOffset = nItemSize;
     668               0 :         nLineOffset = nItemSize * nBands * nCols;
     669               0 :         nBandOffset = nItemSize * nCols;
     670                 :     }
     671                 :     
     672                 : /* -------------------------------------------------------------------- */
     673                 : /*      Create band information objects.                                */
     674                 : /* -------------------------------------------------------------------- */
     675                 :     int i;
     676                 : 
     677              42 :     poDS->nBands = nBands;;
     678             130 :     for( i = 0; i < poDS->nBands; i++ )
     679                 :     {
     680                 :         RawRasterBand *poBand;
     681                 : 
     682                 :         poBand = 
     683                 :             new RawRasterBand( poDS, i+1, poDS->fpImage,
     684                 :                                nSkipBytes + nBandOffset * i, 
     685                 :                                nPixelOffset, nLineOffset, eDataType,
     686                 : #ifdef CPL_LSB                               
     687                 :                                chByteOrder == 'I' || chByteOrder == 'L',
     688                 : #else
     689                 :                                chByteOrder == 'M',
     690                 : #endif        
     691              88 :                                TRUE );
     692                 : 
     693              88 :         if( bNoDataSet )
     694              88 :             poBand->SetNoDataValue( dfNoData );
     695                 : 
     696              88 :         poDS->SetBand( i+1, poBand );
     697                 : 
     698                 :         // Set offset/scale values at the PAM level.
     699                 :         poBand->SetOffset( 
     700              88 :             CPLAtofM(poDS->GetKeyword("QUBE.CORE_BASE","0.0")));
     701                 :         poBand->SetScale( 
     702              88 :             CPLAtofM(poDS->GetKeyword("QUBE.CORE_MULTIPLIER","1.0")));
     703                 :     }
     704                 : 
     705                 : /* -------------------------------------------------------------------- */
     706                 : /*      Check for a .prj file. For isis2 I would like to keep this in   */
     707                 : /* -------------------------------------------------------------------- */
     708              42 :     CPLString osPath, osName;
     709                 : 
     710              42 :     osPath = CPLGetPath( poOpenInfo->pszFilename );
     711              42 :     osName = CPLGetBasename(poOpenInfo->pszFilename);
     712              42 :     const char  *pszPrjFile = CPLFormCIFilename( osPath, osName, "prj" );
     713                 : 
     714              42 :     fp = VSIFOpenL( pszPrjFile, "r" );
     715              42 :     if( fp != NULL )
     716                 :     {
     717                 :         char  **papszLines;
     718               0 :         OGRSpatialReference oSRS;
     719                 : 
     720               0 :         VSIFCloseL( fp );
     721                 :         
     722               0 :         papszLines = CSLLoad( pszPrjFile );
     723                 : 
     724               0 :         if( oSRS.importFromESRI( papszLines ) == OGRERR_NONE )
     725                 :         {
     726               0 :             char *pszResult = NULL;
     727               0 :             oSRS.exportToWkt( &pszResult );
     728               0 :             poDS->osProjection = pszResult;
     729               0 :             CPLFree( pszResult );
     730                 :         }
     731                 : 
     732               0 :         CSLDestroy( papszLines );
     733                 :     }
     734                 : 
     735                 :     
     736              42 :     if( dfULYMap != 0.5 || dfULYMap != 0.5 || dfXDim != 1.0 || dfYDim != 1.0 )
     737                 :     {
     738               1 :         poDS->bGotTransform = TRUE;
     739               1 :         poDS->adfGeoTransform[0] = dfULXMap;
     740               1 :         poDS->adfGeoTransform[1] = dfXDim;
     741               1 :         poDS->adfGeoTransform[2] = 0.0;
     742               1 :         poDS->adfGeoTransform[3] = dfULYMap;
     743               1 :         poDS->adfGeoTransform[4] = 0.0;
     744               1 :         poDS->adfGeoTransform[5] = dfYDim;
     745                 :     }
     746                 :     
     747              42 :     if( !poDS->bGotTransform )
     748                 :         poDS->bGotTransform = 
     749                 :             GDALReadWorldFile( poOpenInfo->pszFilename, "cbw", 
     750              41 :                                poDS->adfGeoTransform );
     751                 : 
     752              42 :     if( !poDS->bGotTransform )
     753                 :         poDS->bGotTransform = 
     754                 :             GDALReadWorldFile( poOpenInfo->pszFilename, "wld", 
     755              41 :                                poDS->adfGeoTransform );
     756                 : 
     757                 : /* -------------------------------------------------------------------- */ 
     758                 : /*      Initialize any PAM information.                                 */ 
     759                 : /* -------------------------------------------------------------------- */ 
     760              42 :     poDS->SetDescription( poOpenInfo->pszFilename ); 
     761              42 :     poDS->TryLoadXML(); 
     762                 : 
     763                 : /* -------------------------------------------------------------------- */
     764                 : /*      Check for overviews.                                            */
     765                 : /* -------------------------------------------------------------------- */
     766              42 :     poDS->oOvManager.Initialize( poDS, poOpenInfo->pszFilename );
     767                 : 
     768              42 :     return( poDS );
     769                 : }
     770                 : 
     771                 : /************************************************************************/
     772                 : /*                             GetKeyword()                             */
     773                 : /************************************************************************/
     774                 : 
     775             924 : const char *ISIS2Dataset::GetKeyword( const char *pszPath, 
     776                 :                                       const char *pszDefault )
     777                 : 
     778                 : {
     779             924 :     return oKeywords.GetKeyword( pszPath, pszDefault );
     780                 : }
     781                 : 
     782                 : /************************************************************************/
     783                 : /*                            GetKeywordSub()                           */
     784                 : /************************************************************************/
     785                 : 
     786             273 : const char *ISIS2Dataset::GetKeywordSub( const char *pszPath, 
     787                 :                                          int iSubscript,
     788                 :                                          const char *pszDefault )
     789                 : 
     790                 : {
     791             273 :     const char *pszResult = oKeywords.GetKeyword( pszPath, NULL );
     792                 :     
     793             273 :     if( pszResult == NULL )
     794               0 :         return pszDefault;
     795                 : 
     796             273 :     if( pszResult[0] != '(' )
     797               0 :         return pszDefault;
     798                 : 
     799                 :     char **papszTokens = CSLTokenizeString2( pszResult, "(,)", 
     800             273 :                                              CSLT_HONOURSTRINGS );
     801                 : 
     802             273 :     if( iSubscript <= CSLCount(papszTokens) )
     803                 :     {
     804             273 :         oTempResult = papszTokens[iSubscript-1];
     805             273 :         CSLDestroy( papszTokens );
     806             273 :         return oTempResult.c_str();
     807                 :     }
     808                 :     else
     809                 :     {
     810               0 :         CSLDestroy( papszTokens );
     811               0 :         return pszDefault;
     812                 :     }
     813                 : }
     814                 : 
     815                 : /************************************************************************/
     816                 : /*                            CleanString()                             */
     817                 : /*                                                                      */
     818                 : /* Removes single or double quotes, and converts spaces to underscores. */
     819                 : /* The change is made in-place to CPLString.                            */
     820                 : /************************************************************************/
     821                 : 
     822              47 : void ISIS2Dataset::CleanString( CPLString &osInput )
     823                 : 
     824                 : {
     825              47 :    if(  ( osInput.size() < 2 ) ||
     826                 :         ((osInput.at(0) != '"'   || osInput.at(osInput.size()-1) != '"' ) &&
     827                 :         ( osInput.at(0) != '\'' || osInput.at(osInput.size()-1) != '\'')) )
     828              47 :         return;
     829                 : 
     830               0 :     char *pszWrk = CPLStrdup(osInput.c_str() + 1);
     831                 :     int i;
     832                 : 
     833               0 :     pszWrk[strlen(pszWrk)-1] = '\0';
     834                 :     
     835               0 :     for( i = 0; pszWrk[i] != '\0'; i++ )
     836                 :     {
     837               0 :         if( pszWrk[i] == ' ' )
     838               0 :             pszWrk[i] = '_';
     839                 :     }
     840                 : 
     841               0 :     osInput = pszWrk;
     842               0 :     CPLFree( pszWrk );
     843                 : }
     844                 : 
     845                 : /************************************************************************/
     846                 : /*                           Create()                                   */
     847                 : /************************************************************************/
     848                 : /**
     849                 :  * Hidden Creation Options:
     850                 :  * INTERLEAVE=BSQ/BIP/BIL: Force the generation specified type of interleaving.
     851                 :  *  BSQ --- band sequental (default),
     852                 :  *  BIP --- band interleaved by pixel,
     853                 :  *  BIL --- band interleaved by line.
     854                 :  * OBJECT=QUBE/IMAGE/SPECTRAL_QUBE, if null default is QUBE
     855                 :  */
     856                 : 
     857              45 : GDALDataset *ISIS2Dataset::Create(const char* pszFilename,
     858                 :                                   int nXSize, int nYSize, int nBands,
     859                 :                                   GDALDataType eType, char** papszParmList) {
     860                 : 
     861                 :     /* Verify settings. In Isis 2 core pixel values can be represented in
     862                 :      * three different ways : 1, 2 4, or 8 Bytes */
     863              45 :     if( eType != GDT_Byte && eType != GDT_Int16 && eType != GDT_Float32
     864                 :         && eType != GDT_UInt16 && eType != GDT_Float64 ){
     865                 :         CPLError(CE_Failure, CPLE_AppDefined,
     866                 :                  "The ISIS2 driver does not supporting creating files of type %s.",
     867              18 :                  GDALGetDataTypeName( eType ) );
     868              18 :         return NULL;
     869                 :     }
     870                 : 
     871                 :     /*  (SAMPLE, LINE, BAND) - Band Sequential (BSQ) - default choice
     872                 :         (SAMPLE, BAND, LINE) - Band Interleaved by Line (BIL)
     873                 :         (BAND, SAMPLE, LINE) - Band Interleaved by Pixel (BIP) */
     874              27 :     const char *pszInterleaving = "(SAMPLE,LINE,BAND)";
     875              27 :     const char *pszInterleavingParam = CSLFetchNameValue( papszParmList, "INTERLEAVE" );
     876              27 :     if ( pszInterleavingParam ) {
     877               0 :         if ( EQUALN( pszInterleavingParam, "bip", 3 ) )
     878               0 :             pszInterleaving = "(BAND,SAMPLE,LINE)";
     879               0 :         else if ( EQUALN( pszInterleavingParam, "bil", 3 ) )
     880               0 :             pszInterleaving = "(SAMPLE,BAND,LINE)";
     881                 :         else
     882               0 :             pszInterleaving = "(SAMPLE,LINE,BAND)";
     883                 :     }
     884                 : 
     885                 :     /* default labeling method is attached */
     886              27 :     bool bAttachedLabelingMethod = true;
     887                 :     /* check if labeling method is set : check the all three first chars */
     888              27 :     const char *pszLabelingMethod = CSLFetchNameValue( papszParmList, "LABELING_METHOD" );
     889              27 :     if ( pszLabelingMethod ){
     890               1 :         if ( EQUALN( pszLabelingMethod, "detached", 3 ) ){
     891               1 :             bAttachedLabelingMethod = false;
     892                 :         }
     893               1 :         if ( EQUALN( pszLabelingMethod, "attached", 3 ) ){
     894               0 :             bAttachedLabelingMethod = true;
     895                 :         }
     896                 :     }
     897                 : 
     898                 :     /*  set the label and data files */
     899              27 :     CPLString osLabelFile, osRasterFile, osOutFile;
     900              27 :     if( bAttachedLabelingMethod ) {
     901              26 :         osLabelFile = "";
     902              26 :         osRasterFile = pszFilename;
     903              26 :         osOutFile = osRasterFile;
     904                 :     }
     905                 :     else
     906                 :     {
     907               1 :         CPLString sExtension = "cub";
     908               1 :         const char* pszExtension = CSLFetchNameValue( papszParmList, "IMAGE_EXTENSION" );
     909               1 :         if( pszExtension ){
     910               1 :             sExtension = pszExtension;
     911                 :         }
     912                 : 
     913               1 :         if( EQUAL(CPLGetExtension( pszFilename ), sExtension) )
     914                 :         {
     915                 :             CPLError( CE_Failure, CPLE_AppDefined,
     916                 :                       "IMAGE_EXTENSION (%s) cannot match LABEL file extension.",
     917               0 :                       sExtension.c_str() );
     918               0 :             return NULL;
     919                 :         }
     920                 : 
     921               1 :         osLabelFile = pszFilename;
     922               1 :         osRasterFile = CPLResetExtension( osLabelFile, sExtension );
     923               1 :         osOutFile = osLabelFile;
     924                 :     }
     925                 : 
     926              27 :     const char *pszObject = CSLFetchNameValue( papszParmList, "OBJECT" );
     927              27 :     CPLString sObject = "QUBE"; // default choice
     928              27 :     if (pszObject) {
     929               0 :         if ( EQUAL( pszObject, "IMAGE") ){
     930               0 :             sObject = "IMAGE";
     931                 :         }
     932               0 :         if ( EQUAL( pszObject, "SPECTRAL_QUBE")){
     933               0 :             sObject = "SPECTRAL_QUBE";
     934                 :         }
     935                 :     }
     936                 : 
     937              27 :     GUIntBig iRecords = ISIS2Dataset::RecordSizeCalculation(nXSize, nYSize, nBands, eType);
     938              27 :     GUIntBig iLabelRecords(2);
     939                 : 
     940              27 :     CPLDebug("ISIS2","irecord = %i",static_cast<int>(iRecords));
     941                 : 
     942              27 :     if( bAttachedLabelingMethod ) 
     943                 :     {
     944              26 :         ISIS2Dataset::WriteLabel(osRasterFile, "", sObject, nXSize, nYSize, nBands, eType, iRecords, pszInterleaving, iLabelRecords, true);
     945                 :     }
     946                 :     else
     947                 :     {
     948               1 :         ISIS2Dataset::WriteLabel(osLabelFile, osRasterFile, sObject, nXSize, nYSize, nBands, eType, iRecords, pszInterleaving, iLabelRecords);
     949                 :     }
     950                 : 
     951              27 :     if( !ISIS2Dataset::WriteRaster(osRasterFile, bAttachedLabelingMethod, 
     952                 :                                    iRecords, iLabelRecords, eType, 
     953                 :                                    pszInterleaving) )
     954               2 :         return NULL;
     955                 : 
     956              25 :     return (GDALDataset *) GDALOpen( osOutFile, GA_Update );
     957                 : }
     958                 : 
     959                 : 
     960                 : /************************************************************************/
     961                 : /*                            WriteRaster()                             */
     962                 : /************************************************************************/
     963                 : 
     964              27 : int ISIS2Dataset::WriteRaster(CPLString osFilename, bool includeLabel, GUIntBig iRecords, GUIntBig iLabelRecords, GDALDataType eType, const char * pszInterleaving) {
     965                 :     GUIntBig nSize;
     966              27 :     GByte byZero(0);
     967              27 :     CPLString pszAccess("wb");
     968              27 :     if(includeLabel)
     969              26 :         pszAccess = "ab";
     970                 : 
     971              27 :     VSILFILE *fpBin = VSIFOpenL( osFilename, pszAccess.c_str() );
     972              27 :     if( fpBin == NULL ) {
     973                 :         CPLError( CE_Failure, CPLE_FileIO,
     974                 :                   "Failed to create %s:\n%s",
     975               2 :                   osFilename.c_str(), VSIStrerror( errno ) );
     976               2 :         return FALSE;
     977                 :     }
     978                 : 
     979              25 :     nSize = iRecords * RECORD_SIZE;
     980              25 :     CPLDebug("ISIS2","nSize = %i", static_cast<int>(nSize));
     981                 : 
     982              25 :     if(includeLabel)
     983              24 :         nSize = iLabelRecords * RECORD_SIZE + nSize;
     984                 : 
     985                 :     // write last byte
     986              25 :     if(VSIFSeekL( fpBin, nSize-1, SEEK_SET ) != 0 ||
     987                 :        VSIFWriteL( &byZero, 1, 1, fpBin ) != 1){
     988                 :         CPLError( CE_Failure, CPLE_FileIO,
     989                 :                   "Failed to write %s:\n%s",
     990               0 :                   osFilename.c_str(), VSIStrerror( errno ) );
     991               0 :         VSIFCloseL( fpBin );
     992               0 :         return FALSE;
     993                 :     }
     994              25 :     VSIFCloseL( fpBin );
     995                 : 
     996              25 :     return TRUE;
     997                 : }
     998                 : 
     999                 : /************************************************************************/
    1000                 : /*                       RecordSizeCalculation()                        */
    1001                 : /************************************************************************/
    1002              27 : GUIntBig ISIS2Dataset::RecordSizeCalculation(unsigned int nXSize, unsigned int nYSize, unsigned int nBands, GDALDataType eType )
    1003                 : 
    1004                 : {
    1005              27 :     GUIntBig n = nXSize * nYSize * nBands * (  GDALGetDataTypeSize(eType) / 8);
    1006                 :     // size of pds file is a multiple of RECORD_SIZE Bytes.
    1007              27 :     CPLDebug("ISIS2","n = %i", static_cast<int>(n));
    1008              27 :     CPLDebug("ISIS2","RECORD SIZE = %i", RECORD_SIZE);
    1009              27 :     CPLDebug("ISIS2","nXSize = %i", nXSize);
    1010              27 :     CPLDebug("ISIS2","nYSize = %i", nYSize);
    1011              27 :     CPLDebug("ISIS2","nBands = %i", nBands);
    1012              27 :     CPLDebug("ISIS2","DataTypeSize = %i", GDALGetDataTypeSize(eType));
    1013              27 :     return (GUIntBig) ceil(static_cast<float>(n)/RECORD_SIZE);
    1014                 : }
    1015                 : 
    1016                 : /************************************************************************/
    1017                 : /*                       WriteQUBE_Information()                        */
    1018                 : /************************************************************************/
    1019                 : 
    1020              25 : int ISIS2Dataset::WriteQUBE_Information(
    1021                 :     VSILFILE *fpLabel, unsigned int iLevel, unsigned int & nWritingBytes,
    1022                 :     unsigned int nXSize,  unsigned int nYSize, unsigned int nBands, 
    1023                 :     GDALDataType eType, const char * pszInterleaving)
    1024                 : 
    1025                 : {
    1026              25 :     nWritingBytes += ISIS2Dataset::WriteFormatting( fpLabel, "");
    1027              50 :     nWritingBytes += ISIS2Dataset::WriteFormatting( fpLabel, "/* Qube structure */");
    1028              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "OBJECT", "QUBE");
    1029              25 :     iLevel++;
    1030              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "AXES", "3");
    1031              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "AXIS_NAME", pszInterleaving);
    1032              50 :     nWritingBytes += ISIS2Dataset::WriteFormatting( fpLabel, "/* Core description */");
    1033                 : 
    1034              25 :     CPLDebug("ISIS2","%d,%d,%d",nXSize,nYSize,nBands);
    1035                 : 
    1036              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_ITEMS",CPLString().Printf("(%d,%d,%d)",nXSize,nYSize,nBands));
    1037              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_NAME", "\"RAW DATA NUMBER\"");
    1038              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_UNIT", "\"N/A\"");
    1039                 :     // TODO change for eType
    1040                 : 
    1041              25 :     if( eType == GDT_Byte ) 
    1042                 :     {
    1043              12 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_ITEM_TYPE", "PC_UNSIGNED_INTEGER");
    1044              24 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_ITEM_BYTES", "1");
    1045                 :     }
    1046              13 :     else if( eType == GDT_UInt16 ) 
    1047                 :     {
    1048               3 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_ITEM_TYPE", "PC_UNSIGNED_INTEGER");
    1049               6 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_ITEM_BYTES", "2");
    1050                 :     }
    1051              10 :     else if( eType == GDT_Int16 ) 
    1052                 :     {
    1053               3 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_ITEM_TYPE", "PC_INTEGER");
    1054               6 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_ITEM_BYTES", "2");
    1055                 :     }
    1056               7 :     else if( eType == GDT_Float32 )
    1057                 :     {
    1058               4 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_ITEM_TYPE", "PC_REAL");
    1059               8 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_ITEM_BYTES", "4");
    1060                 :     }
    1061               3 :     else if( eType == GDT_Float64 )
    1062                 :     {
    1063               3 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_ITEM_TYPE", "PC_REAL");
    1064               6 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_ITEM_BYTES", "8");
    1065                 :     }
    1066                 : 
    1067                 :     // TODO add core null value
    1068                 : 
    1069              25 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_BASE", "0.0");
    1070              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "CORE_MULTIPLIER", "1.0");
    1071              50 :     nWritingBytes += ISIS2Dataset::WriteFormatting( fpLabel, "/* Suffix description */");
    1072              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "SUFFIX_BYTES", "4");
    1073              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "SUFFIX_ITEMS", "( 0, 0, 0)");
    1074              25 :     iLevel--;
    1075              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "END_OBJECT", "QUBE");
    1076                 :     
    1077              25 :     return TRUE;
    1078                 : }
    1079                 : 
    1080                 : /************************************************************************/
    1081                 : /*                             WriteLabel()                             */
    1082                 : /*                                                                      */
    1083                 : /*      osRasterFile : name of raster file but if it is empty we        */
    1084                 : /*                     have only one file with an attached label        */
    1085                 : /*      sObjectTag : QUBE, IMAGE or SPECTRAL_QUBE                       */
    1086                 : /*      bRelaunch : flag to allow recursiv call                         */
    1087                 : /************************************************************************/
    1088                 : 
    1089              27 : int ISIS2Dataset::WriteLabel(
    1090                 :     CPLString osFilename, CPLString osRasterFile, CPLString sObjectTag,
    1091                 :     unsigned int nXSize, unsigned int nYSize, unsigned int nBands, 
    1092                 :     GDALDataType eType,
    1093                 :     GUIntBig iRecords, const char * pszInterleaving, 
    1094                 :     GUIntBig &iLabelRecords, bool bRelaunch)
    1095                 : 
    1096                 : {
    1097              27 :     CPLDebug("ISIS2", "Write Label filename = %s, rasterfile = %s",osFilename.c_str(),osRasterFile.c_str());
    1098              27 :     bool bAttachedLabel = EQUAL(osRasterFile, "");
    1099                 : 
    1100              27 :     VSILFILE *fpLabel = VSIFOpenL( osFilename, "w" );
    1101                 : 
    1102              27 :     if( fpLabel == NULL ){
    1103                 :         CPLError( CE_Failure, CPLE_FileIO,
    1104                 :                   "Failed to create %s:\n%s",
    1105               2 :                   osFilename.c_str(), VSIStrerror( errno ) );
    1106               2 :         return FALSE;
    1107                 :     }
    1108                 : 
    1109              25 :     unsigned int iLevel(0);
    1110              25 :     unsigned int nWritingBytes(0);
    1111                 : 
    1112                 :     /* write common header */
    1113              25 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "PDS_VERSION_ID", "PDS3" );
    1114              50 :     nWritingBytes += ISIS2Dataset::WriteFormatting( fpLabel, "");
    1115              50 :     nWritingBytes += ISIS2Dataset::WriteFormatting( fpLabel, "/* File identification and structure */");
    1116              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "RECORD_TYPE", "FIXED_LENGTH" );
    1117              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "RECORD_BYTES", CPLString().Printf("%d",RECORD_SIZE));
    1118              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "FILE_RECORDS", CPLString().Printf(CPL_FRMT_GUIB,iRecords));
    1119              50 :     nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "LABEL_RECORDS", CPLString().Printf(CPL_FRMT_GUIB,iLabelRecords));
    1120              25 :     if(!bAttachedLabel){
    1121               1 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, "FILE_NAME", CPLGetFilename(osRasterFile));
    1122                 :     }
    1123              25 :     nWritingBytes += ISIS2Dataset::WriteFormatting( fpLabel, "");
    1124                 : 
    1125              50 :     nWritingBytes += ISIS2Dataset::WriteFormatting( fpLabel, "/* Pointers to Data Objects */");
    1126                 : 
    1127              25 :     if(bAttachedLabel){
    1128              24 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, CPLString().Printf("^%s",sObjectTag.c_str()), CPLString().Printf(CPL_FRMT_GUIB,iLabelRecords+1));
    1129                 :     }else{
    1130               1 :         nWritingBytes += ISIS2Dataset::WriteKeyword( fpLabel, iLevel, CPLString().Printf("^%s",sObjectTag.c_str()), CPLString().Printf("(\"%s\",1)",CPLGetFilename(osRasterFile)));
    1131                 :     }
    1132                 : 
    1133              25 :     if(EQUAL(sObjectTag, "QUBE")){
    1134              25 :         ISIS2Dataset::WriteQUBE_Information(fpLabel, iLevel, nWritingBytes, nXSize, nYSize, nBands, eType, pszInterleaving);
    1135                 :     }
    1136                 : 
    1137              25 :     nWritingBytes += ISIS2Dataset::WriteFormatting( fpLabel, "END");
    1138                 : 
    1139                 :     // check if file record is correct
    1140              25 :     unsigned int q = nWritingBytes/RECORD_SIZE;
    1141              25 :     if( q <= iLabelRecords){
    1142                 :         // correct we add space after the label end for complete from iLabelRecords
    1143              25 :         unsigned int nSpaceBytesToWrite = (unsigned int) (iLabelRecords * RECORD_SIZE - nWritingBytes);
    1144              25 :         VSIFPrintfL(fpLabel,"%*c", nSpaceBytesToWrite, ' ');
    1145                 :     }else{
    1146               0 :         iLabelRecords = q+1;
    1147               0 :         ISIS2Dataset::WriteLabel(osFilename, osRasterFile, sObjectTag, nXSize, nYSize, nBands, eType, iRecords, pszInterleaving, iLabelRecords);
    1148                 :     }
    1149              25 :     VSIFCloseL( fpLabel );
    1150                 :   
    1151              25 :     return TRUE;
    1152                 : }
    1153                 : 
    1154                 : /************************************************************************/
    1155                 : /*                            WriteKeyword()                            */
    1156                 : /************************************************************************/
    1157                 : 
    1158             476 :     unsigned int ISIS2Dataset::WriteKeyword(
    1159                 :         VSILFILE *fpLabel, unsigned int iLevel, CPLString key, CPLString value)
    1160                 : 
    1161                 :     {
    1162             476 :         CPLString tab = "";
    1163             476 :         iLevel *= 4; // each struct is idented by 4 spaces
    1164             476 :         int ret = VSIFPrintfL(fpLabel,"%*s%s=%s\n", iLevel, tab.c_str(), key.c_str(), value.c_str());
    1165             476 :         return ret;
    1166                 :     }
    1167                 : 
    1168                 : /************************************************************************/
    1169                 : /*                          WriteFormatting()                           */
    1170                 : /************************************************************************/
    1171                 : 
    1172             225 :     unsigned int ISIS2Dataset::WriteFormatting(VSILFILE *fpLabel, CPLString data)
    1173                 : 
    1174                 :     {
    1175             225 :         int ret = VSIFPrintfL(fpLabel,"%s\n", data.c_str());
    1176             225 :         return ret;
    1177                 :     }
    1178                 : 
    1179                 : /************************************************************************/
    1180                 : /*                         GDALRegister_ISIS2()                         */
    1181                 : /************************************************************************/
    1182                 : 
    1183             558 :     void GDALRegister_ISIS2()
    1184                 : 
    1185                 :     {
    1186                 :         GDALDriver  *poDriver;
    1187                 : 
    1188             558 :         if( GDALGetDriverByName( "ISIS2" ) == NULL )
    1189                 :         {
    1190             537 :             poDriver = new GDALDriver();
    1191                 :         
    1192             537 :             poDriver->SetDescription( "ISIS2" );
    1193                 :             poDriver->SetMetadataItem( GDAL_DMD_LONGNAME, 
    1194             537 :                                        "USGS Astrogeology ISIS cube (Version 2)" );
    1195             537 :             poDriver->SetMetadataItem( GDAL_DMD_HELPTOPIC, "frmt_isis2.html" );
    1196             537 :             poDriver->SetMetadataItem( GDAL_DCAP_VIRTUALIO, "YES" );
    1197             537 :             poDriver->SetMetadataItem( GDAL_DMD_CREATIONDATATYPES, "Byte Int16 UInt16 Float32 Float64");
    1198                 : 
    1199                 :             poDriver->SetMetadataItem( GDAL_DMD_CREATIONOPTIONLIST, 
    1200                 : "<CreationOptionList>\n"
    1201                 : "   <Option name='LABELING_METHOD' type='string-select' default='ATTACHED'>\n"
    1202                 : "     <Value>ATTACHED</Value>"
    1203                 : "     <Value>DETACHED</Value>"
    1204                 : "   </Option>"
    1205                 : "   <Option name='IMAGE_EXTENSION' type='string' default='cub'/>\n"
    1206             537 : "</CreationOptionList>\n" );
    1207                 : 
    1208             537 :             poDriver->pfnIdentify = ISIS2Dataset::Identify;
    1209             537 :             poDriver->pfnOpen = ISIS2Dataset::Open;
    1210             537 :             poDriver->pfnCreate = ISIS2Dataset::Create;
    1211                 : 
    1212             537 :             GetGDALDriverManager()->RegisterDriver( poDriver );
    1213                 :         }
    1214             558 :     }

Generated by: LCOV version 1.7