/****************************************************************************** * $Id: mffdataset.cpp,v 1.29 2006/01/05 16:20:34 gwalter Exp $ * * Project: GView * Purpose: Implementation of Atlantis MFF Support * Author: Frank Warmerdam, warmerda@home.com * ****************************************************************************** * Copyright (c) 2000, Frank Warmerdam * * Permission is hereby granted, free of charge, to any person obtaining a * copy of this software and associated documentation files (the "Software"), * to deal in the Software without restriction, including without limitation * the rights to use, copy, modify, merge, publish, distribute, sublicense, * and/or sell copies of the Software, and to permit persons to whom the * Software is furnished to do so, subject to the following conditions: * * The above copyright notice and this permission notice shall be included * in all copies or substantial portions of the Software. * * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER * DEALINGS IN THE SOFTWARE. ****************************************************************************** * * $Log: mffdataset.cpp,v $ * Revision 1.29 2006/01/05 16:20:34 gwalter * Atlantis->Vexcel in long name. * * Revision 1.28 2005/05/05 13:55:42 fwarmerdam * PAM Enable * * Revision 1.27 2003/04/02 19:05:59 dron * Large file support. * * Revision 1.26 2003/03/18 19:03:44 gwalter * Add flushing, fix coordinate interpretation. * * Revision 1.25 2003/03/13 18:34:02 gwalter * Fix data type determination in CreateCopy. * * Revision 1.24 2003/03/12 16:02:28 gwalter * Fix spheroid support. * * Revision 1.22 2003/03/03 20:10:06 gwalter * Updated MFF and HKV (MFF2) georeferencing support. * * Revision 1.21 2002/11/23 18:54:17 warmerda * added CREATIONDATATYPES metadata for drivers * * Revision 1.20 2002/09/04 06:50:37 warmerda * avoid static driver pointers * * Revision 1.19 2002/06/12 21:12:25 warmerda * update to metadata based driver info * * Revision 1.18 2001/11/11 23:51:00 warmerda * added required class keyword to friend declarations * * Revision 1.17 2001/07/18 04:51:57 warmerda * added CPL_CVSID * * Revision 1.16 2001/07/11 18:08:14 warmerda * Fixed pixel interpretation again. TOP_LEFT is really 0.0,0.0 but * BOTTOM_RIGHT is the top left corner of the bottom right pixel! * * Revision 1.15 2001/07/11 18:00:57 warmerda * Fixed off-by-half-pixel issues with MFF GCPs. * * Revision 1.14 2000/10/11 16:06:51 warmerda * added read support for GCPn and NUM_GCPs keywords * * Revision 1.13 2000/09/25 21:20:13 warmerda * avoid initialization warnings * * Revision 1.12 2000/08/24 14:26:09 warmerda * improved handling of read-only file sets * * Revision 1.11 2000/08/15 19:28:26 warmerda * added help topic * * Revision 1.10 2000/07/26 21:05:21 warmerda * added support for reading tiled MFF files * * Revision 1.9 2000/07/12 19:21:34 warmerda * fixed cleanup of GCPs * * Revision 1.8 2000/06/05 17:24:06 warmerda * added real complex support * * Revision 1.7 2000/05/18 22:06:03 warmerda * added gcp and metadata support * * Revision 1.6 2000/05/15 14:18:27 warmerda * added COMPLEX_INTERPRETATION metadata * * Revision 1.5 2000/04/21 21:59:26 warmerda * added overview support * * Revision 1.4 2000/03/24 20:00:46 warmerda * Don't require IMAGE_FILE_FORMAT. * * Revision 1.3 2000/03/13 14:34:42 warmerda * avoid const problem on write * * Revision 1.2 2000/03/06 21:51:32 warmerda * fixed docs * * Revision 1.1 2000/03/06 19:23:20 warmerda * New * */ #include "rawdataset.h" #include "cpl_string.h" #include #include "ogr_spatialref.h" #include "atlsci_spheroid.h" CPL_CVSID("$Id: mffdataset.cpp,v 1.29 2006/01/05 16:20:34 gwalter Exp $"); CPL_C_START void GDALRegister_MFF(void); CPL_C_END typedef enum { MFFPRJ_NONE, MFFPRJ_LL, MFFPRJ_UTM, MFFPRJ_UNRECOGNIZED } ; static int GetMFFProjectionType(const char * pszNewProjection); /************************************************************************/ /* ==================================================================== */ /* MFFDataset */ /* ==================================================================== */ /************************************************************************/ class MFFDataset : public RawDataset { FILE *fpImage; // image data file. int nGCPCount; GDAL_GCP *pasGCPList; char *pszProjection; char *pszGCPProjection; double adfGeoTransform[6]; void ScanForGCPs(); void ScanForProjectionInfo(); public: MFFDataset(); ~MFFDataset(); char **papszHdrLines; FILE **pafpBandFiles; virtual int GetGCPCount(); virtual const char *GetGCPProjection(); virtual const GDAL_GCP *GetGCPs(); virtual const char *GetProjectionRef(); virtual CPLErr GetGeoTransform( double * ); static GDALDataset *Open( GDALOpenInfo * ); static GDALDataset *Create( const char * pszFilename, int nXSize, int nYSize, int nBands, GDALDataType eType, char ** papszParmList ); static GDALDataset *CreateCopy( const char * pszFilename, GDALDataset *poSrcDS, int bStrict, char ** papszOptions, GDALProgressFunc pfnProgress, void * pProgressData ); }; /************************************************************************/ /* ==================================================================== */ /* MFFTiledBand */ /* ==================================================================== */ /************************************************************************/ class MFFTiledBand : public GDALRasterBand { friend class MFFDataset; FILE *fpRaw; int bNative; public: MFFTiledBand( MFFDataset *, int, FILE *, int, int, GDALDataType, int ); ~MFFTiledBand(); virtual CPLErr IReadBlock( int, int, void * ); }; /************************************************************************/ /* MFFTiledBand() */ /************************************************************************/ MFFTiledBand::MFFTiledBand( MFFDataset *poDS, int nBand, FILE *fp, int nTileXSize, int nTileYSize, GDALDataType eDataType, int bNative ) { this->poDS = poDS; this->nBand = nBand; this->eDataType = eDataType; this->bNative = bNative; this->nBlockXSize = nTileXSize; this->nBlockYSize = nTileYSize; this->fpRaw = fp; } /************************************************************************/ /* ~MFFTiledBand() */ /************************************************************************/ MFFTiledBand::~MFFTiledBand() { VSIFCloseL( fpRaw ); } /************************************************************************/ /* IReadBlock() */ /************************************************************************/ CPLErr MFFTiledBand::IReadBlock( int nBlockXOff, int nBlockYOff, void * pImage ) { long nOffset; int nTilesPerRow; int nWordSize, nBlockSize; nTilesPerRow = (nRasterXSize + nBlockXSize - 1) / nBlockXSize; nWordSize = GDALGetDataTypeSize( eDataType ) / 8; nBlockSize = nWordSize * nBlockXSize * nBlockYSize; nOffset = nBlockSize * (nBlockXOff + nBlockYOff*nTilesPerRow); if( VSIFSeekL( fpRaw, nOffset, SEEK_SET ) == -1 || VSIFReadL( pImage, 1, nBlockSize, fpRaw ) < 1 ) { CPLError( CE_Failure, CPLE_FileIO, "Read of tile %d/%d failed with fseek or fread error.", nBlockXOff, nBlockYOff ); return CE_Failure; } if( !bNative && nWordSize > 1 ) { if( GDALDataTypeIsComplex( eDataType ) ) { GDALSwapWords( pImage, nWordSize/2, nBlockXSize*nBlockYSize, nWordSize ); GDALSwapWords( ((GByte *) pImage)+nWordSize/2, nWordSize/2, nBlockXSize*nBlockYSize, nWordSize ); } else GDALSwapWords( pImage, nWordSize, nBlockXSize * nBlockYSize, nWordSize ); } return CE_None; } /************************************************************************/ /* MFF Spheroids */ /************************************************************************/ class MFFSpheroidList : public SpheroidList { public: MFFSpheroidList(); ~MFFSpheroidList(); }; MFFSpheroidList :: MFFSpheroidList() { num_spheroids = 18; epsilonR = 0.1; epsilonI = 0.000001; spheroids[0].SetValuesByRadii("SPHERE",6371007.0,6371007.0); spheroids[1].SetValuesByRadii("EVEREST",6377304.0,6356103.0); spheroids[2].SetValuesByRadii("BESSEL",6377397.0,6356082.0); spheroids[3].SetValuesByRadii("AIRY",6377563.0,6356300.0); spheroids[4].SetValuesByRadii("CLARKE_1858",6378294.0,6356621.0); spheroids[5].SetValuesByRadii("CLARKE_1866",6378206.4,6356583.8); spheroids[6].SetValuesByRadii("CLARKE_1880",6378249.0,6356517.0); spheroids[7].SetValuesByRadii("HAYFORD",6378388.0,6356915.0); spheroids[8].SetValuesByRadii("KRASOVSKI",6378245.0,6356863.0); spheroids[9].SetValuesByRadii("HOUGH",6378270.0,6356794.0); spheroids[10].SetValuesByRadii("FISHER_60",6378166.0,6356784.0); spheroids[11].SetValuesByRadii("KAULA",6378165.0,6356345.0); spheroids[12].SetValuesByRadii("IUGG_67",6378160.0,6356775.0); spheroids[13].SetValuesByRadii("FISHER_68",6378150.0,6356330.0); spheroids[14].SetValuesByRadii("WGS_72",6378135.0,6356751.0); spheroids[15].SetValuesByRadii("IUGG_75",6378140.0,6356755.0); spheroids[16].SetValuesByRadii("WGS_84",6378137.0,6356752.0); spheroids[17].SetValuesByRadii("HUGHES",6378273.0,6356889.4); } MFFSpheroidList::~MFFSpheroidList() { } /************************************************************************/ /* ==================================================================== */ /* MFFDataset */ /* ==================================================================== */ /************************************************************************/ /************************************************************************/ /* MFFDataset() */ /************************************************************************/ MFFDataset::MFFDataset() { papszHdrLines = NULL; pafpBandFiles = NULL; nGCPCount = 0; pasGCPList = NULL; pszProjection = CPLStrdup(""); pszGCPProjection = CPLStrdup(""); adfGeoTransform[0] = 0.0; adfGeoTransform[1] = 1.0; adfGeoTransform[2] = 0.0; adfGeoTransform[3] = 0.0; adfGeoTransform[4] = 0.0; adfGeoTransform[5] = 1.0; } /************************************************************************/ /* ~MFFDataset() */ /************************************************************************/ MFFDataset::~MFFDataset() { FlushCache(); CSLDestroy( papszHdrLines ); if( pafpBandFiles != NULL ) { for( int i = 0; i < GetRasterCount(); i++ ) { if( pafpBandFiles[i] != NULL ) VSIFCloseL( pafpBandFiles[i] ); } CPLFree( pafpBandFiles ); } if( nGCPCount > 0 ) { GDALDeinitGCPs( nGCPCount, pasGCPList ); CPLFree( pasGCPList ); } CPLFree( pszProjection ); CPLFree( pszGCPProjection ); } /************************************************************************/ /* GetGCPCount() */ /************************************************************************/ int MFFDataset::GetGCPCount() { return nGCPCount; } /************************************************************************/ /* GetGCPProjection() */ /************************************************************************/ const char *MFFDataset::GetGCPProjection() { if( nGCPCount > 0 ) return pszGCPProjection; else return ""; } /************************************************************************/ /* GetProjectionRef() */ /************************************************************************/ const char *MFFDataset::GetProjectionRef() { return ( pszProjection ); } /************************************************************************/ /* GetGeoTransform() */ /************************************************************************/ CPLErr MFFDataset::GetGeoTransform( double * padfTransform ) { memcpy( padfTransform, adfGeoTransform, sizeof(double) * 6 ); return( CE_None ); } /************************************************************************/ /* GetGCP() */ /************************************************************************/ const GDAL_GCP *MFFDataset::GetGCPs() { return pasGCPList; } /************************************************************************/ /* ScanForGCPs() */ /************************************************************************/ void MFFDataset::ScanForGCPs() { int nCorner; int NUM_GCPS = 0; if( CSLFetchNameValue(papszHdrLines, "NUM_GCPS") != NULL ) NUM_GCPS = atoi(CSLFetchNameValue(papszHdrLines, "NUM_GCPS")); nGCPCount = 0; pasGCPList = (GDAL_GCP *) CPLCalloc(sizeof(GDAL_GCP),5+NUM_GCPS); for( nCorner = 0; nCorner < 5; nCorner++ ) { const char * pszBase=NULL; double dfRasterX=0.0, dfRasterY=0.0; char szLatName[40], szLongName[40]; if( nCorner == 0 ) { dfRasterX = 0.5; dfRasterY = 0.5; pszBase = "TOP_LEFT_CORNER"; } else if( nCorner == 1 ) { dfRasterX = GetRasterXSize()-0.5; dfRasterY = 0.5; pszBase = "TOP_RIGHT_CORNER"; } else if( nCorner == 2 ) { dfRasterX = GetRasterXSize()-0.5; dfRasterY = GetRasterYSize()-0.5; pszBase = "BOTTOM_RIGHT_CORNER"; } else if( nCorner == 3 ) { dfRasterX = 0.5; dfRasterY = GetRasterYSize()-0.5; pszBase = "BOTTOM_LEFT_CORNER"; } else if( nCorner == 4 ) { dfRasterX = GetRasterXSize()/2.0; dfRasterY = GetRasterYSize()/2.0; pszBase = "CENTRE"; } sprintf( szLatName, "%s_LATITUDE", pszBase ); sprintf( szLongName, "%s_LONGITUDE", pszBase ); if( CSLFetchNameValue(papszHdrLines, szLatName) != NULL && CSLFetchNameValue(papszHdrLines, szLongName) != NULL ) { GDALInitGCPs( 1, pasGCPList + nGCPCount ); CPLFree( pasGCPList[nGCPCount].pszId ); pasGCPList[nGCPCount].pszId = CPLStrdup( pszBase ); pasGCPList[nGCPCount].dfGCPX = atof(CSLFetchNameValue(papszHdrLines, szLongName)); pasGCPList[nGCPCount].dfGCPY = atof(CSLFetchNameValue(papszHdrLines, szLatName)); pasGCPList[nGCPCount].dfGCPZ = 0.0; pasGCPList[nGCPCount].dfGCPPixel = dfRasterX; pasGCPList[nGCPCount].dfGCPLine = dfRasterY; nGCPCount++; } } /* -------------------------------------------------------------------- */ /* Collect standalone GCPs. They look like: */ /* */ /* GCPn = row, col, lat, long */ /* GCP1 = 1, 1, 45.0, -75.0 */ /* -------------------------------------------------------------------- */ int i; for( i = 0; i < NUM_GCPS; i++ ) { char szName[25]; char **papszTokens; sprintf( szName, "GCP%d", i+1 ); if( CSLFetchNameValue( papszHdrLines, szName ) == NULL ) continue; papszTokens = CSLTokenizeStringComplex( CSLFetchNameValue( papszHdrLines, szName ), ",", FALSE, FALSE ); if( CSLCount(papszTokens) == 4 ) { GDALInitGCPs( 1, pasGCPList + nGCPCount ); CPLFree( pasGCPList[nGCPCount].pszId ); pasGCPList[nGCPCount].pszId = CPLStrdup( szName ); pasGCPList[nGCPCount].dfGCPX = atof(papszTokens[3]); pasGCPList[nGCPCount].dfGCPY = atof(papszTokens[2]); pasGCPList[nGCPCount].dfGCPZ = 0.0; pasGCPList[nGCPCount].dfGCPPixel = atof(papszTokens[1])+0.5; pasGCPList[nGCPCount].dfGCPLine = atof(papszTokens[0])+0.5; nGCPCount++; } } } /************************************************************************/ /* ScanForProjectionInfo */ /************************************************************************/ void MFFDataset::ScanForProjectionInfo() { const char *pszProjName, *pszOriginLong, *pszSpheroidName; const char *pszSpheroidEqRadius, *pszSpheroidPolarRadius; double eq_radius, polar_radius; OGRSpatialReference oProj; OGRSpatialReference oLL; MFFSpheroidList *mffEllipsoids; pszProjName = CSLFetchNameValue(papszHdrLines, "PROJECTION_NAME"); pszOriginLong = CSLFetchNameValue(papszHdrLines, "PROJECTION_ORIGIN_LONGITUDE"); pszSpheroidName = CSLFetchNameValue(papszHdrLines, "SPHEROID_NAME"); if (pszProjName == NULL) { CPLFree( pszProjection ); CPLFree( pszGCPProjection ); pszProjection=CPLStrdup(""); pszGCPProjection=CPLStrdup(""); return; } else if ((!EQUAL(pszProjName,"utm")) && (!EQUAL(pszProjName,"ll"))) { CPLError(CE_Warning,CPLE_AppDefined, "Warning- only utm and lat/long projections are currently supported."); CPLFree( pszProjection ); CPLFree( pszGCPProjection ); pszProjection=CPLStrdup(""); pszGCPProjection=CPLStrdup(""); return; } mffEllipsoids = new MFFSpheroidList; if( EQUAL(pszProjName,"utm") ) { int nZone; if (pszOriginLong == NULL) { /* If origin not specified, assume 0.0 */ CPLError(CE_Warning,CPLE_AppDefined, "Warning- no projection origin longitude specified. Assuming 0.0."); nZone = 31; } else nZone = 31 + (int) floor(atof(pszOriginLong)/6.0); if( pasGCPList[4].dfGCPY < 0 ) oProj.SetUTM( nZone, 0 ); else oProj.SetUTM( nZone, 1 ); if (pszOriginLong != NULL) oProj.SetProjParm(SRS_PP_CENTRAL_MERIDIAN,atof(pszOriginLong)); } if (pszOriginLong != NULL) oLL.SetProjParm(SRS_PP_LONGITUDE_OF_ORIGIN,atof(pszOriginLong)); if ((pszSpheroidName == NULL)) { CPLError(CE_Warning,CPLE_AppDefined, "Warning- unspecified ellipsoid. Using wgs-84 parameters.\n"); oProj.SetWellKnownGeogCS( "WGS84" ); oLL.SetWellKnownGeogCS( "WGS84" ); } else { if (mffEllipsoids->SpheroidInList(pszSpheroidName)) { oProj.SetGeogCS( "unknown","unknown",pszSpheroidName, mffEllipsoids->GetSpheroidEqRadius(pszSpheroidName), mffEllipsoids->GetSpheroidInverseFlattening(pszSpheroidName) ); oLL.SetGeogCS( "unknown","unknown",pszSpheroidName, mffEllipsoids->GetSpheroidEqRadius(pszSpheroidName), mffEllipsoids->GetSpheroidInverseFlattening(pszSpheroidName) ); } else if (EQUAL(pszSpheroidName,"USER_DEFINED")) { pszSpheroidEqRadius = CSLFetchNameValue(papszHdrLines, "SPHEROID_EQUATORIAL_RADIUS"); pszSpheroidPolarRadius = CSLFetchNameValue(papszHdrLines, "SPHEROID_POLAR_RADIUS"); if ((pszSpheroidEqRadius != NULL) && (pszSpheroidPolarRadius != NULL)) { eq_radius = atof( pszSpheroidEqRadius ); polar_radius = atof( pszSpheroidPolarRadius ); oProj.SetGeogCS( "unknown","unknown","unknown", eq_radius, eq_radius/(eq_radius - polar_radius)); oLL.SetGeogCS( "unknown","unknown","unknown", eq_radius, eq_radius/(eq_radius - polar_radius)); } else { CPLError(CE_Warning,CPLE_AppDefined, "Warning- radii not specified for user-defined ellipsoid. Using wgs-84 parameters. \n"); oProj.SetWellKnownGeogCS( "WGS84" ); oLL.SetWellKnownGeogCS( "WGS84" ); } } else { CPLError(CE_Warning,CPLE_AppDefined, "Warning- unrecognized ellipsoid. Using wgs-84 parameters.\n"); oProj.SetWellKnownGeogCS( "WGS84" ); oLL.SetWellKnownGeogCS( "WGS84" ); } } /* If a geotransform is sufficient to represent the GCP's (ie. each */ /* estimated gcp is within 0.25*pixel size of the actual value- this */ /* is the test applied by GDALGCPsToGeoTransform), store the */ /* geotransform. */ int transform_ok = FALSE; if (EQUAL(pszProjName,"LL")) { transform_ok = GDALGCPsToGeoTransform(nGCPCount,pasGCPList,adfGeoTransform,0); } else { OGRCoordinateTransformation *poTransform = NULL; double *dfPrjX, *dfPrjY; int gcp_index; int bSuccess = TRUE; dfPrjX = (double *) CPLMalloc(nGCPCount*sizeof(double)); dfPrjY = (double *) CPLMalloc(nGCPCount*sizeof(double)); poTransform = OGRCreateCoordinateTransformation( &oLL, &oProj ); if( poTransform == NULL ) bSuccess = FALSE; for(gcp_index=0;gcp_indexTransform( 1, &(dfPrjX[gcp_index]), &(dfPrjY[gcp_index]) ) ) bSuccess = FALSE; } if( bSuccess ) { for(gcp_index=0;gcp_indexnHeaderBytes < 17 || poOpenInfo->fp == NULL ) return NULL; if( !EQUAL(CPLGetExtension(poOpenInfo->pszFilename),"hdr") ) return NULL; /* -------------------------------------------------------------------- */ /* Load the .hdr file, and compress white space out around the */ /* equal sign. */ /* -------------------------------------------------------------------- */ papszHdrLines = CSLLoad( poOpenInfo->pszFilename ); if( papszHdrLines == NULL ) return NULL; for( i = 0; papszHdrLines[i] != NULL; i++ ) { int bAfterEqual = FALSE; int iSrc, iDst; char *pszLine = papszHdrLines[i]; for( iSrc=0, iDst=0; pszLine[iSrc] != '\0'; iSrc++ ) { if( bAfterEqual || pszLine[iSrc] != ' ' ) { pszLine[iDst++] = pszLine[iSrc]; } if( iDst > 0 && pszLine[iDst-1] == '=' ) bAfterEqual = FALSE; } pszLine[iDst] = '\0'; } /* -------------------------------------------------------------------- */ /* Verify it is an MFF file. */ /* -------------------------------------------------------------------- */ if( CSLFetchNameValue( papszHdrLines, "IMAGE_FILE_FORMAT" ) != NULL && !EQUAL(CSLFetchNameValue(papszHdrLines,"IMAGE_FILE_FORMAT"),"MFF") ) { CSLDestroy( papszHdrLines ); return NULL; } if( (CSLFetchNameValue( papszHdrLines, "IMAGE_LINES" ) == NULL || CSLFetchNameValue(papszHdrLines,"LINE_SAMPLES") == NULL) && (CSLFetchNameValue( papszHdrLines, "no_rows" ) == NULL || CSLFetchNameValue(papszHdrLines,"no_columns") == NULL) ) { CSLDestroy( papszHdrLines ); return NULL; } /* -------------------------------------------------------------------- */ /* Create a corresponding GDALDataset. */ /* -------------------------------------------------------------------- */ MFFDataset *poDS; poDS = new MFFDataset(); poDS->papszHdrLines = papszHdrLines; /* -------------------------------------------------------------------- */ /* Set some dataset wide information. */ /* -------------------------------------------------------------------- */ if( CSLFetchNameValue(papszHdrLines,"no_rows") != NULL && CSLFetchNameValue(papszHdrLines,"no_columns") != NULL ) { poDS->RasterInitialize( atoi(CSLFetchNameValue(papszHdrLines,"no_columns")), atoi(CSLFetchNameValue(papszHdrLines,"no_rows")) ); } else { poDS->RasterInitialize( atoi(CSLFetchNameValue(papszHdrLines,"LINE_SAMPLES")), atoi(CSLFetchNameValue(papszHdrLines,"IMAGE_LINES")) ); } if( CSLFetchNameValue( papszHdrLines, "BYTE_ORDER" ) != NULL ) { #ifdef CPL_MSB bNative = EQUAL(CSLFetchNameValue(papszHdrLines,"BYTE_ORDER"),"MSB"); #else bNative = EQUAL(CSLFetchNameValue(papszHdrLines,"BYTE_ORDER"),"LSB"); #endif } /* -------------------------------------------------------------------- */ /* Get some information specific to APP tiled files. */ /* -------------------------------------------------------------------- */ int bTiled, nTileXSize=0, nTileYSize=0; const char *pszRefinedType = NULL; pszRefinedType = CSLFetchNameValue(papszHdrLines, "type" ); bTiled = CSLFetchNameValue(papszHdrLines,"no_rows") != NULL; if( bTiled ) { if( CSLFetchNameValue(papszHdrLines,"tile_size_rows") ) nTileYSize = atoi(CSLFetchNameValue(papszHdrLines,"tile_size_rows")); if( CSLFetchNameValue(papszHdrLines,"tile_size_columns") ) nTileXSize = atoi(CSLFetchNameValue(papszHdrLines,"tile_size_columns")); } /* -------------------------------------------------------------------- */ /* Read the directory to find matching band files. */ /* -------------------------------------------------------------------- */ char **papszDirFiles; char *pszTargetBase, *pszTargetPath; int nRawBand, nSkipped=0; pszTargetPath = CPLStrdup(CPLGetPath(poOpenInfo->pszFilename)); pszTargetBase = CPLStrdup(CPLGetBasename( poOpenInfo->pszFilename )); papszDirFiles = CPLReadDir( CPLGetPath( poOpenInfo->pszFilename ) ); if( papszDirFiles == NULL ) return NULL; for( nRawBand = 0; TRUE; nRawBand++ ) { const char *pszExtension; int nBand; GDALDataType eDataType; /* Find the next raw band file. */ for( i = 0; papszDirFiles[i] != NULL; i++ ) { if( !EQUAL(CPLGetBasename(papszDirFiles[i]),pszTargetBase) ) continue; pszExtension = CPLGetExtension(papszDirFiles[i]); if( isdigit(pszExtension[1]) && atoi(pszExtension+1) == nRawBand && strchr("bBcCiIjJrRxXzZ",pszExtension[0]) != NULL ) break; } if( papszDirFiles[i] == NULL ) break; /* open the file for required level of access */ FILE *fpRaw; const char *pszRawFilename = CPLFormFilename(pszTargetPath, papszDirFiles[i], NULL ); if( poOpenInfo->eAccess == GA_Update ) fpRaw = VSIFOpenL( pszRawFilename, "rb+" ); else fpRaw = VSIFOpenL( pszRawFilename, "rb" ); if( fpRaw == NULL ) { CPLError( CE_Warning, CPLE_OpenFailed, "Unable to open %s ... skipping.\n", pszRawFilename ); nSkipped++; continue; } pszExtension = CPLGetExtension(papszDirFiles[i]); if( pszRefinedType != NULL ) { if( EQUAL(pszRefinedType,"C*4") ) eDataType = GDT_CFloat32; else if( EQUAL(pszRefinedType,"C*8") ) eDataType = GDT_CFloat64; else if( EQUAL(pszRefinedType,"R*4") ) eDataType = GDT_Float32; else if( EQUAL(pszRefinedType,"R*8") ) eDataType = GDT_Float64; else if( EQUAL(pszRefinedType,"I*1") ) eDataType = GDT_Byte; else if( EQUAL(pszRefinedType,"I*2") ) eDataType = GDT_Int16; else if( EQUAL(pszRefinedType,"I*4") ) eDataType = GDT_Int32; else if( EQUAL(pszRefinedType,"U*2") ) eDataType = GDT_UInt16; else if( EQUAL(pszRefinedType,"U*4") ) eDataType = GDT_UInt32; else if( EQUAL(pszRefinedType,"J*1") ) continue; /* we don't support 1 byte complex */ else if( EQUAL(pszRefinedType,"J*2") ) eDataType = GDT_CInt16; else if( EQUAL(pszRefinedType,"K*4") ) eDataType = GDT_CInt32; else continue; } else if( EQUALN(pszExtension,"b",1) ) { eDataType = GDT_Byte; } else if( EQUALN(pszExtension,"i",1) ) { eDataType = GDT_UInt16; } else if( EQUALN(pszExtension,"j",1) ) { eDataType = GDT_CInt16; } else if( EQUALN(pszExtension,"r",1) ) { eDataType = GDT_Float32; } else if( EQUALN(pszExtension,"x",1) ) { eDataType = GDT_CFloat32; } else continue; nBand = poDS->GetRasterCount() + 1; int nPixelOffset = GDALGetDataTypeSize(eDataType)/8; GDALRasterBand *poBand = NULL; if( bTiled ) { poBand = new MFFTiledBand( poDS, nBand, fpRaw, nTileXSize, nTileYSize, eDataType, bNative ); } else { poBand = new RawRasterBand( poDS, nBand, fpRaw, 0, nPixelOffset, nPixelOffset * poDS->GetRasterXSize(), eDataType, bNative, TRUE ); } poDS->SetBand( nBand, poBand ); } CPLFree(pszTargetPath); CPLFree(pszTargetBase); CSLDestroy(papszDirFiles); /* -------------------------------------------------------------------- */ /* Check if we have bands. */ /* -------------------------------------------------------------------- */ if( poDS->GetRasterCount() == 0 ) { if( nSkipped > 0 && poOpenInfo->eAccess ) { CPLError( CE_Failure, CPLE_OpenFailed, "Failed to open %d files that were apparently bands.\n" "Perhaps this dataset is readonly?\n", nSkipped ); delete poDS; return NULL; } else { CPLError( CE_Failure, CPLE_OpenFailed, "MFF header file read successfully, but no bands\n" "were successfully found and opened." ); delete poDS; return NULL; } } /* -------------------------------------------------------------------- */ /* Check for overviews. */ /* -------------------------------------------------------------------- */ poDS->oOvManager.Initialize( poDS, poOpenInfo->pszFilename ); /* -------------------------------------------------------------------- */ /* Initialize any PAM information. */ /* -------------------------------------------------------------------- */ poDS->SetDescription( poOpenInfo->pszFilename ); poDS->TryLoadXML(); /* -------------------------------------------------------------------- */ /* Set all information from the .hdr that isn't well know to be */ /* metadata. */ /* -------------------------------------------------------------------- */ for( i = 0; papszHdrLines[i] != NULL; i++ ) { const char *pszValue; char *pszName; pszValue = CPLParseNameValue(papszHdrLines[i], &pszName); if( pszName == NULL || pszValue == NULL ) continue; if( !EQUAL(pszName,"END") && !EQUAL(pszName,"FILE_TYPE") && !EQUAL(pszName,"BYTE_ORDER") && !EQUAL(pszName,"no_columns") && !EQUAL(pszName,"no_rows") && !EQUAL(pszName,"type") && !EQUAL(pszName,"tile_size_rows") && !EQUAL(pszName,"tile_size_columns") && !EQUAL(pszName,"IMAGE_FILE_FORMAT") && !EQUAL(pszName,"IMAGE_LINES") && !EQUAL(pszName,"LINE_SAMPLES") ) { poDS->SetMetadataItem( pszName, pszValue ); } CPLFree( pszName ); } /* -------------------------------------------------------------------- */ /* Any GCPs in header file? */ /* -------------------------------------------------------------------- */ poDS->ScanForGCPs(); poDS->ScanForProjectionInfo(); return( poDS ); } int GetMFFProjectionType(const char *pszNewProjection) { OGRSpatialReference *oSRS; char *modifiableProjection = NULL; if( !EQUALN(pszNewProjection,"GEOGCS",6) && !EQUALN(pszNewProjection,"PROJCS",6) && !EQUAL(pszNewProjection,"") ) { return MFFPRJ_UNRECOGNIZED; } else if (EQUAL(pszNewProjection,"")) { return MFFPRJ_NONE; } else { /* importFromWkt updates the pointer, so don't use pszNewProjection directly */ modifiableProjection=CPLStrdup(pszNewProjection); oSRS = new OGRSpatialReference; oSRS->importFromWkt(&modifiableProjection); if ((oSRS->GetAttrValue("PROJECTION") != NULL) && (EQUAL(oSRS->GetAttrValue("PROJECTION"),SRS_PT_TRANSVERSE_MERCATOR))) { return MFFPRJ_UTM; } else if ((oSRS->GetAttrValue("PROJECTION") == NULL) && (oSRS->IsGeographic())) { return MFFPRJ_LL; } else { return MFFPRJ_UNRECOGNIZED; } } } /************************************************************************/ /* Create() */ /************************************************************************/ GDALDataset *MFFDataset::Create( const char * pszFilenameIn, int nXSize, int nYSize, int nBands, GDALDataType eType, char ** papszParmList ) { /* -------------------------------------------------------------------- */ /* Verify input options. */ /* -------------------------------------------------------------------- */ if( eType != GDT_Byte && eType != GDT_Float32 && eType != GDT_UInt16 && eType != GDT_CInt16 && eType != GDT_CFloat32 ) { CPLError( CE_Failure, CPLE_AppDefined, "Attempt to create MFF file with currently unsupported\n" "data type (%s).\n", GDALGetDataTypeName(eType) ); return NULL; } /* -------------------------------------------------------------------- */ /* Establish the base filename (path+filename, less extension). */ /* -------------------------------------------------------------------- */ char *pszBaseFilename; int i; pszBaseFilename = (char *) CPLMalloc(strlen(pszFilenameIn)+5); strcpy( pszBaseFilename, pszFilenameIn ); for( i = strlen(pszBaseFilename)-1; i > 0; i-- ) { if( pszBaseFilename[i] == '.' ) { pszBaseFilename[i] = '\0'; break; } if( pszBaseFilename[i] == '/' || pszBaseFilename[i] == '\\' ) break; } /* -------------------------------------------------------------------- */ /* Create the header file. */ /* -------------------------------------------------------------------- */ FILE *fp; const char *pszFilename; pszFilename = CPLFormFilename( NULL, pszBaseFilename, "hdr" ); fp = VSIFOpen( pszFilename, "wt" ); if( fp == NULL ) { CPLError( CE_Failure, CPLE_OpenFailed, "Couldn't create %s.\n", pszFilename ); return NULL; } fprintf( fp, "IMAGE_FILE_FORMAT = MFF\n" ); fprintf( fp, "FILE_TYPE = IMAGE\n" ); fprintf( fp, "IMAGE_LINES = %d\n", nYSize ); fprintf( fp, "LINE_SAMPLES = %d\n", nXSize ); #ifdef CPL_MSB fprintf( fp, "BYTE_ORDER = MSB\n" ); #else fprintf( fp, "BYTE_ORDER = LSB\n" ); #endif if (CSLFetchNameValue(papszParmList,"NO_END") == NULL) fprintf( fp, "END\n" ); VSIFClose( fp ); /* -------------------------------------------------------------------- */ /* Create the data files, but don't bother writing any data to them.*/ /* -------------------------------------------------------------------- */ for( int iBand = 0; iBand < nBands; iBand++ ) { char szExtension[4]; if( eType == GDT_Byte ) sprintf( szExtension, "b%02d", iBand ); else if( eType == GDT_UInt16 ) sprintf( szExtension, "i%02d", iBand ); else if( eType == GDT_Float32 ) sprintf( szExtension, "r%02d", iBand ); else if( eType == GDT_CInt16 ) sprintf( szExtension, "j%02d", iBand ); else if( eType == GDT_CFloat32 ) sprintf( szExtension, "x%02d", iBand ); pszFilename = CPLFormFilename( NULL, pszBaseFilename, szExtension ); fp = VSIFOpen( pszFilename, "wb" ); if( fp == NULL ) { CPLError( CE_Failure, CPLE_OpenFailed, "Couldn't create %s.\n", pszFilename ); return NULL; } VSIFWrite( (void *) "", 1, 1, fp ); VSIFClose( fp ); } /* -------------------------------------------------------------------- */ /* Open the dataset normally. */ /* -------------------------------------------------------------------- */ GDALDataset *poDS; strcat( pszBaseFilename, ".hdr" ); poDS = (GDALDataset *) GDALOpen( pszBaseFilename, GA_Update ); CPLFree( pszBaseFilename ); return poDS; } /************************************************************************/ /* CreateCopy() */ /************************************************************************/ GDALDataset * MFFDataset::CreateCopy( const char * pszFilename, GDALDataset *poSrcDS, int bStrict, char ** papszOptions, GDALProgressFunc pfnProgress, void * pProgressData ) { MFFDataset *poDS; GDALDataType eType = poSrcDS->GetRasterBand(1)->GetRasterDataType(); int iBand; char **newpapszOptions=NULL; if( !pfnProgress( 0.0, NULL, pProgressData ) ) return NULL; /* check that other bands match type- sets type */ /* to unknown if they differ. */ for( iBand = 1; iBand < poSrcDS->GetRasterCount(); iBand++ ) { GDALRasterBand *poBand = poSrcDS->GetRasterBand( iBand+1 ); eType = GDALDataTypeUnion( eType, poBand->GetRasterDataType() ); } newpapszOptions=CSLDuplicate(papszOptions); newpapszOptions=CSLSetNameValue(newpapszOptions,"NO_END","TRUE"); poDS = (MFFDataset *) Create( pszFilename, poSrcDS->GetRasterXSize(), poSrcDS->GetRasterYSize(), poSrcDS->GetRasterCount(), eType, newpapszOptions ); CSLDestroy(newpapszOptions); /* Check that Create worked- return Null if it didn't */ if (poDS == NULL) return NULL; /* -------------------------------------------------------------------- */ /* Copy the image data. */ /* -------------------------------------------------------------------- */ int nXSize = poDS->GetRasterXSize(); int nYSize = poDS->GetRasterYSize(); int nBlockXSize, nBlockYSize, nBlockTotal, nBlocksDone; poDS->GetRasterBand(1)->GetBlockSize( &nBlockXSize, &nBlockYSize ); nBlockTotal = ((nXSize + nBlockXSize - 1) / nBlockXSize) * ((nYSize + nBlockYSize - 1) / nBlockYSize) * poSrcDS->GetRasterCount(); nBlocksDone = 0; for( iBand = 0; iBand < poSrcDS->GetRasterCount(); iBand++ ) { GDALRasterBand *poSrcBand = poSrcDS->GetRasterBand( iBand+1 ); GDALRasterBand *poDstBand = poDS->GetRasterBand( iBand+1 ); int iYOffset, iXOffset; void *pData; CPLErr eErr; pData = CPLMalloc(nBlockXSize * nBlockYSize * GDALGetDataTypeSize(eType) / 8); for( iYOffset = 0; iYOffset < nYSize; iYOffset += nBlockYSize ) { for( iXOffset = 0; iXOffset < nXSize; iXOffset += nBlockXSize ) { int nTBXSize, nTBYSize; if( !pfnProgress( (nBlocksDone++) / (float) nBlockTotal, NULL, pProgressData ) ) { CPLError( CE_Failure, CPLE_UserInterrupt, "User terminated" ); delete poDS; GDALDriver *poMFFDriver = (GDALDriver *) GDALGetDriverByName( "MFF" ); poMFFDriver->Delete( pszFilename ); return NULL; } nTBXSize = MIN(nBlockXSize,nXSize-iXOffset); nTBYSize = MIN(nBlockYSize,nYSize-iYOffset); eErr = poSrcBand->RasterIO( GF_Read, iXOffset, iYOffset, nTBXSize, nTBYSize, pData, nTBXSize, nTBYSize, eType, 0, 0 ); if( eErr != CE_None ) { return NULL; } eErr = poDstBand->RasterIO( GF_Write, iXOffset, iYOffset, nTBXSize, nTBYSize, pData, nTBXSize, nTBYSize, eType, 0, 0 ); if( eErr != CE_None ) { return NULL; } } } CPLFree( pData ); } /* -------------------------------------------------------------------- */ /* Copy georeferencing information, if enough is available. */ /* -------------------------------------------------------------------- */ /* -------------------------------------------------------------------- */ /* Establish the base filename (path+filename, less extension). */ /* -------------------------------------------------------------------- */ char *pszBaseFilename; int i; FILE *fp; const char *pszFilenameGEO; pszBaseFilename = (char *) CPLMalloc(strlen(pszFilename)+5); strcpy( pszBaseFilename, pszFilename ); for( i = strlen(pszBaseFilename)-1; i > 0; i-- ) { if( pszBaseFilename[i] == '.' ) { pszBaseFilename[i] = '\0'; break; } if( pszBaseFilename[i] == '/' || pszBaseFilename[i] == '\\' ) break; } pszFilenameGEO = CPLFormFilename( NULL, pszBaseFilename, "hdr" ); fp = VSIFOpen( pszFilenameGEO, "at" ); if( fp == NULL ) { CPLError( CE_Failure, CPLE_OpenFailed, "Couldn't open %s for appending.\n", pszFilenameGEO ); return NULL; } /* MFF requires corner and center gcps */ double *padfTiepoints; int src_prj; int georef_created = FALSE; padfTiepoints = (double *) CPLMalloc(2*sizeof(double)*5); src_prj = GetMFFProjectionType(poSrcDS->GetProjectionRef()); if ((src_prj != MFFPRJ_NONE) && (src_prj != MFFPRJ_UNRECOGNIZED)) { double *tempGeoTransform = NULL; tempGeoTransform = (double *) CPLMalloc(6*sizeof(double)); if (( poSrcDS->GetGeoTransform( tempGeoTransform ) == CE_None) && (tempGeoTransform[0] != 0.0 || tempGeoTransform[1] != 1.0 || tempGeoTransform[2] != 0.0 || tempGeoTransform[3] != 0.0 || tempGeoTransform[4] != 0.0 || ABS(tempGeoTransform[5]) != 1.0 )) { OGRSpatialReference oUTMorLL; OGRSpatialReference oLL; OGRCoordinateTransformation *poTransform = NULL; char *srcProjection=NULL; char *newGCPProjection=NULL; padfTiepoints[0]=tempGeoTransform[0] + tempGeoTransform[1]*0.5 +\ tempGeoTransform[2]*0.5; padfTiepoints[1]=tempGeoTransform[3] + tempGeoTransform[4]*0.5 +\ tempGeoTransform[5]*0.5; padfTiepoints[2]=tempGeoTransform[0] + tempGeoTransform[2]*0.5 +\ tempGeoTransform[1]*(poSrcDS->GetRasterXSize()-0.5); padfTiepoints[3]=tempGeoTransform[3] + tempGeoTransform[5]*0.5 +\ tempGeoTransform[4]*(poSrcDS->GetRasterXSize()-0.5); padfTiepoints[4]=tempGeoTransform[0] + tempGeoTransform[1]*0.5 +\ tempGeoTransform[2]*(poSrcDS->GetRasterYSize()-0.5); padfTiepoints[5]=tempGeoTransform[3] + tempGeoTransform[4]*0.5 +\ tempGeoTransform[5]*(poSrcDS->GetRasterYSize()-0.5); padfTiepoints[6]=tempGeoTransform[0] +\ tempGeoTransform[1]*(poSrcDS->GetRasterXSize()-0.5) +\ tempGeoTransform[2]*(poSrcDS->GetRasterYSize()-0.5); padfTiepoints[7]=tempGeoTransform[3]+\ tempGeoTransform[4]*(poSrcDS->GetRasterXSize()-0.5)+\ tempGeoTransform[5]*(poSrcDS->GetRasterYSize()-0.5); padfTiepoints[8]=tempGeoTransform[0]+\ tempGeoTransform[1]*(poSrcDS->GetRasterXSize())/2.0+\ tempGeoTransform[2]*(poSrcDS->GetRasterYSize())/2.0; padfTiepoints[9]=tempGeoTransform[3]+\ tempGeoTransform[4]*(poSrcDS->GetRasterXSize())/2.0+\ tempGeoTransform[5]*(poSrcDS->GetRasterYSize())/2.0; srcProjection = CPLStrdup(poSrcDS->GetProjectionRef()); oUTMorLL.importFromWkt(&srcProjection); (oUTMorLL.GetAttrNode("GEOGCS"))->exportToWkt(&newGCPProjection); oLL.importFromWkt(&newGCPProjection); if EQUALN(poSrcDS->GetProjectionRef(),"PROJCS",6) { // projected coordinate system- need to translate gcps */ int bSuccess=TRUE; int index; poTransform = OGRCreateCoordinateTransformation( &oUTMorLL, &oLL ); if( poTransform == NULL ) bSuccess = FALSE; for (index=0;index<5;index++) { if( !bSuccess || !poTransform->Transform( 1, &(padfTiepoints[index*2]), &(padfTiepoints[index*2+1]) ) ) bSuccess = FALSE; } if (bSuccess == TRUE) georef_created = TRUE; } else { georef_created = TRUE; } } CPLFree(tempGeoTransform); } if (georef_created == TRUE) { char *szValue; szValue = (char *) CPLMalloc(255); /* -------------------------------------------------------------------- */ /* top left */ /* -------------------------------------------------------------------- */ sprintf( szValue, "TOP_LEFT_CORNER_LATITUDE = %.10f\n", padfTiepoints[1] ); fprintf( fp, szValue ); sprintf( szValue, "TOP_LEFT_CORNER_LONGITUDE = %.10f\n", padfTiepoints[0] ); fprintf( fp, szValue ); /* -------------------------------------------------------------------- */ /* top_right */ /* -------------------------------------------------------------------- */ sprintf( szValue, "TOP_RIGHT_CORNER_LATITUDE = %.10f\n", padfTiepoints[3] ); fprintf( fp, szValue ); sprintf( szValue, "TOP_RIGHT_CORNER_LONGITUDE = %.10f\n", padfTiepoints[2] ); fprintf( fp, szValue ); /* -------------------------------------------------------------------- */ /* bottom_left */ /* -------------------------------------------------------------------- */ sprintf( szValue, "BOTTOM_LEFT_CORNER_LATITUDE = %.10f\n", padfTiepoints[5] ); fprintf( fp, szValue ); sprintf( szValue, "BOTTOM_LEFT_CORNER_LONGITUDE = %.10f\n", padfTiepoints[4] ); fprintf( fp, szValue ); /* -------------------------------------------------------------------- */ /* bottom_right */ /* -------------------------------------------------------------------- */ sprintf( szValue, "BOTTOM_RIGHT_CORNER_LATITUDE = %.10f\n", padfTiepoints[7] ); fprintf( fp, szValue ); sprintf( szValue, "BOTTOM_RIGHT_CORNER_LONGITUDE = %.10f\n", padfTiepoints[6] ); fprintf( fp, szValue ); /* -------------------------------------------------------------------- */ /* Center */ /* -------------------------------------------------------------------- */ sprintf( szValue, "CENTRE_LATITUDE = %.10f\n", padfTiepoints[9] ); fprintf( fp, szValue ); sprintf( szValue, "CENTRE_LONGITUDE = %.10f\n", padfTiepoints[8] ); fprintf( fp, szValue ); CPLFree(szValue); /* ------------------------------------------------------------------- */ /* Ellipsoid/projection */ /* --------------------------------------------------------------------*/ OGRSpatialReference *oSRS; MFFSpheroidList *mffEllipsoids; double eq_radius, inv_flattening; OGRErr ogrerrorEq=OGRERR_NONE; OGRErr ogrerrorInvf=OGRERR_NONE; OGRErr ogrerrorOl=OGRERR_NONE; char *modifiableProjection = NULL; char *pszNewProjection; char *spheroid_name = NULL; pszNewProjection = CPLStrdup( poSrcDS->GetProjectionRef() ); if( !EQUALN(pszNewProjection,"GEOGCS",6) && !EQUALN(pszNewProjection,"PROJCS",6) && !EQUAL(pszNewProjection,"") ) { CPLError( CE_Warning, CPLE_AppDefined, "Only OGC WKT Projections supported for writing to MFF.\n" "%s not supported.", pszNewProjection ); } else if (!EQUAL(pszNewProjection,"")) { /* importFromWkt updates the pointer, so don't use pszNewProjection directly */ modifiableProjection=CPLStrdup(pszNewProjection); oSRS = new OGRSpatialReference; oSRS->importFromWkt(&modifiableProjection); if ((oSRS->GetAttrValue("PROJECTION") != NULL) && (EQUAL(oSRS->GetAttrValue("PROJECTION"),SRS_PT_TRANSVERSE_MERCATOR))) { char *ol_txt; ol_txt=(char *) CPLMalloc(255); fprintf(fp,"PROJECTION_NAME = UTM\n"); sprintf(ol_txt,"PROJECTION_ORIGIN_LONGITUDE = %f\n", oSRS->GetProjParm(SRS_PP_CENTRAL_MERIDIAN,0.0,&ogrerrorOl)); fprintf(fp,ol_txt); CPLFree(ol_txt); } else if ((oSRS->GetAttrValue("PROJECTION") == NULL) && (oSRS->IsGeographic())) { fprintf(fp,"PROJECTION_NAME = LL\n"); } else { CPLError( CE_Warning, CPLE_AppDefined, "Unrecognized projection- no georeferencing information transferred."); fprintf(fp,"PROJECTION_NAME = LL\n"); } eq_radius = oSRS->GetSemiMajor(&ogrerrorEq); inv_flattening = oSRS->GetInvFlattening(&ogrerrorInvf); if ((ogrerrorEq == OGRERR_NONE) && (ogrerrorInvf == OGRERR_NONE)) { char *ol_txt; ol_txt=(char *) CPLMalloc(255); mffEllipsoids = new MFFSpheroidList; spheroid_name = mffEllipsoids->GetSpheroidNameByEqRadiusAndInvFlattening(eq_radius,inv_flattening); if (spheroid_name != NULL) { sprintf(ol_txt,"SPHEROID_NAME = %s\n",spheroid_name ); fprintf(fp,ol_txt); } else { sprintf(ol_txt, "SPHEROID_NAME = USER_DEFINED\nSPHEROID_EQUATORIAL_RADIUS = %.10f\nSPHEROID_POLAR_RADIUS = %.10f\n", eq_radius,eq_radius*(1-1.0/inv_flattening) ); fprintf(fp,ol_txt); } delete mffEllipsoids; CPLFree(ol_txt); CPLFree(spheroid_name); } } CPLFree(pszNewProjection); } CPLFree( padfTiepoints ); fprintf( fp, "END\n" ); VSIFClose( fp ); /* End of georeferencing stuff */ /* Make sure image data gets flushed */ for( iBand = 0; iBand < poDS->GetRasterCount(); iBand++ ) { RawRasterBand *poDstBand = (RawRasterBand *) poDS->GetRasterBand( iBand+1 ); poDstBand->FlushCache(); } if( !pfnProgress( 1.0, NULL, pProgressData ) ) { CPLError( CE_Failure, CPLE_UserInterrupt, "User terminated" ); delete poDS; GDALDriver *poMFFDriver = (GDALDriver *) GDALGetDriverByName( "MFF" ); poMFFDriver->Delete( pszFilename ); return NULL; } poDS->CloneInfo( poSrcDS, GCIF_PAM_DEFAULT ); return poDS; } /************************************************************************/ /* GDALRegister_MFF() */ /************************************************************************/ void GDALRegister_MFF() { GDALDriver *poDriver; if( GDALGetDriverByName( "MFF" ) == NULL ) { poDriver = new GDALDriver(); poDriver->SetDescription( "MFF" ); poDriver->SetMetadataItem( GDAL_DMD_LONGNAME, "Vexcel MFF Raster" ); poDriver->SetMetadataItem( GDAL_DMD_HELPTOPIC, "frmt_various.html#MFF" ); poDriver->SetMetadataItem( GDAL_DMD_EXTENSION, "hdr" ); poDriver->SetMetadataItem( GDAL_DMD_CREATIONDATATYPES, "Byte UInt16 Float32 CInt16 CFloat32" ); poDriver->pfnOpen = MFFDataset::Open; poDriver->pfnCreate = MFFDataset::Create; poDriver->pfnCreateCopy = MFFDataset::CreateCopy; GetGDALDriverManager()->RegisterDriver( poDriver ); } }