自学教程：C++ GDALOpen函数代码示例

// -------------------------------------------------------------------int main(int argc, char **argv){	   char *ERROR_PRM = "Ungueltiger numerischer Wert fuer %s: %s/n!";	   // Alle GDAL Treiber registrieren      GDALAllRegister();	   // Dateiname der Geotiff-Datei        const char *format = "GTiff";   GDALDriverH h_drv = GDALGetDriverByName( format );   if( h_drv == NULL ) {	   error_exit(10,"Treiber %s nicht vorhanden!" ,format);   }   // Test ob Geotiffdateien erzeugt werden koennen   char **test_meta;   test_meta = GDALGetMetadata( h_drv, NULL );   if( ! CSLFetchBoolean( test_meta, GDAL_DCAP_CREATE, FALSE ) ) {	   error_exit(10,"Das Format %s kann nicht erzeugt werden" ,format);	}         // 3 Kommandozeilenparameter   if (argc<6) {	   error_exit(10,	   "Fehlende Parameter/nUsage %s IN OUT EXT SZ ID X Y ID X Y ID X Y....!/n",	   argv[0]);   }      // Eingabemuster einlesen   char *ifile = argv[1];   // Ausgabemuster einlesen   char *ofile = argv[2];   // zusammengesetzte Ausgabedatei    char cfile[512];   // Dateierweiterung setzen   char *ext   = argv[3];   // Fenstergroesse   int size = 64;    if (! sscanf(argv[4],"%d",&size) ) {        error_exit(1000+3,ERROR_PRM,"SZ",argv[4]);   }   double trfm[] ={0,0,0,0,0,0};   // Vektoren fuer die Positionen und ID   int_vector_t id;   int_vector_init(&id, 10);   dbl_vector_t pos_x;   dbl_vector_init(&pos_x, 10);   dbl_vector_t pos_y;   dbl_vector_init(&pos_y, 10);      // Positionen einlesen    int a = 5; double dbl; int pk;   while( a < argc-2 ) {	   // X Koordinate parsen       if (! sscanf(argv[a],"%d",&pk) ) {          error_exit(1000+a,ERROR_PRM,"ID", argv[a]);       }	   int_vector_add(&id,pk);       	   // X Koordinate parsen       if (! sscanf(argv[a+1],"%lf",&dbl) ) {          error_exit(1000+a+1,ERROR_PRM,"X", argv[a+1]);       }	   dbl_vector_add(&pos_x,dbl);	   // Y Koordinate parsen       if (! sscanf(argv[a+2],"%lf",&dbl) ) {          error_exit(1000+a+2,ERROR_PRM,"Y", argv[a+2]);       }	   dbl_vector_add(&pos_y,dbl);	   a+=3;   }         // Alle GDAL Treiber registrieren      GDALAllRegister();      // Geotiff oeffnen   printf("# IN FILE:  %s/n", ifile);   GDALDatasetH h_dset = GDALOpen( ifile, GA_ReadOnly);   printf("# OUT FILE: %s/n", ofile);   printf("# EXTENTION: %s/n", ext);      // Transformation holen    if( GDALGetGeoTransform( h_dset, trfm ) == CE_None ) {        printf("# TRANSFORM: /n");        printf("#  X = %.6f + %.6f * COL + %.6f * ROW/n",                   trfm[0], trfm[1], trfm[2] );        printf("#  Y = %.6f + %.6f * COL + %.6f * ROW/n# EOF:/n",                   trfm[3], trfm[4], trfm[5] );//.........这里部分代码省略.........

//.........这里部分代码省略.........                                     char ** papszOptions ){    if( eType != GDT_Float32 )    {        CPLError(CE_Failure, CPLE_AppDefined,                 "Attempt to create CTable2 file with unsupported data type '%s'.",                 GDALGetDataTypeName( eType ) );        return NULL;    }/* -------------------------------------------------------------------- *//*      Try to open or create file.                                     *//* -------------------------------------------------------------------- */    VSILFILE	*fp;    fp = VSIFOpenL( pszFilename, "wb" );        if( fp == NULL )    {        CPLError( CE_Failure, CPLE_OpenFailed,                  "Attempt to create file `%s' failed./n",                  pszFilename );        return NULL;    }/* -------------------------------------------------------------------- *//*      Create a file header, with a defaulted georeferencing.          *//* -------------------------------------------------------------------- */    char achHeader[160];    int nValue32;    double dfValue;    memset( achHeader, 0, sizeof(achHeader));    memcpy( achHeader+0, "CTABLE V2.0     ", 16 );        if( CSLFetchNameValue( papszOptions, "DESCRIPTION" ) != NULL )        strncpy( achHeader + 16,                  CSLFetchNameValue( papszOptions, "DESCRIPTION" ),                  80 );        // lower left origin (longitude, center of pixel, radians)    dfValue = 0;    CPL_LSBPTR64( &dfValue );    memcpy( achHeader + 96, &dfValue, 8 );    // lower left origin (latitude, center of pixel, radians)    dfValue = 0;    CPL_LSBPTR64( &dfValue );    memcpy( achHeader + 104, &dfValue, 8 );    // pixel width (radians)    dfValue = 0.01 * M_PI / 180.0;    CPL_LSBPTR64( &dfValue );    memcpy( achHeader + 112, &dfValue, 8 );    // pixel height (radians)    dfValue = 0.01 * M_PI / 180.0;    CPL_LSBPTR64( &dfValue );    memcpy( achHeader + 120, &dfValue, 8 );    // raster width in pixels    nValue32 = nXSize;    CPL_LSBPTR32( &nValue32 );    memcpy( achHeader + 128, &nValue32, 4 );    // raster width in pixels    nValue32 = nYSize;    CPL_LSBPTR32( &nValue32 );    memcpy( achHeader + 132, &nValue32, 4 );    VSIFWriteL( achHeader, 1, sizeof(achHeader), fp );/* -------------------------------------------------------------------- *//*      Write zeroed grid data.                                         *//* -------------------------------------------------------------------- */    float *pafLine = (float *) CPLCalloc(sizeof(float)*2,nXSize);    int i;    for( i = 0; i < nYSize; i++ )    {        if( (int)VSIFWriteL( pafLine, sizeof(float)*2, nXSize, fp ) != nXSize )         {            CPLError( CE_Failure, CPLE_FileIO,                       "Write failed at line %d, perhaps the disk is full?",                      i );            return NULL;        }    }    /* -------------------------------------------------------------------- *//*      Cleanup and return.                                             *//* -------------------------------------------------------------------- */    CPLFree( pafLine );    VSIFCloseL( fp );    return (GDALDataset *) GDALOpen( pszFilename, GA_Update );}

