GDALDriverH hDriver; GDALDataType eDT; GDALDatasetH hDstDS; GDALDatasetH hSrcDS; GDALRasterBandH hSrcbandDS; // Open the source file. GDALAllRegister(); hSrcDS = GDALOpen( "D://7.tif", GA_ReadOnly ); CPLAssert( hSrcDS != NULL ); // Create output with same datatype as first input band. hSrcbandDS = GDALGetRasterBand(hSrcDS,1); eDT = GDALGetRasterDataType(hSrcbandDS); // Get output driver (GeoTIFF format) hDriver = GDALGetDriverByName( "GTiff" ); CPLAssert( hDriver != NULL ); // Get Source coordinate system. const char *pszSrcWKT; char*pszDstWKT = NULL; pszSrcWKT = GDALGetProjectionRef( hSrcDS ); CPLAssert( pszSrcWKT != NULL && strlen(pszSrcWKT) > 0 ); // Get approximate output georeferenced bounds and resolution for file. double adfDstGeoTransform[6]; int nPixels=0, nLines=0; CPLErr eErr; nPixels=6000; nLines=6000; // Create the output file. hDstDS = GDALCreate( hDriver, "D://11.tif", nPixels, nLines, GDALGetRasterCount(hSrcDS), eDT, NULL ); CPLAssert( hDstDS != NULL ); // Write out the projection definition. // Setup output coordinate system that is UTM 11 WGS84. OGRSpatialReference oSRS; oSRS.SetUTM( 50, TRUE ); oSRS.SetWellKnownGeogCS( "WGS84" ); oSRS.exportToWkt(&pszDstWKT ); // Create a transformer that maps from source pixel/line coordinates // to destination georeferenced coordinates (not destination // pixel line). We do that by omitting the destination dataset // handle (setting it to NULL). void *hTransformArg; hTransformArg = GDALCreateGenImgProjTransformer( hSrcDS, pszSrcWKT, NULL, pszDstWKT, FALSE, 1000, 1 ); eErr = GDALSuggestedWarpOutput( hSrcDS, GDALGenImgProjTransform, hTransformArg, adfDstGeoTransform, &nPixels, &nLines ); CPLAssert( hTransformArg != NULL ); GDALDestroyGenImgProjTransformer( hTransformArg ); adfDstGeoTransform[1]=adfDstGeoTransform[1]/2; adfDstGeoTransform[5]=adfDstGeoTransform[5]/2; CPLAssert( eErr == CE_None ); GDALSetProjection( hDstDS, pszDstWKT ); GDALSetGeoTransform( hDstDS, adfDstGeoTransform ); // Copy the color table, if required. GDALColorTableH hCT; hCT = GDALGetRasterColorTable( GDALGetRasterBand(hSrcDS,1) ); if( hCT != NULL ) GDALSetRasterColorTable( GDALGetRasterBand(hDstDS,1), hCT ); // Setup warp options. GDALWarpOptions *psWarpOptions = GDALCreateWarpOptions(); psWarpOptions->hSrcDS = hSrcDS; psWarpOptions->hDstDS = hDstDS; psWarpOptions->nBandCount = 3; psWarpOptions->panSrcBands = (int *) CPLMalloc(sizeof(int) * psWarpOptions->nBandCount ); psWarpOptions->panSrcBands[0] = 1; psWarpOptions->panSrcBands[1] = 2; psWarpOptions->panSrcBands[2] = 3; psWarpOptions->panDstBands = (int *) CPLMalloc(sizeof(int) * psWarpOptions->nBandCount ); psWarpOptions->panDstBands[0] = 1; psWarpOptions->panDstBands[1] = 2; psWarpOptions->panDstBands[2] = 3; psWarpOptions->pfnProgress = GDALTermProgress; // Establish reprojection transformer. psWarpOptions->eResampleAlg=GRA_NearestNeighbour; // psWarpOptions->pTransformerArg = // GDALCreateGenImgProjTransformer( hSrcDS, // GDALGetProjectionRef(hSrcDS), // hDstDS, // GDALGetProjectionRef(hDstDS), // FALSE, 1000.0, 3 ); // psWarpOptions->pfnTransformer = GDALGenImgProjTransform; void * m_hGenTransformArg=NULL; void * m_hApproxArg=NULL; m_hGenTransformArg = GDALCreateGenImgProjTransformer( hSrcDS, GDALGetProjectionRef(hSrcDS), hDstDS, GDALGetProjectionRef(hDstDS), FALSE, 1000.0, 3 ); if( m_hGenTransformArg == NULL ) return CE_Failure; m_hApproxArg = GDALCreateApproxTransformer( GDALGenImgProjTransform, m_hGenTransformArg, 0.125 ); if ( m_hApproxArg == NULL) return CE_Failure; psWarpOptions->pfnTransformer = GDALApproxTransform; psWarpOptions->pTransformerArg = m_hApproxArg; // Initialize and execute the warp operation. GDALWarpOperation oOperation; oOperation.Initialize( psWarpOptions ); // int Typesize=GDALGetDataTypeSize(GDT_Byte); // int Dstwidth,DstHeight; // int Srcwidth,Srcheight; // Dstwidth=500; // DstHeight=500; // Srcheight=1000; // Srcwidth=1000; // void *pDstBuffer = new byte[ Typesize* Srcheight* Srcwidth*GDALGetRasterCount(hDstDS) ]; // memset(pDstBuffer, 0,Typesize* Srcheight* Srcwidth*GDALGetRasterCount(hDstDS)); // int nJpgBandCount = GDALGetRasterCount(hDstDS); // int *panBands = new int[nJpgBandCount]; // for (int i = 0; i < nJpgBandCount; i++) // { // panBands[i] = i+1; // } // // //AfxMessageBox("123"); // CPLErr err = GDALDatasetRasterIO( hSrcDS, GF_Read, // 0, 0, Srcwidth, Srcheight, // pDstBuffer, Srcwidth, Srcheight, // GDT_Byte, nJpgBandCount, panBands, // 0, 0, 0 ); // // oOperation.WarpRegionToBuffer(0,0,Dstwidth,DstHeight,pDstBuffer,GDT_Byte); // // eErr = GDALDatasetRasterIO( hDstDS, GF_Write, // 0, 0, Dstwidth, DstHeight, // pDstBuffer, Srcwidth, Srcheight, // GDT_Byte, // 3, // panBands, // 0, 0, 0 ); // delete[] panBands; oOperation.ChunkAndWarpImage( 0, 0, GDALGetRasterXSize(hDstDS), GDALGetRasterYSize(hDstDS)); GDALDestroyGenImgProjTransformer( m_hGenTransformArg ); GDALDestroyWarpOptions( psWarpOptions ); GDALDestroyApproxTransformer(m_hApproxArg); GDALClose( hDstDS ); GDALClose( hSrcDS ); return TRUE; BOOL CMapDivisionDoc::Interpolation(CString strFileTemp,CString strFile) { GDALDriverH hDriver; GDALDataType eDT; GDALDatasetH hDstDS; GDALDatasetH hSrcDS; GDALRasterBandH hSrcbandDS; // Open the source file. GDALAllRegister(); hSrcDS = GDALOpen(strFileTemp, GA_ReadOnly ); CPLAssert( hSrcDS != NULL ); // Create output with same datatype as first input band. hSrcbandDS = GDALGetRasterBand(hSrcDS,1); eDT = GDALGetRasterDataType(hSrcbandDS); // Get output driver (GeoTIFF format) char *pszFormat=GetOutimgInfo(strFileTemp); hDriver = GDALGetDriverByName(pszFormat); CPLAssert( hDriver != NULL ); // Get Source coordinate system. const char *pszSrcWKT; char*pszDstWKT = NULL; pszSrcWKT = GDALGetProjectionRef( hSrcDS ); CPLAssert( pszSrcWKT != NULL && strlen(pszSrcWKT) > 0 ); // Get approximate output georeferenced bounds and resolution for file. double adfDstGeoTransform[6]; int nPixels=0, nLines=0; CPLErr eErr; nPixels=m_nFileSize; nLines=m_nFileSize; // Create the output file. hDstDS = GDALCreate( hDriver, strFile, nPixels, nLines, GDALGetRasterCount(hSrcDS), eDT, NULL ); CPLAssert( hDstDS != NULL ); // Write out the projection definition. // Setup output coordinate system that is UTM 11 WGS84. // OGRSpatialReference oSRS; // oSRS.SetUTM( 50, TRUE ); // oSRS.SetWellKnownGeogCS( "WGS84" ); // oSRS.exportToWkt(&pszDstWKT ); // Create a transformer that maps from source pixel/line coordinates // to destination georeferenced coordinates (not destination // pixel line). We do that by omitting the destination dataset // handle (setting it to NULL). void *hTransformArg; hTransformArg = GDALCreateGenImgProjTransformer( hSrcDS, pszSrcWKT, NULL, pszSrcWKT, FALSE, 1000, 1 ); eErr = GDALSuggestedWarpOutput( hSrcDS, GDALGenImgProjTransform, hTransformArg, adfDstGeoTransform, &nPixels, &nLines ); CPLAssert( hTransformArg != NULL ); GDALDestroyGenImgProjTransformer( hTransformArg ); adfDstGeoTransform[1]=adfDstGeoTransform[1]*GDALGetRasterXSize(hSrcDS)/m_nFileSize; adfDstGeoTransform[5]=adfDstGeoTransform[5]*GDALGetRasterXSize(hSrcDS)/m_nFileSize; CPLAssert( eErr == CE_None ); GDALSetProjection( hDstDS, pszSrcWKT ); GDALSetGeoTransform( hDstDS, adfDstGeoTransform ); // Copy the color table, if required. GDALColorTableH hCT; hCT = GDALGetRasterColorTable( GDALGetRasterBand(hSrcDS,1) ); if( hCT != NULL ) GDALSetRasterColorTable( GDALGetRasterBand(hDstDS,1), hCT ); // Setup warp options. GDALWarpOptions *psWarpOptions = GDALCreateWarpOptions(); psWarpOptions->hSrcDS = hSrcDS; psWarpOptions->hDstDS = hDstDS; psWarpOptions->nBandCount = GDALGetRasterCount(hDstDS) ; psWarpOptions->panSrcBands = (int *) CPLMalloc(sizeof(int) * psWarpOptions->nBandCount ); for (int i=0;i<GDALGetRasterCount(hDstDS);i++) { psWarpOptions->panSrcBands[i] = i+1; } // psWarpOptions->panSrcBands[0] = 1; // psWarpOptions->panSrcBands[1] = 2; // psWarpOptions->panSrcBands[2] = 3; psWarpOptions->panDstBands = (int *) CPLMalloc(sizeof(int) * psWarpOptions->nBandCount ); for (int m=0;m<GDALGetRasterCount(hDstDS);m++) { psWarpOptions->panDstBands[m] = m+1; } // psWarpOptions->panDstBands[0] = 1; // psWarpOptions->panDstBands[1] = 2; // psWarpOptions->panDstBands[2] = 3; psWarpOptions->pfnProgress = GDALTermProgress; // Establish reprojection transformer. psWarpOptions->eResampleAlg=GRA_Bilinear; void * m_hGenTransformArg=NULL; void * m_hApproxArg=NULL; m_hGenTransformArg = GDALCreateGenImgProjTransformer( hSrcDS, GDALGetProjectionRef(hSrcDS), hDstDS, GDALGetProjectionRef(hDstDS), FALSE, 1000.0, 3 ); if( m_hGenTransformArg == NULL ) return CE_Failure; m_hApproxArg = GDALCreateApproxTransformer( GDALGenImgProjTransform, m_hGenTransformArg, 0.125 ); if ( m_hApproxArg == NULL) return CE_Failure; psWarpOptions->pfnTransformer = GDALApproxTransform; psWarpOptions->pTransformerArg = m_hApproxArg; // Initialize and execute the warp operation. GDALWarpOperation oOperation; oOperation.Initialize( psWarpOptions ); oOperation.ChunkAndWarpMulti( 0, 0, GDALGetRasterXSize(hDstDS), GDALGetRasterYSize(hDstDS)); GDALDestroyGenImgProjTransformer( m_hGenTransformArg ); GDALDestroyWarpOptions( psWarpOptions ); GDALDestroyApproxTransformer(m_hApproxArg); GDALClose( hDstDS ); GDALClose( hSrcDS ); remove(strFileTemp); return TRUE; }