Please disregard I found what I needed to know in here: http://www.gdal.org/warptut.html
and in here: http://trac.osgeo.org/gdal/browser/trunk/gdal/alg/gdalwarper.cpp and I now have it coming out correctly tilted. Thank You, Andy Canfield On Mon, Mar 25, 2013 at 5:51 PM, Andy Canfield <[email protected]>wrote: > I am attempting to tranform a GeoTiff from UTM14N to WGS84 using GDAL from > VC++ 2008 and I am brand new to using GDAL. > > When I do the transformation everything seems to go off correctly but when > I compare my transformed image in Quantum GIS to the source image projected > on the fly by Quantum GIS their image is slightly rotated and mine is > straight. Theirs is correct, mine is not(I believe), yet I am not sure what > I am doing incorrectly. Is it my affine translation parameters, is it that > I am warping it improperly, did I calculate corners improperly or something > else entirely? > > The source image is a 3 band GeoTiff if that helps at all. > > Here is my complete code as I have it right now: > > #include "stdafx.h" > #include "gdalwarper.h" > #include "cpl_conv.h" > #include "gdal_priv.h" > #include <ogr_spatialref.h> > #include "ogr_srs_api.h" > > int _tmain(int argc, _TCHAR* argv[]) > { > GDALAllRegister(); > > GDALDataset *poDataset, *poDstDS; > const char *pszFormat = "GTiff"; > const char *pszSrcFilename = "C:\\Example\\SampleData\\In.tif"; > const char *pszDstFilename = "C:\\Example\\SampleData\\Out.tif"; > GDALDriver *poDriver; > int rasterCount; > > OGRSpatialReference oSrcSRS,oDstSRS; > char *pszSrcSRSWKT = NULL; > char *pszDstSRSWKT = NULL; > > OGRCoordinateTransformation *transform = NULL; > > double adfSrcGeoTransform[6]; > double adfDstGeoTransform[6]; > double xS[2]; > double yS[2]; > > //Set destination spatial reference to WGS84 and get its WKT > oDstSRS.SetWellKnownGeogCS("WGS84"); > oDstSRS.exportToWkt(&pszDstSRSWKT); > > //Now set the source spatial reference to UTM14N WGS84 and get its WKT > oSrcSRS.SetProjCS( "UTM 14 (WGS84) in northern hemisphere." ); > oSrcSRS.SetWellKnownGeogCS( "WGS84" ); > oSrcSRS.SetUTM( 14, TRUE ); > oSrcSRS.exportToWkt(&pszSrcSRSWKT); > > //Set the coordinate transformation > transform = OGRCreateCoordinateTransformation(&oSrcSRS,&oDstSRS); > > //Set the GeoTiff driver > poDriver = GetGDALDriverManager()->GetDriverByName(pszFormat); > > // Open input and output files. > poDataset = (GDALDataset *) GDALOpen( pszSrcFilename, GA_ReadOnly ); > rasterCount = poDataset->GetRasterCount(); > > poDstDS = poDriver->CreateCopy( pszDstFilename, poDataset, FALSE, NULL, > NULL, NULL ); > > // Setup warp options. > GDALWarpOptions *psWarpOptions = GDALCreateWarpOptions(); > > psWarpOptions->hSrcDS = poDataset; > psWarpOptions->hDstDS = poDstDS; > > > psWarpOptions->nBandCount = rasterCount; > psWarpOptions->panSrcBands = (int *) CPLMalloc(sizeof(int) * > psWarpOptions->nBandCount ); > psWarpOptions->panDstBands = (int *) CPLMalloc(sizeof(int) * > psWarpOptions->nBandCount ); > > for(int i=0;i<rasterCount;i++) > { > psWarpOptions->panSrcBands[i] = i+1;//src band number > psWarpOptions->panDstBands[i] = i+1;//dst band number > } > > psWarpOptions->pfnProgress = GDALTermProgress; > > // Establish reprojection transformer. > > psWarpOptions->pTransformerArg = > GDALCreateGenImgProjTransformer( poDataset, > pszSrcSRSWKT, > poDstDS, > pszDstSRSWKT, > FALSE, 0.0, 1 ); > > psWarpOptions->pfnTransformer = GDALGenImgProjTransform; > psWarpOptions->eResampleAlg = GRA_Cubic; > > // Initialize and execute the warp operation. > GDALWarpOperation oOperation; > > oOperation.Initialize( psWarpOptions ); > oOperation.ChunkAndWarpImage( 0, 0, GDALGetRasterXSize( poDstDS ), > GDALGetRasterYSize( poDstDS ) ); > > GDALDestroyGenImgProjTransformer( psWarpOptions->pTransformerArg ); > GDALDestroyWarpOptions( psWarpOptions ); > > //Set the new image's projection tag info > poDstDS->SetProjection(pszDstSRSWKT); > > //Get the original affine parameters > poDataset->GetGeoTransform(adfSrcGeoTransform); > > double srcWidth = poDataset->GetRasterXSize(); > double srcHeight = poDataset->GetRasterYSize(); > > > //Get the corners > xS[0] = adfSrcGeoTransform[0];//minx > xS[1] = adfSrcGeoTransform[0] + srcWidth*adfSrcGeoTransform[1] + > srcHeight*adfSrcGeoTransform[2];//??? maxx > > yS[0] = adfSrcGeoTransform[3];//maxy > yS[1] = adfSrcGeoTransform[3] + srcWidth*adfSrcGeoTransform[4] + > srcHeight*adfSrcGeoTransform[5];//??? miny > > > //Reproject its corners coordinate info > transform->Transform(2,xS,yS); > > //Set all the new affine parameters > adfDstGeoTransform[0] = xS[0]; > adfDstGeoTransform[1] = (xS[1] - xS[0]) / srcWidth;//(xMax - xMin) / > imageWidth > adfDstGeoTransform[2] = adfSrcGeoTransform[2]; > adfDstGeoTransform[3] = yS[0]; > adfDstGeoTransform[4] = adfSrcGeoTransform[4]; > adfDstGeoTransform[5] = ((yS[0] - yS[1]) / srcHeight)*(-1);//((yMax - > yMin) / imageHeight)*(-1) > > //Set the new images affine information > poDstDS->SetGeoTransform(adfDstGeoTransform); > > GDALClose( poDstDS ); > GDALClose( poDataset ); > > return 0; > } > > Thank you in advance for any help you folks may be able to give me on this. > > -Andy >
_______________________________________________ gdal-dev mailing list [email protected] http://lists.osgeo.org/mailman/listinfo/gdal-dev