int main( int argc, char ** argv ){    /* Check that we are running against at least GDAL 1.4 (probably older in fact !) */    /* Note to developers : if we use newer API, please change the requirement */    if (atoi(GDALVersionInfo("VERSION_NUM")) < 1400)    {        fprintf(stderr, "At least, GDAL >= 1.4.0 is required for this version of %s, "                "which was compiled against GDAL %s/n", argv[0], GDAL_RELEASE_NAME);        exit(1);    }/* -------------------------------------------------------------------- *//*      Generic arg processing.                                         *//* -------------------------------------------------------------------- */    GDALAllRegister();    GDALSetCacheMax( 100000000 );    argc = GDALGeneralCmdLineProcessor( argc, &argv, 0 );    if( argc < 1 )        exit( -argc );    /* -------------------------------------------------------------------- *//*      Parse arguments.                                                *//* -------------------------------------------------------------------- */    int i;    const char *pszOutFile = NULL;    const char *pszInFile = NULL;    int nMaxNonBlack = 2;    int nNearDist = 15;    int bNearWhite = FALSE;        for( i = 1; i < argc; i++ )    {        if( EQUAL(argv[i], "--utility_version") )        {            printf("%s was compiled against GDAL %s and is running against GDAL %s/n",                   argv[0], GDAL_RELEASE_NAME, GDALVersionInfo("RELEASE_NAME"));            return 0;        }        else if( EQUAL(argv[i], "-o") && i < argc-1 )            pszOutFile = argv[++i];        else if( EQUAL(argv[i], "-white") )            bNearWhite = TRUE;        else if( EQUAL(argv[i], "-nb") && i < argc-1 )            nMaxNonBlack = atoi(argv[++i]);        else if( EQUAL(argv[i], "-near") && i < argc-1 )            nNearDist = atoi(argv[++i]);        else if( argv[i][0] == '-' )            Usage();        else if( pszInFile == NULL )            pszInFile = argv[i];        else            Usage();    }    if( pszInFile == NULL )        Usage();    if( pszOutFile == NULL )        pszOutFile = pszInFile;/* -------------------------------------------------------------------- *//*      Open input file.                                                *//* -------------------------------------------------------------------- */    GDALDatasetH hInDS, hOutDS = NULL;    int nXSize, nYSize, nBands;    if( pszOutFile == pszInFile )        hInDS = hOutDS = GDALOpen( pszInFile, GA_Update );    else        hInDS = GDALOpen( pszInFile, GA_ReadOnly );    if( hInDS == NULL )        exit( 1 );    nXSize = GDALGetRasterXSize( hInDS );    nYSize = GDALGetRasterYSize( hInDS );    nBands = GDALGetRasterCount( hInDS );/* -------------------------------------------------------------------- *//*      Do we need to create output file?                               *//* -------------------------------------------------------------------- */    if( hOutDS == NULL )    {        GDALDriverH hDriver = GDALGetDriverByName( "HFA" );                hOutDS = GDALCreate( hDriver, pszOutFile,                              nXSize, nYSize, nBands, GDT_Byte,                              NULL );        if( hOutDS == NULL )            exit( 1 );        double adfGeoTransform[6];        if( GDALGetGeoTransform( hInDS, adfGeoTransform ) == CE_None )        {            GDALSetGeoTransform( hOutDS, adfGeoTransform );            GDALSetProjection( hOutDS, GDALGetProjectionRef( hInDS ) );        }    }//.........这里部分代码省略.........

//.........这里部分代码省略.........    double dfMinY = adfGeoTransform[5] * (nYSize - 0.5) + adfGeoTransform[3];    double dfMaxY = adfGeoTransform[3] + adfGeoTransform[5] / 2;    CPLErr eErr = WriteHeader( fp, nXSize, nYSize,			       dfMinX, dfMaxX, dfMinY, dfMaxY, 0.0, 0.0 );    if( eErr != CE_None )    {	VSIFCloseL( fp );        return NULL;    }/* -------------------------------------------------------------------- *//*      Copy band data.							*//* -------------------------------------------------------------------- */    float *pfData = (float *)VSIMalloc2( nXSize, sizeof( float ) );    if( pfData == NULL )    {	VSIFCloseL( fp );	CPLError( CE_Failure, CPLE_OutOfMemory,		  "Unable to create copy, unable to allocate line buffer./n" );	return NULL;    }    int     bSrcHasNDValue;    float   fSrcNoDataValue = (float) poSrcBand->GetNoDataValue( &bSrcHasNDValue );    double  dfMinZ = DBL_MAX;    double  dfMaxZ = -DBL_MAX;    for( GInt16 iRow = nYSize - 1; iRow >= 0; iRow-- )    {	eErr = poSrcBand->RasterIO( GF_Read, 0, iRow,				    nXSize, 1, pfData,				    nXSize, 1, GDT_Float32, 0, 0 );	if( eErr != CE_None )	{	    VSIFCloseL( fp );	    VSIFree( pfData );	    return NULL;	}	for( int iCol=0; iCol<nXSize; iCol++ )	{	    if( bSrcHasNDValue && pfData[iCol] == fSrcNoDataValue )	    {		pfData[iCol] = fNODATA_VALUE;	    }	    else	    {		if( pfData[iCol] > dfMaxZ )		    dfMaxZ = pfData[iCol];		if( pfData[iCol] < dfMinZ )		    dfMinZ = pfData[iCol];	    }	    CPL_LSBPTR32( pfData+iCol );	}	if( VSIFWriteL( (void *)pfData, 4, nXSize,			fp ) != static_cast<unsigned>(nXSize) )	{	    VSIFCloseL( fp );	    VSIFree( pfData );	    CPLError( CE_Failure, CPLE_FileIO,		      "Unable to write grid row. Disk full?/n" );	    return NULL;	}	if( !pfnProgress( static_cast<double>(nYSize - iRow)/nYSize,			  NULL, pProgressData ) )	{	    VSIFCloseL( fp );	    VSIFree( pfData );	    CPLError( CE_Failure, CPLE_UserInterrupt, "User terminated" );	    return NULL;	}    }    VSIFree( pfData );    /* write out the min and max values */    eErr = WriteHeader( fp, nXSize, nYSize,			dfMinX, dfMaxX, dfMinY, dfMaxY, dfMinZ, dfMaxZ );    if( eErr != CE_None )    {	VSIFCloseL( fp );        return NULL;    }    VSIFCloseL( fp );    GDALPamDataset *poDS = (GDALPamDataset *)GDALOpen( pszFilename,                                                GA_Update );    if (poDS)    {        poDS->CloneInfo( poSrcDS, GCIF_PAM_DEFAULT );    }    return poDS;}

//.........这里部分代码省略.........            {                char** papszMEMDSOptions = NULL;                char szTmp[64];                memset(szTmp, 0, sizeof(szTmp));                CPLPrintPointer(szTmp,                                pabyMEMDSBuffer + iBand * nDataTypeSize *                                nReqXSize * nReqYSize, sizeof(szTmp));                papszMEMDSOptions = CSLSetNameValue(papszMEMDSOptions, "DATAPOINTER", szTmp);                GDALAddBand(hMemDS, eDataType, papszMEMDSOptions);                CSLDestroy(papszMEMDSOptions);            }                        GDALDatasetH hOutDS = GDALCreateCopy(hTileDriver,                                        osTempFileName.c_str(), hMemDS, FALSE,                                        papszTileDriverOptions, NULL, NULL);            GDALClose(hMemDS);            if (hOutDS)                GDALClose(hOutDS);            else            {                eErr = CE_Failure;                break;            }/* -------------------------------------------------------------------- *//*      Insert new entry into raster table                              *//* -------------------------------------------------------------------- */            vsi_l_offset nDataLength;            GByte *pabyData = VSIGetMemFileBuffer( osTempFileName.c_str(),                                                   &nDataLength, FALSE);            OGRFeatureH hFeat = OGR_F_Create( OGR_L_GetLayerDefn(hRasterLayer) );            OGR_F_SetFieldBinary(hFeat, 0, (int)nDataLength, pabyData);                        OGR_L_CreateFeature(hRasterLayer, hFeat);            /* Query raster ID to set it as the ID of the associated metadata */            int nRasterID = (int)OGR_F_GetFID(hFeat);                        OGR_F_Destroy(hFeat);                        VSIUnlink(osTempFileName.c_str());            /* -------------------------------------------------------------------- *//*      Insert new entry into metadata table                            *//* -------------------------------------------------------------------- */                        hFeat = OGR_F_Create( OGR_L_GetLayerDefn(hMetadataLayer) );            OGR_F_SetFID(hFeat, nRasterID);            OGR_F_SetFieldString(hFeat, 0, GDALGetDescription(poSrcDS));            OGR_F_SetFieldInteger(hFeat, 1, nTileId ++);            OGR_F_SetFieldInteger(hFeat, 2, nReqXSize);            OGR_F_SetFieldInteger(hFeat, 3, nReqYSize);            OGR_F_SetFieldDouble(hFeat, 4, adfGeoTransform[1]);            OGR_F_SetFieldDouble(hFeat, 5, -adfGeoTransform[5]);                        minx = adfGeoTransform[0] +                (nBlockXSize * nBlockXOff) * adfGeoTransform[1];            maxx = adfGeoTransform[0] +                (nBlockXSize * nBlockXOff + nReqXSize) * adfGeoTransform[1];            maxy = adfGeoTransform[3] +                (nBlockYSize * nBlockYOff) * adfGeoTransform[5];            miny = adfGeoTransform[3] +                (nBlockYSize * nBlockYOff + nReqYSize) * adfGeoTransform[5];                        OGRGeometryH hRectangle = OGR_G_CreateGeometry(wkbPolygon);            OGRGeometryH hLinearRing = OGR_G_CreateGeometry(wkbLinearRing);            OGR_G_AddPoint_2D(hLinearRing, minx, miny);            OGR_G_AddPoint_2D(hLinearRing, minx, maxy);            OGR_G_AddPoint_2D(hLinearRing, maxx, maxy);            OGR_G_AddPoint_2D(hLinearRing, maxx, miny);            OGR_G_AddPoint_2D(hLinearRing, minx, miny);            OGR_G_AddGeometryDirectly(hRectangle, hLinearRing);                        OGR_F_SetGeometryDirectly(hFeat, hRectangle);                        OGR_L_CreateFeature(hMetadataLayer, hFeat);            OGR_F_Destroy(hFeat);                        nBlocks++;            if (pfnProgress && !pfnProgress(1.0 * nBlocks / nTotalBlocks,                                            NULL, pProgressData))                eErr = CE_Failure;        }    }        if (eErr == CE_None)        OGR_DS_ExecuteSQL(hDS, "COMMIT", NULL, NULL);    else        OGR_DS_ExecuteSQL(hDS, "ROLLBACK", NULL, NULL);        CSLDestroy(papszTileDriverOptions);        VSIFree(pabyMEMDSBuffer);        OGRReleaseDataSource(hDS);            return (GDALDataset*) GDALOpen(pszFilename, GA_Update);}

//.........这里部分代码省略.........		    if( dfValue > dfMax )			dfMax = dfValue;		    if( dfValue < dfMin )			dfMin = dfValue;		}		std::ostringstream ssOut;		ssOut.precision(nFIELD_PRECISION);		ssOut.setf( std::ios::uppercase );		ssOut << dfValue << " ";		CPLString sOut = ssOut.str();		if( VSIFWriteL( sOut.c_str(), 1, sOut.length(), fp )		    != sOut.length() )		{		    VSIFCloseL( fp );		    VSIFree( pdfData );		    CPLError( CE_Failure, CPLE_FileIO,			      "Unable to write grid cell.  Disk full?/n" );		    return NULL;		}	    }	    if( VSIFWriteL( (void *)"/x0D/x0A", 1, 2, fp ) != 2 )	    {		VSIFCloseL( fp );		VSIFree( pdfData );		CPLError( CE_Failure, CPLE_FileIO,			  "Unable to finish write of grid line. Disk full?/n" );		return NULL;	    }	}	if( VSIFWriteL( (void *)"/x0D/x0A", 1, 2, fp ) != 2 )	{	    VSIFCloseL( fp );	    VSIFree( pdfData );	    CPLError( CE_Failure, CPLE_FileIO,		      "Unable to finish write of grid row. Disk full?/n" );	    return NULL;	}	if( !pfnProgress( static_cast<double>(iRow + 1)/nYSize,			  NULL, pProgressData ) )	{	    VSIFCloseL( fp );	    VSIFree( pdfData );	    CPLError( CE_Failure, CPLE_UserInterrupt, "User terminated" );	    return NULL;	}    }    VSIFree( pdfData );    /* write out the min and max values */    std::ostringstream ssRange;    ssRange.precision( nFIELD_PRECISION );    ssRange.setf( std::ios::uppercase );    ssRange << dfMin << " " << dfMax << "/x0D/x0A";    if( ssRange.str().length() != nDummyRangeLen )    {	int nShiftSize = ssRange.str().length() - nDummyRangeLen;	if( ShiftFileContents( fp, nRangeStart + nDummyRangeLen,			       nShiftSize, "/x0D/x0A" ) != CE_None )	{	    VSIFCloseL( fp );	    CPLError( CE_Failure, CPLE_FileIO,		      "Unable to shift file contents./n" );	    return NULL;	}    }    if( VSIFSeekL( fp, nRangeStart, SEEK_SET ) != 0 )    {	VSIFCloseL( fp );	CPLError( CE_Failure, CPLE_FileIO,		  "Unable to seek to start of grid file copy./n" );	return NULL;    }    if( VSIFWriteL( (void *)ssRange.str().c_str(), 1, ssRange.str().length(),		    fp ) != ssRange.str().length() )    {	VSIFCloseL( fp );        CPLError( CE_Failure, CPLE_FileIO,                  "Unable to write range information./n" );        return NULL;    }    VSIFCloseL( fp );    GDALPamDataset *poDS = (GDALPamDataset *)GDALOpen( pszFilename,                                                GA_Update );    if (poDS)    {        poDS->CloneInfo( poSrcDS, GCIF_PAM_DEFAULT );    }    return poDS;}

int main( int nArgc, char ** papszArgv ){    GDALDatasetH     hDataset;    const char      *pszResampling = "nearest";    const char      *pszFilename = NULL;    int              anLevels[1024];    int              nLevelCount = 0;    int              nResultStatus = 0;    int              bReadOnly = FALSE;    int              bClean = FALSE;    GDALProgressFunc pfnProgress = GDALTermProgress;     /* Check that we are running against at least GDAL 1.7 */    /* Note to developers : if we use newer API, please change the requirement */    if (atoi(GDALVersionInfo("VERSION_NUM")) < 1700)    {        fprintf(stderr, "At least, GDAL >= 1.7.0 is required for this version of %s, "                        "which was compiled against GDAL %s/n", papszArgv[0], GDAL_RELEASE_NAME);        exit(1);    }    GDALAllRegister();    nArgc = GDALGeneralCmdLineProcessor( nArgc, &papszArgv, 0 );    if( nArgc < 1 )        exit( -nArgc );/* -------------------------------------------------------------------- *//*      Parse commandline.                                              *//* -------------------------------------------------------------------- */    for( int iArg = 1; iArg < nArgc; iArg++ )    {        if( EQUAL(papszArgv[iArg], "--utility_version") )        {            printf("%s was compiled against GDAL %s and is running against GDAL %s/n",                   papszArgv[0], GDAL_RELEASE_NAME, GDALVersionInfo("RELEASE_NAME"));            return 0;        }        else if( EQUAL(papszArgv[iArg],"-r") && iArg < nArgc-1 )            pszResampling = papszArgv[++iArg];        else if( EQUAL(papszArgv[iArg],"-ro"))            bReadOnly = TRUE;        else if( EQUAL(papszArgv[iArg],"-clean"))            bClean = TRUE;        else if( EQUAL(papszArgv[iArg],"-q") || EQUAL(papszArgv[iArg],"-quiet") )             pfnProgress = GDALDummyProgress;         else if( pszFilename == NULL )            pszFilename = papszArgv[iArg];        else if( atoi(papszArgv[iArg]) > 0 )            anLevels[nLevelCount++] = atoi(papszArgv[iArg]);        else            Usage();    }    if( pszFilename == NULL || (nLevelCount == 0 && !bClean) )        Usage();/* -------------------------------------------------------------------- *//*      Open data file.                                                 *//* -------------------------------------------------------------------- */    if (bReadOnly)        hDataset = NULL;    else    {        CPLPushErrorHandler( CPLQuietErrorHandler );        hDataset = GDALOpen( pszFilename, GA_Update );        CPLPopErrorHandler();    }    if( hDataset == NULL )        hDataset = GDALOpen( pszFilename, GA_ReadOnly );    if( hDataset == NULL )        exit( 2 );/* -------------------------------------------------------------------- *//*      Clean overviews.                                                *//* -------------------------------------------------------------------- */    if ( bClean &&        GDALBuildOverviews( hDataset,pszResampling, 0, 0,                              0, NULL, pfnProgress, NULL ) != CE_None )    {        printf( "Cleaning overviews failed./n" );        nResultStatus = 200;    }/* -------------------------------------------------------------------- *//*      Generate overviews.                                             *//* -------------------------------------------------------------------- */    if (nLevelCount > 0 && nResultStatus == 0 &&        GDALBuildOverviews( hDataset,pszResampling, nLevelCount, anLevels,                             0, NULL, pfnProgress, NULL ) != CE_None )    {        printf( "Overview building failed./n" );        nResultStatus = 100;    }/* -------------------------------------------------------------------- *//*      Cleanup                                                         *///.........这里部分代码省略.........

CPLErr GDALWMSRasterBand::ReadBlockFromFile(int x, int y, const char *file_name, int to_buffer_band, void *buffer, int advise_read) {    CPLErr ret = CE_None;    GDALDataset *ds = 0;    GByte *color_table = NULL;    int i;    //CPLDebug("WMS", "ReadBlockFromFile: to_buffer_band=%d, (x,y)=(%d, %d)", to_buffer_band, x, y);    /* expected size */    const int esx = MIN(MAX(0, (x + 1) * nBlockXSize), nRasterXSize) - MIN(MAX(0, x * nBlockXSize), nRasterXSize);    const int esy = MIN(MAX(0, (y + 1) * nBlockYSize), nRasterYSize) - MIN(MAX(0, y * nBlockYSize), nRasterYSize);    ds = reinterpret_cast<GDALDataset*>(GDALOpen(file_name, GA_ReadOnly));    if (ds != NULL) {        int sx = ds->GetRasterXSize();        int sy = ds->GetRasterYSize();        bool accepted_as_no_alpha = false;  // if the request is for 4 bands but the wms returns 3          /* Allow bigger than expected so pre-tiled constant size images work on corners */        if ((sx > nBlockXSize) || (sy > nBlockYSize) || (sx < esx) || (sy < esy)) {            CPLError(CE_Failure, CPLE_AppDefined, "GDALWMS: Incorrect size %d x %d of downloaded block, expected %d x %d, max %d x %d.",                sx, sy, esx, esy, nBlockXSize, nBlockYSize);            ret = CE_Failure;        }        if (ret == CE_None) {            int nDSRasterCount = ds->GetRasterCount();            if (nDSRasterCount != m_parent_dataset->nBands) {                /* Maybe its an image with color table */                bool accepted_as_ct = false;                if ((eDataType == GDT_Byte) && (ds->GetRasterCount() == 1)) {                    GDALRasterBand *rb = ds->GetRasterBand(1);                    if (rb->GetRasterDataType() == GDT_Byte) {                        GDALColorTable *ct = rb->GetColorTable();                        if (ct != NULL) {                            accepted_as_ct = true;                            if (!advise_read) {                                color_table = new GByte[256 * 4];                                const int count = MIN(256, ct->GetColorEntryCount());                                for (i = 0; i < count; ++i) {                                    GDALColorEntry ce;                                    ct->GetColorEntryAsRGB(i, &ce);                                    color_table[i] = static_cast<GByte>(ce.c1);                                    color_table[i + 256] = static_cast<GByte>(ce.c2);                                    color_table[i + 512] = static_cast<GByte>(ce.c3);                                    color_table[i + 768] = static_cast<GByte>(ce.c4);                                }                                for (i = count; i < 256; ++i) {                                    color_table[i] = 0;                                    color_table[i + 256] = 0;                                    color_table[i + 512] = 0;                                    color_table[i + 768] = 0;                                }                            }                        }                    }                }                if (nDSRasterCount == 4 && m_parent_dataset->nBands == 3)                {                    /* metacarta TMS service sometimes return a 4 band PNG instead of the expected 3 band... */                }                else if (!accepted_as_ct) {                   if (ds->GetRasterCount()==3 && m_parent_dataset->nBands == 4 && (eDataType == GDT_Byte))                   { // WMS returned a file with no alpha so we will fill the alpha band with "opaque"                       accepted_as_no_alpha = true;                   }                   else                   {                      CPLError(CE_Failure, CPLE_AppDefined, "GDALWMS: Incorrect bands count %d in downloaded block, expected %d.",                         nDSRasterCount, m_parent_dataset->nBands);                      ret = CE_Failure;                   }                }            }        }        if (!advise_read) {            for (int ib = 1; ib <= m_parent_dataset->nBands; ++ib) {                if (ret == CE_None) {                    void *p = NULL;                    GDALRasterBlock *b = NULL;                    if ((buffer != NULL) && (ib == to_buffer_band)) {                        p = buffer;                    } else {                        GDALWMSRasterBand *band = static_cast<GDALWMSRasterBand *>(m_parent_dataset->GetRasterBand(ib));                        if (m_overview >= 0) band = static_cast<GDALWMSRasterBand *>(band->GetOverview(m_overview));                        if (!band->IsBlockInCache(x, y)) {                            b = band->GetLockedBlockRef(x, y, true);                            if (b != NULL) {                                p = b->GetDataRef();                                if (p == NULL) {                                  CPLError(CE_Failure, CPLE_AppDefined, "GDALWMS: GetDataRef returned NULL.");                                  ret = CE_Failure;                                }                            }                        }                        else                        {                            //CPLDebug("WMS", "Band %d, block (x,y)=(%d, %d) already in cache", band->GetBand(), x, y);                        }                    }                    if (p != NULL) {                        int pixel_space = GDALGetDataTypeSize(eDataType) / 8;//.........这里部分代码省略.........

//		GDALRegister_ADRG();		GDALRegister_COASP();		GDALRegister_BLX();		GDALRegister_LCP();//		GDALRegister_PGCHIP();//		GDALRegister_TMS();		GDALRegister_EIR();//		GDALRegister_GEOR();        MemoryIgnoreLeaksEnd();    } 	Close();		MemoryIgnoreLeaksBegin();	dataset = (GDALDataset *) GDALOpen(fn, GA_ReadOnly);	MemoryIgnoreLeaksEnd();	if(!dataset)		return false;	int nbands = dataset->GetRasterCount();	for(int i = 1; i <= nbands; i++)		bands.Add(new Band(*this, dataset->GetRasterBand(i)));	double affine[6];	dataset->GetGeoTransform(affine);	transform.a.x = affine[0];	transform.x.x = affine[1];	transform.y.x = affine[2];	transform.a.y = affine[3];	transform.x.y = affine[4];	transform.y.y = affine[5];	if(fabs(Determinant(transform)) <= 1e-10)

static void msTransformToGeospatialPDF(imageObj *img, mapObj *map, cairo_renderer *r){  /* We need a GDAL 1.10 PDF driver at runtime, but as far as the C API is concerned, GDAL 1.9 is */  /* largely sufficient. */#if defined(USE_GDAL) && defined(GDAL_VERSION_NUM) && GDAL_VERSION_NUM >= 1900  GDALDatasetH hDS = NULL;  const char* pszGEO_ENCODING = NULL;  GDALDriverH hPDFDriver = NULL;  const char* pszVirtualIO = NULL;  int bVirtualIO = FALSE;  char* pszTmpFilename = NULL;  VSILFILE* fp = NULL;  if (map == NULL)    return;  pszGEO_ENCODING = msGetOutputFormatOption(img->format, "GEO_ENCODING", NULL);  if (pszGEO_ENCODING == NULL)    return;  msGDALInitialize();  hPDFDriver = GDALGetDriverByName("PDF");  if (hPDFDriver == NULL)    return;  /* When compiled against libpoppler, the PDF driver is VirtualIO capable */  /* but not, when it is compiled against libpodofo. */  pszVirtualIO = GDALGetMetadataItem( hPDFDriver, GDAL_DCAP_VIRTUALIO, NULL );  if (pszVirtualIO)    bVirtualIO = CSLTestBoolean(pszVirtualIO);  if (bVirtualIO)    pszTmpFilename = msTmpFile(map, NULL, "/vsimem/mscairopdf/", "pdf");  else    pszTmpFilename = msTmpFile(map, map->mappath, NULL, "pdf");  /* Copy content of outputStream buffer into file */  fp = VSIFOpenL(pszTmpFilename, "wb");  if (fp == NULL) {    msFree(pszTmpFilename);    return;  }  VSIFWriteL(r->outputStream->data, 1, r->outputStream->size, fp);  VSIFCloseL(fp);  fp = NULL;  hDS = GDALOpen(pszTmpFilename, GA_Update);  if ( hDS != NULL ) {    char* pszWKT = msProjectionObj2OGCWKT( &(map->projection) );    if( pszWKT != NULL ) {      double adfGeoTransform[6];      int i;      /* Add user-specified options */      for( i = 0; i < img->format->numformatoptions; i++ ) {        const char* pszOption = img->format->formatoptions[i];        if( strncasecmp(pszOption,"METADATA_ITEM:",14) == 0 ) {          char* pszKey = NULL;          const char* pszValue = CPLParseNameValue(pszOption + 14,                                 &pszKey);          if( pszKey != NULL ) {            GDALSetMetadataItem(hDS, pszKey, pszValue, NULL);            CPLFree(pszKey);          }        }      }      /* We need to rescale the geotransform because GDAL will not necessary */      /* open the PDF with the DPI that was used to generate it */      memcpy(adfGeoTransform, map->gt.geotransform, 6 * sizeof(double));      adfGeoTransform[1] = adfGeoTransform[1] * map->width / GDALGetRasterXSize(hDS);      adfGeoTransform[5] = adfGeoTransform[5] * map->height / GDALGetRasterYSize(hDS);      GDALSetGeoTransform(hDS, adfGeoTransform);      GDALSetProjection(hDS, pszWKT);      msFree( pszWKT );      pszWKT = NULL;      CPLSetThreadLocalConfigOption("GDAL_PDF_GEO_ENCODING", pszGEO_ENCODING);      GDALClose(hDS);      hDS = NULL;      CPLSetThreadLocalConfigOption("GDAL_PDF_GEO_ENCODING", NULL);      /* We need to replace the buffer with the content of the GDAL file */      fp = VSIFOpenL(pszTmpFilename, "rb");      if( fp != NULL ) {        int nFileSize;        VSIFSeekL(fp, 0, SEEK_END);        nFileSize = (int)VSIFTellL(fp);        msBufferResize(r->outputStream, nFileSize);        VSIFSeekL(fp, 0, SEEK_SET);        r->outputStream->size = VSIFReadL(r->outputStream->data, 1, nFileSize, fp);        VSIFCloseL(fp);//.........这里部分代码省略.........

GDALDataset *PNMDataset::Create( const char * pszFilename,                                 int nXSize, int nYSize, int nBands,                                 GDALDataType eType,                                 char ** papszOptions ){/* -------------------------------------------------------------------- *//*      Verify input options.                                           *//* -------------------------------------------------------------------- */    if( eType != GDT_Byte && eType != GDT_UInt16 )    {        CPLError( CE_Failure, CPLE_AppDefined,              "Attempt to create PNM dataset with an illegal/n"              "data type (%s), only Byte and UInt16 supported./n",              GDALGetDataTypeName(eType) );        return NULL;    }    if( nBands != 1 && nBands != 3 )    {        CPLError( CE_Failure, CPLE_AppDefined,                  "Attempt to create PNM dataset with an illegal number/n"                  "of bands (%d).  Must be 1 (greyscale) or 3 (RGB)./n",                  nBands );        return NULL;    }/* -------------------------------------------------------------------- *//*      Try to create the file.                                         *//* -------------------------------------------------------------------- */    VSILFILE        *fp;    fp = VSIFOpenL( pszFilename, "wb" );    if( fp == NULL )    {        CPLError( CE_Failure, CPLE_OpenFailed,                  "Attempt to create file `%s' failed./n",                  pszFilename );        return NULL;    }/* -------------------------------------------------------------------- *//*      Write out the header.                                           *//* -------------------------------------------------------------------- */    char        szHeader[500];    const char  *pszMaxValue = NULL;    int         nMaxValue = 0;    pszMaxValue = CSLFetchNameValue( papszOptions, "MAXVAL" );    if ( pszMaxValue )    {        nMaxValue = atoi( pszMaxValue );        if ( eType == GDT_Byte && (nMaxValue > 255 || nMaxValue < 0) )            nMaxValue = 255;        else if ( nMaxValue > 65535 || nMaxValue < 0 )            nMaxValue = 65535;    }    else    {        if ( eType == GDT_Byte )            nMaxValue = 255;        else            nMaxValue = 65535;    }    memset( szHeader, 0, sizeof(szHeader) );    if( nBands == 3 )        sprintf( szHeader, "P6/n%d %d/n%d/n", nXSize, nYSize, nMaxValue );    else        sprintf( szHeader, "P5/n%d %d/n%d/n", nXSize, nYSize, nMaxValue );    VSIFWriteL( (void *) szHeader, strlen(szHeader) + 2, 1, fp );    VSIFCloseL( fp );    return (GDALDataset *) GDALOpen( pszFilename, GA_Update );}

int LLVMFuzzerTestOneInput(const uint8_t *buf, size_t len){    VSILFILE* fp = VSIFileFromMemBuffer( "/vsimem/test.tar",            reinterpret_cast<GByte*>(const_cast<uint8_t*>(buf)), len, FALSE );    VSIFCloseL(fp);    CPLPushErrorHandler(CPLQuietErrorHandler);    char** papszArgv = nullptr;    // Prevent generating too big output raster. Make sure they are set at    // the beginning to avoid being accidentally eaten by invalid arguments    // afterwards.    papszArgv = CSLAddString(papszArgv, "-limit_outsize");    papszArgv = CSLAddString(papszArgv, "1000000");    fp = VSIFOpenL("/vsitar//vsimem/test.tar/cmd.txt", "rb");    if( fp != nullptr )    {        const char* pszLine = nullptr;        while( (pszLine = CPLReadLineL(fp)) != nullptr )        {            if( !EQUAL(pszLine, "-limit_outsize") )                papszArgv = CSLAddString(papszArgv, pszLine);        }        VSIFCloseL(fp);    }    int nXDim = -1;    int nYDim = -1;    bool bXDimPct = false;    bool bYDimPct = false;    bool bNonNearestResampling = false;    int nBlockXSize = 0;    int nBlockYSize = 0;    bool bStatsEnabled = false;    bool bHFA = false;    if( papszArgv != nullptr )    {        int nCount = CSLCount(papszArgv);        for( int i = 0; i < nCount; i++ )        {            if( EQUAL(papszArgv[i], "-outsize") && i + 2 < nCount )            {                nXDim = atoi(papszArgv[i+1]);                bXDimPct = (papszArgv[i+1][0] != '/0' &&                            papszArgv[i+1][strlen(papszArgv[i+1])-1] == '%');                nYDim = atoi(papszArgv[i+2]);                bYDimPct = (papszArgv[i+2][0] != '/0' &&                            papszArgv[i+2][strlen(papszArgv[i+2])-1] == '%');            }            else if( EQUAL(papszArgv[i], "-r") && i + 1 < nCount )            {                bNonNearestResampling = !STARTS_WITH_CI(papszArgv[i+1], "NEAR");            }            else if( EQUAL(papszArgv[i], "-co") && i + 1 < nCount )            {                if( STARTS_WITH_CI(papszArgv[i+1], "BLOCKSIZE=") )                {                    nBlockXSize = std::max(nBlockXSize,                                atoi(papszArgv[i+1]+strlen("BLOCKSIZE=")));                    nBlockYSize = std::max(nBlockYSize,                                atoi(papszArgv[i+1]+strlen("BLOCKSIZE=")));                }                else if( STARTS_WITH_CI(papszArgv[i+1], "BLOCKXSIZE=") )                {                    nBlockXSize = std::max(nBlockXSize,                                atoi(papszArgv[i+1]+strlen("BLOCKXSIZE=")));                }                else if( STARTS_WITH_CI(papszArgv[i+1], "BLOCKYSIZE=") )                {                    nBlockYSize = std::max(nBlockYSize,                                atoi(papszArgv[i+1]+strlen("BLOCKYSIZE=")));                }            }            else if( EQUAL(papszArgv[i], "-stats") )            {                bStatsEnabled = true;            }            else if( EQUAL(papszArgv[i], "-of") && i + 1 < nCount )            {                bHFA = EQUAL( papszArgv[i+1], "HFA" );            }        }        if( bHFA )        {            // Disable statistics computation for HFA, as it can be time            // consuming.            // See https://bugs.chromium.org/p/oss-fuzz/issues/detail?id=10067            papszArgv = CSLInsertString(papszArgv, 0, "-co");            papszArgv = CSLInsertString(papszArgv, 1, "STATISTICS=NO");        }    }    if( papszArgv != nullptr )    {        GDALTranslateOptions* psOptions = GDALTranslateOptionsNew(papszArgv, nullptr);        if( psOptions )        {            GDALDatasetH hSrcDS = GDALOpen( "/vsitar//vsimem/test.tar/in", GA_ReadOnly );//.........这里部分代码省略.........

//.........这里部分代码省略.........        else if( EQUAL(argv[i], "-nomd") )            bShowMetadata = FALSE;        else if( EQUAL(argv[i], "-norat") )            bShowRAT = FALSE;        else if( EQUAL(argv[i], "-noct") )            bShowColorTable = FALSE;        else if( EQUAL(argv[i], "-mdd") )        {            CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);            papszExtraMDDomains = CSLAddString( papszExtraMDDomains,                                                argv[++i] );        }        else if( EQUAL(argv[i], "-nofl") )            bShowFileList = FALSE;        else if( EQUAL(argv[i], "-sd") )        {            CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);            nSubdataset = atoi(argv[++i]);        }        else if( argv[i][0] == '-' )            Usage(CPLSPrintf("Unkown option name '%s'", argv[i]));        else if( pszFilename == NULL )            pszFilename = argv[i];        else            Usage("Too many command options.");    }    if( pszFilename == NULL )        Usage("No datasource specified.");/* -------------------------------------------------------------------- *//*      Open dataset.                                                   *//* -------------------------------------------------------------------- */    hDataset = GDALOpen( pszFilename, GA_ReadOnly );    if( hDataset == NULL )    {        fprintf( stderr,                 "gdalinfo failed - unable to open '%s'./n",                 pszFilename );/* -------------------------------------------------------------------- *//*      If argument is a VSIFILE, then print its contents               *//* -------------------------------------------------------------------- */        if ( strncmp( pszFilename, "/vsizip/", 8 ) == 0 ||              strncmp( pszFilename, "/vsitar/", 8 ) == 0 )         {            papszFileList = VSIReadDirRecursive( pszFilename );            if ( papszFileList )            {                int nCount = CSLCount( papszFileList );                fprintf( stdout,                          "Unable to open source `%s' directly./n"                         "The archive contains %d files:/n",                          pszFilename, nCount );                for ( i = 0; i < nCount; i++ )                {                    fprintf( stdout, "       %s/%s/n", pszFilename, papszFileList[i] );                }                CSLDestroy( papszFileList );                papszFileList = NULL;            }        }        CSLDestroy( argv );        CSLDestroy( papszExtraMDDomains );

int equalize_density(char *infile, char *outfile, int fast, int accurate) {		int xsize, ysize;				// Size of the density grid.	double *gridx, *gridy;			// Array for grid		double **rho;					// Initial population density	GDALDatasetH hDataset;			// The input density raster file.	GDALRasterBandH hBand;			// The raster band we are going to use.	FILE *outfp;					// The morphing file (a text file).	double adfGeoTransform[6];		// For the georeference of the raster.			// Register all GDAL drivers.    GDALAllRegister();		#if defined (_OPENMP)	omp_set_num_threads(omp_get_num_procs());#endif			hDataset = GDALOpen(infile, GA_ReadOnly);    if (hDataset == NULL) {		fprintf(stderr,"Error. Unable to open file `%s'/n", infile);		exit(1);	}		outfp = fopen(outfile, "w");	if (outfp == NULL) {		fprintf(stderr,"Error. Unable to open file `%s'/n", outfile);		exit(1);	}			// Get the raster band for the dataset; we are using the first band.	hBand = GDALGetRasterBand(hDataset, 1);	if (hBand == NULL) {		fprintf(stderr, "Error. Unable to read band 1 in file `%s'/n", infile);		exit(1);	}		// Determine the raster size	xsize = GDALGetRasterBandXSize(hBand);	ysize = GDALGetRasterBandYSize(hBand);				// Allocate space for the cartogram code to use	cart_makews(xsize, ysize);			// Read in the population data, transform it, then destroy it again	rho = cart_dmalloc(xsize, ysize);	if (readpop(hBand, rho, xsize, ysize)) {		fprintf(stderr,"Error. Density file contains too few or incorrect data/n");		exit(1);	}	cart_transform(rho, xsize, ysize);	cart_dfree(rho);			// Create the grid of points	gridx = malloc((xsize+1)*(ysize+1)*sizeof(double));	gridy = malloc((xsize+1)*(ysize+1)*sizeof(double));	creategrid(gridx, gridy, xsize, ysize);			// Compute the cartogram	cart_makecart(gridx, gridy, (xsize+1)*(ysize+1), xsize, ysize, 0.0);			// Write out the final positions of the grid points	GDALGetGeoTransform(hDataset, adfGeoTransform);	writepoints(outfp, gridx, gridy, xsize, ysize, adfGeoTransform);	//writepoints(outfp, gridx, gridy, (xsize+1)*(ysize+1));			// Free up the allocated memory	cart_freews(xsize, ysize);	free(gridx);	free(gridy);			// Close the input and output files	GDALClose(hDataset);	fclose(outfp);		return 0;}

//.........这里部分代码省略........./* -------------------------------------------------------------------- */    VSILFILE *fp = VSIFOpenL( pszFilename, "w" );    if( fp == NULL )    {        CPLError( CE_Failure, CPLE_OpenFailed,                  "Failed to create file '%s': %s",                  pszFilename, VSIStrerror( errno ) );        return NULL;    }/* -------------------------------------------------------------------- *//*      Prepare and write 512 byte header.                              *//* -------------------------------------------------------------------- */    GByte abyHeader[512];    memset( abyHeader, 0, 512 );    abyHeader[0] = 1;    abyHeader[1] = 218;    abyHeader[2] = 1; // RLE    abyHeader[3] = 1;  // 8bit    GInt16 nShortValue;    if( nBands == 1 )        nShortValue = CPL_MSBWORD16(2);    else        nShortValue = CPL_MSBWORD16(3);    memcpy( abyHeader + 4, &nShortValue, 2 );    nShortValue = CPL_MSBWORD16(nXSize);    memcpy( abyHeader + 6, &nShortValue, 2 );    nShortValue = CPL_MSBWORD16(nYSize);    memcpy( abyHeader + 8, &nShortValue, 2 );    nShortValue = CPL_MSBWORD16(nBands);    memcpy( abyHeader + 10, &nShortValue, 2 );    GInt32 nIntValue = CPL_MSBWORD32(0);    memcpy( abyHeader + 12, &nIntValue, 4 );    GUInt32 nUIntValue = CPL_MSBWORD32(255);    memcpy( abyHeader + 16, &nUIntValue, 4 );    VSIFWriteL( abyHeader, 1, 512, fp );/* -------------------------------------------------------------------- *//*      Create our RLE compressed zero-ed dummy line.                   *//* -------------------------------------------------------------------- */    GByte *pabyRLELine = reinterpret_cast<GByte *>(        CPLMalloc( ( nXSize / 127 ) * 2 + 4 ) );    int nPixelsRemaining = nXSize;    GInt32 nRLEBytes = 0;    while( nPixelsRemaining > 0 )    {        pabyRLELine[nRLEBytes] = static_cast<GByte>(            std::min( 127, nPixelsRemaining  ) );        pabyRLELine[nRLEBytes+1] = 0;        nPixelsRemaining -= pabyRLELine[nRLEBytes];        nRLEBytes += 2;    }/* -------------------------------------------------------------------- *//*      Prepare and write RLE offset/size tables with everything        *//*      zeroed indicating dummy lines.                                  *//* -------------------------------------------------------------------- */    const int nTableLen = nYSize * nBands;    GInt32 nDummyRLEOffset = 512 + 4 * nTableLen * 2;    CPL_MSBPTR32( &nRLEBytes );    CPL_MSBPTR32( &nDummyRLEOffset );    for( int i = 0; i < nTableLen; i++ )        VSIFWriteL( &nDummyRLEOffset, 1, 4, fp );    for( int i = 0; i < nTableLen; i++ )        VSIFWriteL( &nRLEBytes, 1, 4, fp );/* -------------------------------------------------------------------- *//*      write the dummy RLE blank line.                                 *//* -------------------------------------------------------------------- */    CPL_MSBPTR32( &nRLEBytes );    if( static_cast<GInt32>( VSIFWriteL( pabyRLELine, 1, nRLEBytes, fp ) )        != nRLEBytes )    {        CPLError( CE_Failure, CPLE_FileIO,                  "Failure writing SGI file '%s'./n%s",                  pszFilename,                  VSIStrerror( errno ) );        return NULL;    }    VSIFCloseL( fp );    CPLFree( pabyRLELine );    return reinterpret_cast<GDALDataset *>(        GDALOpen( pszFilename, GA_Update ) );}

GDALDatasetH QgsGdalProviderBase::gdalOpen( const char *pszFilename, GDALAccess eAccess ){  GDALDatasetH hDS = GDALOpen( pszFilename, eAccess );  return hDS;}

GDALDataset *GTXDataset::Create( const char * pszFilename,                                 int nXSize, int nYSize, int nBands,                                 GDALDataType eType,                                 char ** papszOptions ){    if( eType != GDT_Float32 )    {        CPLError( CE_Failure, CPLE_AppDefined,                  "Attempt to create gtx file with unsupported data type '%s'.",                  GDALGetDataTypeName( eType ) );        return NULL;    }    if( !EQUAL(CPLGetExtension(pszFilename),"gtx") )    {        CPLError( CE_Failure, CPLE_AppDefined,                  "Attempt to create gtx file with extension other than gtx." );        return NULL;    }/* -------------------------------------------------------------------- *//*      Try to create the file.                                         *//* -------------------------------------------------------------------- */    VSILFILE	*fp;    fp = VSIFOpenL( pszFilename, "wb" );    if( fp == NULL )    {        CPLError( CE_Failure, CPLE_OpenFailed,                  "Attempt to create file `%s' failed./n",                  pszFilename );        return NULL;    }/* -------------------------------------------------------------------- *//*      Write out the header with stub georeferencing.                  *//* -------------------------------------------------------------------- */    unsigned char header[40];    double dfXOrigin=0, dfYOrigin=0, dfXSize=0.01, dfYSize=0.01;    GInt32 nXSize32 = nXSize, nYSize32 = nYSize;    memcpy( header + 0, &dfYOrigin, 8 );    CPL_MSBPTR64( header + 0 );    memcpy( header + 8, &dfXOrigin, 8 );    CPL_MSBPTR64( header + 8 );    memcpy( header + 16, &dfYSize, 8 );    CPL_MSBPTR64( header + 16 );    memcpy( header + 24, &dfXSize, 8 );    CPL_MSBPTR64( header + 24 );    memcpy( header + 32, &nYSize32, 4 );    CPL_MSBPTR32( header + 32 );    memcpy( header + 36, &nXSize32, 4 );    CPL_MSBPTR32( header + 36 );    VSIFWriteL( header, 40, 1, fp );    VSIFCloseL( fp );    return (GDALDataset *) GDALOpen( pszFilename, GA_Update );}

//.........这里部分代码省略........./* -------------------------------------------------------------------- *//*      Open the file.                                                  *//* -------------------------------------------------------------------- */    fp = VSIFOpenL( pszAuxFilename, "wt" );    if( fp == NULL )    {        CPLError( CE_Failure, CPLE_OpenFailed,                  "Attempt to create file `%s' failed./n",                  pszAuxFilename );        return NULL;    }    CPLFree( pszAuxFilename );    /* -------------------------------------------------------------------- *//*      We need to write out the original filename but without any      *//*      path components in the AuxilaryTarget line.  Do so now.         *//* -------------------------------------------------------------------- */    int		iStart;    iStart = strlen(pszFilename)-1;    while( iStart > 0 && pszFilename[iStart-1] != '/'           && pszFilename[iStart-1] != '//' )        iStart--;    VSIFPrintfL( fp, "AuxilaryTarget: %s/n", pszFilename + iStart );/* -------------------------------------------------------------------- *//*      Write out the raw definition for the dataset as a whole.        *//* -------------------------------------------------------------------- */    VSIFPrintfL( fp, "RawDefinition: %d %d %d/n",                nXSize, nYSize, nBands );/* -------------------------------------------------------------------- *//*      Write out a definition for each band.  We always write band     *//*      sequential files for now as these are pretty efficiently        *//*      handled by GDAL.                                                *//* -------------------------------------------------------------------- */    vsi_l_offset    nImgOffset = 0;        for( iBand = 0; iBand < nBands; iBand++ )    {        const char  *pszTypeName;        int         nPixelOffset, nLineOffset;		vsi_l_offset nNextImgOffset;/* -------------------------------------------------------------------- *//*      Establish our file layout based on supplied interleaving.       *//* -------------------------------------------------------------------- */		if( EQUAL(pszInterleave,"LINE") )		{			nPixelOffset = GDALGetDataTypeSize(eType)/8;			nLineOffset = nXSize * nPixelSizeSum;			nNextImgOffset = nImgOffset + nPixelOffset * nXSize;		}		else if( EQUAL(pszInterleave,"PIXEL") )		{			nPixelOffset = nPixelSizeSum;			nLineOffset = nXSize * nPixelOffset;			nNextImgOffset = nImgOffset + (GDALGetDataTypeSize(eType)/8);		}		else /* default to band */		{			nPixelOffset = GDALGetDataTypeSize(eType)/8;			nLineOffset = nXSize * nPixelOffset;			nNextImgOffset = nImgOffset + nYSize * (vsi_l_offset) nLineOffset;		}/* -------------------------------------------------------------------- *//*      Write out line indicating layout.                               *//* -------------------------------------------------------------------- */        if( eType == GDT_Float32 )            pszTypeName = "32R";        else if( eType == GDT_Int16 )            pszTypeName = "16S";        else if( eType == GDT_UInt16 )            pszTypeName = "16U";        else            pszTypeName = "8U";        VSIFPrintfL( fp, "ChanDefinition-%d: %s " CPL_FRMT_GIB " %d %d %s/n", 					 iBand+1,					 pszTypeName, (GIntBig) nImgOffset,					 nPixelOffset, nLineOffset,#ifdef CPL_LSB                    "Swapped"#else                    "Unswapped"#endif                    );        nImgOffset = nNextImgOffset;    }/* -------------------------------------------------------------------- *//*      Cleanup                                                         *//* -------------------------------------------------------------------- */    VSIFCloseL( fp );    return (GDALDataset *) GDALOpen( pszFilename, GA_Update );}

GDALDataset *HDF5Dataset::Open( GDALOpenInfo * poOpenInfo ){    HDF5Dataset *poDS;    if( !Identify( poOpenInfo ) )        return NULL;/* -------------------------------------------------------------------- *//*      Create datasource.                                              *//* -------------------------------------------------------------------- */    poDS = new HDF5Dataset();    poDS->SetDescription( poOpenInfo->pszFilename );/* -------------------------------------------------------------------- *//*      Try opening the dataset.                                        *//* -------------------------------------------------------------------- */    poDS->hHDF5 = H5Fopen( poOpenInfo->pszFilename,                           H5F_ACC_RDONLY,                           H5P_DEFAULT );    if( poDS->hHDF5 < 0 )  {        delete poDS;        return NULL;    }    poDS->hGroupID = H5Gopen( poDS->hHDF5, "/" );    if( poDS->hGroupID < 0 ) {        poDS->bIsHDFEOS=false;        delete poDS;        return NULL;    }    poDS->bIsHDFEOS=true;    poDS->ReadGlobalAttributes( true );    poDS->SetMetadata( poDS->papszMetadata  );    if ( CSLCount( poDS->papszSubDatasets ) / 2 >= 1 )        poDS->SetMetadata( poDS->papszSubDatasets, "SUBDATASETS" );    // Make sure we don't try to do any pam stuff with this dataset.    poDS->nPamFlags |= GPF_NOSAVE;/* -------------------------------------------------------------------- *//*      If we have single subdataset only, open it immediately          *//* -------------------------------------------------------------------- */    int nSubDatasets = CSLCount( poDS->papszSubDatasets ) / 2;    if( nSubDatasets == 1 )    {        CPLString osDSName = CSLFetchNameValue( poDS->papszSubDatasets,                                                "SUBDATASET_1_NAME" );        delete poDS;        return (GDALDataset *) GDALOpen( osDSName, poOpenInfo->eAccess );    }    else    {/* -------------------------------------------------------------------- *//*      Confirm the requested access is supported.                      *//* -------------------------------------------------------------------- */        if( poOpenInfo->eAccess == GA_Update )        {            delete poDS;            CPLError( CE_Failure, CPLE_NotSupported,                      "The HDF5 driver does not support update access to existing"                      " datasets./n" );            return NULL;        }    }    return( poDS );}

//.........这里部分代码省略.........    int nDecimalCount = 7;    int bHasNoDataValue = FALSE;    double dfNoDataValue =        poSrcDS->GetRasterBand(1)->GetNoDataValue(&bHasNoDataValue);    if (!bHasNoDataValue)        dfNoDataValue = 1.e30;    VSIFPrintfL(fp, "!/n");    VSIFPrintfL(fp, "! Created by GDAL./n");    VSIFPrintfL(fp, "!/n");    VSIFPrintfL(fp, "@GRID FILE, GRID, %d/n", nValuesPerLine);    WriteRightJustified(fp, nFieldSize, 10);    VSIFPrintfL(fp, ",");    WriteRightJustified(fp, dfNoDataValue, 10);    VSIFPrintfL(fp, ",");    WriteRightJustified(fp, "", 10);    VSIFPrintfL(fp, ",");    WriteRightJustified(fp, nDecimalCount, 10);    VSIFPrintfL(fp, ",");    WriteRightJustified(fp, 1, 10);    VSIFPrintfL(fp, "/n");    WriteRightJustified(fp, nYSize, 10);    VSIFPrintfL(fp, ",");    WriteRightJustified(fp, nXSize, 10);    VSIFPrintfL(fp, ",");    if (CSLTestBoolean(CPLGetConfigOption("ZMAP_PIXEL_IS_POINT", "FALSE")))    {        WriteRightJustified(fp, adfGeoTransform[0] + adfGeoTransform[1] / 2, 14, 7);        VSIFPrintfL(fp, ",");        WriteRightJustified(fp, adfGeoTransform[0] + adfGeoTransform[1] * nXSize -                                adfGeoTransform[1] / 2, 14, 7);        VSIFPrintfL(fp, ",");        WriteRightJustified(fp, adfGeoTransform[3] + adfGeoTransform[5] * nYSize -                                adfGeoTransform[5] / 2, 14, 7);        VSIFPrintfL(fp, ",");        WriteRightJustified(fp, adfGeoTransform[3] + adfGeoTransform[5] / 2, 14, 7);    }    else    {        WriteRightJustified(fp, adfGeoTransform[0], 14, 7);        VSIFPrintfL(fp, ",");        WriteRightJustified(fp, adfGeoTransform[0] + adfGeoTransform[1] * nXSize, 14, 7);        VSIFPrintfL(fp, ",");        WriteRightJustified(fp, adfGeoTransform[3] + adfGeoTransform[5] * nYSize, 14, 7);        VSIFPrintfL(fp, ",");        WriteRightJustified(fp, adfGeoTransform[3], 14, 7);    }    VSIFPrintfL(fp, "/n");    VSIFPrintfL(fp, "0.0, 0.0, 0.0/n");    VSIFPrintfL(fp, "@/n");/* -------------------------------------------------------------------- *//*      Copy imagery                                                    *//* -------------------------------------------------------------------- */    double* padfLineBuffer = (double*) CPLMalloc(nYSize * sizeof(double));    int i, j;    CPLErr eErr = CE_None;    for(i=0;i<nXSize && eErr == CE_None;i++)    {        eErr = poSrcDS->GetRasterBand(1)->RasterIO(                                            GF_Read, i, 0, 1, nYSize,                                            padfLineBuffer, 1, nYSize,                                            GDT_Float64, 0, 0);        if (eErr != CE_None)            break;        int bEOLPrinted = FALSE;        for(j=0;j<nYSize;j++)        {            WriteRightJustified(fp, padfLineBuffer[j], nFieldSize, nDecimalCount);            if (((j + 1) % nValuesPerLine) == 0)            {                bEOLPrinted = TRUE;                VSIFPrintfL(fp, "/n");            }            else                bEOLPrinted = FALSE;        }        if (!bEOLPrinted)            VSIFPrintfL(fp, "/n");        if (!pfnProgress( (j+1) * 1.0 / nYSize, NULL, pProgressData))        {            eErr = CE_Failure;            break;        }    }    CPLFree(padfLineBuffer);    VSIFCloseL(fp);    if (eErr != CE_None)        return NULL;    return (GDALDataset*) GDALOpen(pszFilename, GA_ReadOnly);}

//.........这里部分代码省略.........            if (bMustMemset)                memset(pabyBuffer, 0, nBands * nBlockXSize * nBlockYSize);                        eErr = poSrcDS->RasterIO(GF_Read,                              nBlockXOff * nBlockXSize,                              nBlockYOff * nBlockYSize,                              nReqXSize, nReqYSize,                              pabyBuffer,                              nReqXSize, nReqYSize,                              GDT_Byte, nBands, NULL,                              1,                              nBlockXSize,                              nBlockXSize * nBlockYSize);                        int nOutBufSize = nTargetBlockSize;            if (eErr == CE_None && nBands == 1)            {                if (EPS_OK != eps_encode_grayscale_block(apapbyRawBuffer[0],                                           nXSize, nYSize,                                           nReqXSize, nReqYSize,                                           nBlockXOff * nBlockXSize,                                           nBlockYOff * nBlockYSize,                                           pabyOutBuf, &nOutBufSize,                                           (char*) pszFilter, eMode))                {                    CPLError(CE_Failure, CPLE_AppDefined,                             "Error occured when encoding block (%d, %d)",                             nBlockXOff, nBlockYOff);                    eErr = CE_Failure;                }            }            else if (eErr == CE_None)            {                if (EPS_OK != eps_encode_truecolor_block(                                           apapbyRawBuffer[0],                                           apapbyRawBuffer[1],                                           apapbyRawBuffer[2],                                           nXSize, nYSize,                                           nReqXSize, nReqYSize,                                           nBlockXOff * nBlockXSize,                                           nBlockYOff * nBlockYSize,                                           eResample,                                           pabyOutBuf, &nOutBufSize,                                           nYRatio, nCbRatio, nCrRatio,                                           (char*) pszFilter, eMode))                {                    CPLError(CE_Failure, CPLE_AppDefined,                             "Error occured when encoding block (%d, %d)",                             nBlockXOff, nBlockYOff);                    eErr = CE_Failure;                }            }                        if (eErr == CE_None)            {                if ((int)VSIFWriteL(pabyOutBuf, 1, nOutBufSize, fp) !=                                                                nOutBufSize)                    eErr = CE_Failure;                char chEPSMarker = EPS_MARKER;                VSIFWriteL(&chEPSMarker, 1, 1, fp);                                if (pfnProgress && !pfnProgress(                      1.0 * (nBlockYOff * nXBlocks + nBlockXOff + 1) / nBlocks,                      NULL, pProgressData))                {                    eErr = CE_Failure;                }            }        }    }        if (bRasterliteOutput)    {        const char* pszFooter = RASTERLITE_WAVELET_FOOTER;        VSIFWriteL(pszFooter, 1, strlen(pszFooter) + 1, fp);    }/* -------------------------------------------------------------------- *//*      Cleanup work buffers                                            *//* -------------------------------------------------------------------- */        for(i=0;i<nBands;i++)    {        VSIFree(apapbyRawBuffer[i]);    }        VSIFree(pabyOutBuf);    VSIFree(pabyBuffer);            VSIFCloseL(fp);        if (eErr != CE_None)        return NULL;/* -------------------------------------------------------------------- *//*      Reopen the dataset, unless asked for not (Rasterlite optim)     *//* -------------------------------------------------------------------- */    return (GDALDataset*) GDALOpen(pszFilename, GA_ReadOnly);}

int QgsZonalStatistics::calculateStatistics( QProgressDialog* p ){  if ( !mPolygonLayer || mPolygonLayer->geometryType() != QGis::Polygon )  {    return 1;  }  QgsVectorDataProvider* vectorProvider = mPolygonLayer->dataProvider();  if ( !vectorProvider )  {    return 2;  }  //open the raster layer and the raster band  GDALAllRegister();  GDALDatasetH inputDataset = GDALOpen( TO8F( mRasterFilePath ), GA_ReadOnly );  if ( inputDataset == NULL )  {    return 3;  }  if ( GDALGetRasterCount( inputDataset ) < ( mRasterBand - 1 ) )  {    GDALClose( inputDataset );    return 4;  }  GDALRasterBandH rasterBand = GDALGetRasterBand( inputDataset, mRasterBand );  if ( rasterBand == NULL )  {    GDALClose( inputDataset );    return 5;  }  mInputNodataValue = GDALGetRasterNoDataValue( rasterBand, NULL );  //get geometry info about raster layer  int nCellsXGDAL = GDALGetRasterXSize( inputDataset );  int nCellsYGDAL = GDALGetRasterYSize( inputDataset );  double geoTransform[6];  if ( GDALGetGeoTransform( inputDataset, geoTransform ) != CE_None )  {    GDALClose( inputDataset );    return 6;  }  double cellsizeX = geoTransform[1];  if ( cellsizeX < 0 )  {    cellsizeX = -cellsizeX;  }  double cellsizeY = geoTransform[5];  if ( cellsizeY < 0 )  {    cellsizeY = -cellsizeY;  }  QgsRectangle rasterBBox( geoTransform[0], geoTransform[3] - ( nCellsYGDAL * cellsizeY ),                           geoTransform[0] + ( nCellsXGDAL * cellsizeX ), geoTransform[3] );  //add the new count, sum, mean fields to the provider  QList<QgsField> newFieldList;  QString countFieldName = getUniqueFieldName( mAttributePrefix + "count" );  QString sumFieldName = getUniqueFieldName( mAttributePrefix + "sum" );  QString meanFieldName = getUniqueFieldName( mAttributePrefix + "mean" );  QgsField countField( countFieldName, QVariant::Double, "double precision" );  QgsField sumField( sumFieldName, QVariant::Double, "double precision" );  QgsField meanField( meanFieldName, QVariant::Double, "double precision" );  newFieldList.push_back( countField );  newFieldList.push_back( sumField );  newFieldList.push_back( meanField );  vectorProvider->addAttributes( newFieldList );  //index of the new fields  int countIndex = vectorProvider->fieldNameIndex( countFieldName );  int sumIndex = vectorProvider->fieldNameIndex( sumFieldName );  int meanIndex = vectorProvider->fieldNameIndex( meanFieldName );  if ( countIndex == -1 || sumIndex == -1 || meanIndex == -1 )  {    return 8;  }  //progress dialog  long featureCount = vectorProvider->featureCount();  if ( p )  {    p->setMaximum( featureCount );  }  //iterate over each polygon  QgsFeatureRequest request;  request.setSubsetOfAttributes( QgsAttributeList() );  QgsFeatureIterator fi = vectorProvider->getFeatures( request );  QgsFeature f;  double count = 0;  double sum = 0;  double mean = 0;  int featureCounter = 0;  while ( fi.nextFeature( f ) )  {//.........这里部分代码省略.........

int main(int argc, char ** argv) {	program_options options;	if (!options.parse_args(argc, argv)) return 1;	static const size_t memory = options.memory * 1024 * 1024;	/// Initialize GDAL	GDALAllRegister();	std::unique_ptr<GDALDataset> in((GDALDataset*)GDALOpen(options.input_file.c_str(), GA_ReadOnly));	int xsize = in->GetRasterXSize();	int ysize = in->GetRasterYSize();	if (options.outputxsize == -1) options.outputxsize=xsize;	if (options.outputysize == -1) options.outputysize=ysize;	GDALDriver * driver = GetGDALDriverManager()->GetDriverByName("ENVI");	std::unique_ptr<GDALDataset> out(driver->Create(options.output_file.c_str(),													options.outputxsize, options.outputysize,													1, GDT_Float32 , nullptr));	double geoCoords[6];  //Geographic metadata	if (in->GetGeoTransform(geoCoords) == CE_None) out->SetGeoTransform(geoCoords);	const char * sref = in->GetProjectionRef();	if (sref != NULL) out->SetProjection(sref);	GDALRasterBand * in_band = in->GetRasterBand(1);	GDALRasterBand * out_band = out->GetRasterBand(1);	float nodata;	int has_nodata = false;	nodata = in_band->GetNoDataValue(&has_nodata);	if (!has_nodata) nodata=-9999;	out_band->SetNoDataValue(nodata);	/// Done initializing GDAL	// std::vector memory usage is at most twice the size of the elements.	const size_t row_memory = sizeof(float) * xsize * 2;	// Only one sorter in phase 1.	const size_t phase1_sort_memory = memory;	// Memory per sorter for the two sorters in phase 2.	const size_t phase2_sort_memory = (memory - row_memory) / 2;	// Row from previous phase is not deallocated automatically,	// so we have two rows and one sorter in memory in phase 3.	const size_t phase3_sort_memory = memory - 2*row_memory;	point_generator pg(options.outputxsize, options.outputysize);	typedef point_mapper<point_generator, matrix_transform> mapped_pg_type;	mapped_pg_type mapped_pg(pg, options.transform, xsize, ysize);	typedef sxs::runs_creator<mapped_pg_type, map_point_value_less> rc_type;	rc_type rc(mapped_pg, map_point_value_less(), phase1_sort_memory);	typedef sxs::runs_merger<rc_type::sorted_runs_type, map_point_value_less> rm_type;	rm_type rm(rc.result(), map_point_value_less(), phase2_sort_memory);	std::cout << "Initialized runs merger 1" << std::endl;	typedef RasterReader raster_reader_type;	raster_reader_type raster_reader(in_band, xsize, ysize);	typedef Filler<rm_type, raster_reader_type> filler_type;	filler_type filler(rm, raster_reader, nodata);	typedef sxs::runs_creator<filler_type, value_point_less> rc_filler_type;	rc_filler_type rc_filler(filler, value_point_less(), phase2_sort_memory);	typedef sxs::runs_merger<rc_filler_type::sorted_runs_type, value_point_less> rm_filler_type;	rm_filler_type rm_filler(rc_filler.result(), value_point_less(), phase3_sort_memory);	std::cout << "Initialized runs merger 2" << std::endl;	typedef PointToRaster<rm_filler_type> point_to_raster_type;	point_to_raster_type point_to_raster(rm_filler, options.outputxsize, options.outputysize, nodata);	std::cout << "Initialized point_to_raster" << std::endl;	write_raster(point_to_raster, out_band, options.outputxsize);	std::cout << "I am now done" << std::endl;	return 0;}

int main( int argc, char ** argv ){    EarlySetConfigOptions(argc, argv);    GDALAllRegister();    argc = GDALGeneralCmdLineProcessor( argc, &argv, 0 );    if( argc < 1 )        exit( -argc );    for( int i = 0; argv != NULL && argv[i] != NULL; i++ )    {        if( EQUAL(argv[i], "--utility_version") )        {            printf("%s was compiled against GDAL %s and is running against GDAL %s/n",                   argv[0], GDAL_RELEASE_NAME, GDALVersionInfo("RELEASE_NAME"));            CSLDestroy( argv );            return 0;        }        else if( EQUAL(argv[i],"--help") )        {            Usage();        }    }    argv = CSLAddString(argv, "-stdout");    GDALInfoOptionsForBinary* psOptionsForBinary = GDALInfoOptionsForBinaryNew();    GDALInfoOptions *psOptions        = GDALInfoOptionsNew(argv + 1, psOptionsForBinary);    if( psOptions == NULL )        Usage();    if( psOptionsForBinary->pszFilename == NULL )        Usage("No datasource specified.");/* -------------------------------------------------------------------- *//*      Open dataset.                                                   *//* -------------------------------------------------------------------- */#ifdef __AFL_HAVE_MANUAL_CONTROL    int iIter = 0;    while (__AFL_LOOP(1000)) {        iIter ++;#endif    GDALDatasetH hDataset        = GDALOpenEx( psOptionsForBinary->pszFilename, GDAL_OF_READONLY | GDAL_OF_RASTER | GDAL_OF_VERBOSE_ERROR, NULL,                      (const char* const* )psOptionsForBinary->papszOpenOptions, NULL );    if( hDataset == NULL )    {#ifdef __AFL_HAVE_MANUAL_CONTROL        continue;#else        fprintf( stderr,                 "gdalinfo failed - unable to open '%s'./n",                 psOptionsForBinary->pszFilename );/* -------------------------------------------------------------------- *//*      If argument is a VSIFILE, then print its contents               *//* -------------------------------------------------------------------- */        if ( STARTS_WITH(psOptionsForBinary->pszFilename, "/vsizip/") ||             STARTS_WITH(psOptionsForBinary->pszFilename, "/vsitar/") )        {            char** papszFileList = VSIReadDirRecursive( psOptionsForBinary->pszFilename );            if ( papszFileList )            {                int nCount = CSLCount( papszFileList );                fprintf( stdout,                         "Unable to open source `%s' directly./n"                         "The archive contains %d files:/n",                         psOptionsForBinary->pszFilename, nCount );                for ( int i = 0; i < nCount; i++ )                {                    fprintf( stdout, "       %s/%s/n", psOptionsForBinary->pszFilename, papszFileList[i] );                }                CSLDestroy( papszFileList );            }        }        CSLDestroy( argv );        GDALInfoOptionsForBinaryFree(psOptionsForBinary);        GDALInfoOptionsFree( psOptions );        GDALDumpOpenDatasets( stderr );        GDALDestroyDriverManager();        CPLDumpSharedList( NULL );        exit( 1 );#endif    }/* -------------------------------------------------------------------- *//*      Read specified subdataset if requested.                         *//* -------------------------------------------------------------------- *///.........这里部分代码省略.........

GDALDataset *ISIS2Dataset::Create(const char* pszFilename,                                  int nXSize, int nYSize, int nBands,                                  GDALDataType eType, char** papszParmList) {    /* Verify settings. In Isis 2 core pixel values can be represented in     * three different ways : 1, 2 4, or 8 Bytes */    if( eType != GDT_Byte && eType != GDT_Int16 && eType != GDT_Float32        && eType != GDT_UInt16 && eType != GDT_Float64 ){        CPLError(CE_Failure, CPLE_AppDefined,                 "The ISIS2 driver does not supporting creating files of type %s.",                 GDALGetDataTypeName( eType ) );        return nullptr;    }    /*  (SAMPLE, LINE, BAND) - Band Sequential (BSQ) - default choice        (SAMPLE, BAND, LINE) - Band Interleaved by Line (BIL)        (BAND, SAMPLE, LINE) - Band Interleaved by Pixel (BIP) */    const char *pszInterleaving = "(SAMPLE,LINE,BAND)";    const char *pszInterleavingParam = CSLFetchNameValue( papszParmList, "INTERLEAVE" );    if ( pszInterleavingParam ) {        if ( STARTS_WITH_CI(pszInterleavingParam, "bip") )            pszInterleaving = "(BAND,SAMPLE,LINE)";        else if ( STARTS_WITH_CI(pszInterleavingParam, "bil") )            pszInterleaving = "(SAMPLE,BAND,LINE)";        else            pszInterleaving = "(SAMPLE,LINE,BAND)";    }    /* default labeling method is attached */    bool bAttachedLabelingMethod = true;    /* check if labeling method is set : check the all three first chars */    const char *pszLabelingMethod = CSLFetchNameValue( papszParmList, "LABELING_METHOD" );    if ( pszLabelingMethod ){        if ( STARTS_WITH_CI( pszLabelingMethod, "det" /* "detached" */ ) ){            bAttachedLabelingMethod = false;        }        if ( STARTS_WITH_CI( pszLabelingMethod, "att" /* attached" */ ) ){            bAttachedLabelingMethod = true;        }    }    /*  set the label and data files */    CPLString osLabelFile, osRasterFile, osOutFile;    if( bAttachedLabelingMethod ) {        osLabelFile = "";        osRasterFile = pszFilename;        osOutFile = osRasterFile;    }    else    {        CPLString sExtension = "cub";        const char* pszExtension = CSLFetchNameValue( papszParmList, "IMAGE_EXTENSION" );        if( pszExtension ){            sExtension = pszExtension;        }        if( EQUAL(CPLGetExtension( pszFilename ), sExtension) )        {            CPLError( CE_Failure, CPLE_AppDefined,                      "IMAGE_EXTENSION (%s) cannot match LABEL file extension.",                      sExtension.c_str() );            return nullptr;        }        osLabelFile = pszFilename;        osRasterFile = CPLResetExtension( osLabelFile, sExtension );        osOutFile = osLabelFile;    }    const char *pszObject = CSLFetchNameValue( papszParmList, "OBJECT" );    CPLString sObject = "QUBE"; // default choice    if (pszObject) {        if ( EQUAL( pszObject, "IMAGE") ){            sObject = "IMAGE";        }        if ( EQUAL( pszObject, "SPECTRAL_QUBE")){            sObject = "SPECTRAL_QUBE";        }    }    GUIntBig iRecords = ISIS2Dataset::RecordSizeCalculation(nXSize, nYSize, nBands, eType);    GUIntBig iLabelRecords(2);    CPLDebug("ISIS2","irecord = %i",static_cast<int>(iRecords));    if( bAttachedLabelingMethod )    {        ISIS2Dataset::WriteLabel(osRasterFile, "", sObject, nXSize, nYSize, nBands, eType, iRecords, pszInterleaving, iLabelRecords, true);    }    else    {        ISIS2Dataset::WriteLabel(osLabelFile, osRasterFile, sObject, nXSize, nYSize, nBands, eType, iRecords, pszInterleaving, iLabelRecords);    }    if( !ISIS2Dataset::WriteRaster(osRasterFile, bAttachedLabelingMethod,                                   iRecords, iLabelRecords, eType,                                   pszInterleaving) )        return nullptr;    return reinterpret_cast<GDALDataset *>( GDALOpen( osOutFile, GA_Update ) );//.........这里部分代码省略.........

//.........这里部分代码省略.........    src_yoff = MAX(0,(int) floor(ury+0.5));    src_xsize = MIN(MAX(0,(int) (urx - llx + 0.5)),                    GDALGetRasterXSize(clinfo->hOrigDS) - src_xoff);    src_ysize = MIN(MAX(0,(int) (lly - ury + 0.5)),                    GDALGetRasterYSize(clinfo->hOrigDS) - src_yoff);    /* Update the geographic extent (buffer added) */    /* TODO: a better way to go the geo_trans */    copyRect.minx = GEO_TRANS(adfGeoTransform+0,src_xoff,0);    copyRect.maxx = GEO_TRANS(adfGeoTransform+0,src_xoff+src_xsize,0);    copyRect.miny = GEO_TRANS(adfGeoTransform+3,0,src_yoff+src_ysize);    copyRect.maxy = GEO_TRANS(adfGeoTransform+3,0,src_yoff);        /*      * If input window is to small then stop here     */    if (src_xsize < 2 || src_ysize < 2)    {      if (layer->debug)        msDebug("msContourLayerReadRaster(): input window too small, or no apparent overlap between map view and this window(1)./n");      return MS_SUCCESS;    }    /* Target buffer size */    dst_xsize = (int)ceil((copyRect.maxx - copyRect.minx) / dst_cellsize_x);    dst_ysize = (int)ceil((copyRect.maxy - copyRect.miny) / dst_cellsize_y);    if (dst_xsize == 0 || dst_ysize == 0) {      if (layer->debug)        msDebug("msContourLayerReadRaster(): no apparent overlap between map view and this window(2)./n");      return MS_SUCCESS;    }    if (layer->debug)      msDebug( "msContourLayerReadRaster(): src=%d,%d,%d,%d, dst=%d,%d,%d,%d/n",               src_xoff, src_yoff, src_xsize, src_ysize,               0, 0, dst_xsize, dst_ysize );  } else {    src_xoff = 0;    src_yoff = 0;    dst_xsize = src_xsize = MIN(map->width,src_xsize);    dst_ysize = src_ysize = MIN(map->height,src_ysize);    copyRect.minx = copyRect.miny = 0;    copyRect.maxx = map->width;    copyRect.maxy = map->height;    dst_cellsize_y = dst_cellsize_x = 1;  }  /* -------------------------------------------------------------------- */  /*      Allocate buffer, and read data into it.                         */  /* -------------------------------------------------------------------- */  clinfo->buffer = (double *) malloc(sizeof(double) * dst_xsize * dst_ysize);  if (clinfo->buffer == NULL) {    msSetError(MS_MEMERR, "Malloc(): Out of memory.", "msContourLayerReadRaster()");    return MS_FAILURE;  }  eErr = GDALRasterIO(hBand, GF_Read,                      src_xoff, src_yoff, src_xsize, src_ysize,                      clinfo->buffer, dst_xsize, dst_ysize, GDT_Float64,                      0, 0);  if (eErr != CE_None) {    msSetError( MS_IOERR, "GDALRasterIO() failed: %s",                "msContourLayerReadRaster()", CPLGetLastErrorMsg() );    free(clinfo->buffer);    return MS_FAILURE;  }  memset(pointer, 0, sizeof(pointer));  CPLPrintPointer(pointer, clinfo->buffer, sizeof(pointer));  sprintf(memDSPointer,"MEM:::DATAPOINTER=%s,PIXELS=%d,LINES=%d,BANDS=1,DATATYPE=Float64",          pointer, dst_xsize, dst_ysize);  clinfo->hDS = GDALOpen(memDSPointer,  GA_ReadOnly);  if (clinfo->hDS == NULL) {    msSetError(MS_IMGERR,               "Unable to open GDAL Memory dataset.",               "msContourLayerReadRaster()");    free(clinfo->buffer);    return MS_FAILURE;  }  adfGeoTransform[0] = copyRect.minx;  adfGeoTransform[1] = dst_cellsize_x;  adfGeoTransform[2] = 0;  adfGeoTransform[3] = copyRect.maxy;  adfGeoTransform[4] = 0;  adfGeoTransform[5] = -dst_cellsize_y;  clinfo->cellsize = MAX(dst_cellsize_x, dst_cellsize_y);  {    char buf[64];    sprintf(buf, "%lf", clinfo->cellsize);    msInsertHashTable(&layer->metadata, "__data_cellsize__", buf);  }        GDALSetGeoTransform(clinfo->hDS, adfGeoTransform);  return MS_SUCCESS;}

//.........这里部分代码省略.........               "The corner coordinates of the source are not properly "               "aligned on plain latitude/longitude boundaries.");    }/* -------------------------------------------------------------------- *//*      Check image dimensions.                                         *//* -------------------------------------------------------------------- */    if (!((nXSize == 1201 && nYSize == 1201) || (nXSize == 3601 && nYSize == 3601)))    {        CPLError( CE_Failure, CPLE_AppDefined,                   "Image dimensions should be 1201x1201 or 3601x3601.");        return NULL;    }/* -------------------------------------------------------------------- *//*      Check filename.                                                 *//* -------------------------------------------------------------------- */    char expectedFileName[12];    snprintf(expectedFileName, sizeof(expectedFileName), "%c%02d%c%03d.HGT",             (nLLOriginLat >= 0) ? 'N' : 'S',             (nLLOriginLat >= 0) ? nLLOriginLat : -nLLOriginLat,             (nLLOriginLong >= 0) ? 'E' : 'W',             (nLLOriginLong >= 0) ? nLLOriginLong : -nLLOriginLong);    if (!EQUAL(expectedFileName, CPLGetFilename(pszFilename)))    {        CPLError( CE_Warning, CPLE_AppDefined,                   "Expected output filename is %s.", expectedFileName);    }/* -------------------------------------------------------------------- *//*      Write output file.                                              *//* -------------------------------------------------------------------- */    VSILFILE* fp = VSIFOpenL(pszFilename, "wb");    if (fp == NULL)    {        CPLError( CE_Failure, CPLE_FileIO,                  "Cannot create file %s", pszFilename );        return NULL;    }    GInt16* panData = (GInt16*) CPLMalloc(sizeof(GInt16) * nXSize);    GDALRasterBand* poSrcBand = poSrcDS->GetRasterBand(1);    int bSrcBandHasNoData;    double srcBandNoData = poSrcBand->GetNoDataValue(&bSrcBandHasNoData);    for( int iY = 0; iY < nYSize; iY++ )    {        poSrcBand->RasterIO( GF_Read, 0, iY, nXSize, 1,                            (void *) panData, nXSize, 1,                            GDT_Int16, 0, 0 );        /* Translate nodata values */        if (bSrcBandHasNoData && srcBandNoData != SRTMHG_NODATA_VALUE)        {            for( int iX = 0; iX < nXSize; iX++ )            {                if (panData[iX] == srcBandNoData)                    panData[iX] = SRTMHG_NODATA_VALUE;            }        }#ifdef CPL_LSB        GDALSwapWords(panData, 2, nXSize, 2);#endif        if( VSIFWriteL( panData,sizeof(GInt16) * nXSize,1,fp ) != 1)        {            CPLError( CE_Failure, CPLE_FileIO,                      "Failed to write line %d in SRTMHGT dataset./n",                      iY );            VSIFCloseL(fp);            CPLFree( panData );            return NULL;        }        if( pfnProgress && !pfnProgress((iY+1) / (double) nYSize, NULL, pProgressData ) )        {            CPLError( CE_Failure, CPLE_UserInterrupt,                         "User terminated CreateCopy()" );            VSIFCloseL(fp);            CPLFree( panData );            return NULL;        }    }    CPLFree( panData );    VSIFCloseL(fp);/* -------------------------------------------------------------------- *//*      Reopen and copy missing information into a PAM file.            *//* -------------------------------------------------------------------- */    GDALPamDataset *poDS = (GDALPamDataset *)         GDALOpen( pszFilename, GA_ReadOnly );    if( poDS )        poDS->CloneInfo( poSrcDS, GCIF_PAM_DEFAULT);    return poDS;}

GDALDataset *PNMDataset::Create( const char * pszFilename,                                 int nXSize, int nYSize, int nBands,                                 GDALDataType eType,                                 char ** papszOptions ){/* -------------------------------------------------------------------- *//*      Verify input options.                                           *//* -------------------------------------------------------------------- */    if( eType != GDT_Byte && eType != GDT_UInt16 )    {        CPLError( CE_Failure, CPLE_AppDefined,              "Attempt to create PNM dataset with an illegal "              "data type (%s), only Byte and UInt16 supported.",              GDALGetDataTypeName(eType) );        return nullptr;    }    if( nBands != 1 && nBands != 3 )    {        CPLError( CE_Failure, CPLE_AppDefined,                  "Attempt to create PNM dataset with an illegal number"                  "of bands (%d).  Must be 1 (greyscale) or 3 (RGB).",                  nBands );        return nullptr;    }/* -------------------------------------------------------------------- *//*      Try to create the file.                                         *//* -------------------------------------------------------------------- */    VSILFILE *fp = VSIFOpenL( pszFilename, "wb" );    if( fp == nullptr )    {        CPLError( CE_Failure, CPLE_OpenFailed,                  "Attempt to create file `%s' failed.",                  pszFilename );        return nullptr;    }/* -------------------------------------------------------------------- *//*      Write out the header.                                           *//* -------------------------------------------------------------------- */    int nMaxValue = 0;    const char *pszMaxValue = CSLFetchNameValue( papszOptions, "MAXVAL" );    if ( pszMaxValue )    {        nMaxValue = atoi( pszMaxValue );        if ( eType == GDT_Byte && (nMaxValue > 255 || nMaxValue < 0) )            nMaxValue = 255;        else if ( nMaxValue > 65535 || nMaxValue < 0 )            nMaxValue = 65535;    }    else    {        if ( eType == GDT_Byte )            nMaxValue = 255;        else            nMaxValue = 65535;    }    char szHeader[500] = { '/0' };    if( nBands == 3 )        snprintf( szHeader, sizeof(szHeader),                  "P6/n%d %d/n%d/n", nXSize, nYSize, nMaxValue );    else        snprintf( szHeader, sizeof(szHeader),                  "P5/n%d %d/n%d/n", nXSize, nYSize, nMaxValue );    bool bOK = VSIFWriteL( szHeader, strlen(szHeader) + 2, 1, fp ) == 1;    if( VSIFCloseL( fp ) != 0 )        bOK = false;    if( !bOK )        return nullptr;    return        reinterpret_cast<GDALDataset *>( GDALOpen( pszFilename, GA_Update ) );}