<div dir="ltr">When attempting to use a polygon shapefile (in custom web mercator projection) with the gdalwarp cutline operation, operating on a GeoTIFF (with same projection), the following error is output in both versions 1.10.0 and 1.10.1:<div>
<br></div><div> ERROR 6: No translation for 'Mercator_Auxiliary_Sphere' to PROJ.4 format is known.</div><div><br></div><div>With GDAL 1.10.0, the operation still completes successfully with a zero exit status, but with GDAL 1.10.1 the operation fails with a nonzero exit status, and the resulting GeoTIFF has zeros for all pixels.</div>
<div><br></div><div>Here is the command I'm using:</div><div><br></div><div> gdalwarp -co COMPRESS=LZW -cutline oceans_glakes_50m_3857.shp -srcnodata 9999.0 -dstnodata 9999.0 -overwrite sst_global.tif sst_global_cut.tif</div>
<div><br></div><div><br></div><div>The input shapefile was produced by reprojecting with ogr2ogr from EPSG:4326 to my custom projection (WKT is included below), which is stored in a *.prf file. The source GeoTIFF is single-band with Float32 pixel type.</div>
<div><br></div><div>The custom WKT can be found here: <a href="http://jgreenlaw.org/gdal/esri-3857.prf">http://jgreenlaw.org/gdal/esri-3857.prf</a></div><div><br></div><div>The cutline shapefile can be downloaded at: <a href="http://jgreenlaw.org/gdal/oceans_glakes_50m_3857.tar.gz">http://jgreenlaw.org/gdal/oceans_glakes_50m_3857.tar.gz</a> (900kb)</div>
<div><br></div><div>The source GeoTIFF can be generated from this FLT/HDR raster: <a href="http://jgreenlaw.org/gdal/sst_flt.tar.gz">http://jgreenlaw.org/gdal/sst_flt.tar.gz</a> (11mb) by running these commands:</div><div>
<br></div><div> gdalwarp -s_srs EPSG:4326 -t_srs esri-3857.prf -r bilinear -srcnodata 9999.0 -dstnodata 9999.0 -tr 4000 9000 -overwrite -co COMPRESS=LZW -co TILED=yes sst.flt sst_projected.tif</div><div><br></div><div>
gdalwarp -co COMPRESS=LZW -te -20037507.0671618 -19971868.8804086 20037507.0671618 19971868.8804086 -srcnodata 9999.0 -dstnodata 9999.0 -tr 4000 9000 sst_projected.tif sst_global.tif</div><div><br></div><div>Does anyone know if this behavior is to be expected when the above error is output? Is there a better method for reprojecting to Web Mercator that will align properly in ArcGIS?</div>
<div><br></div><div>Thanks</div><div>-Jason</div><div><br></div><div><br></div><div>Background:</div><div>The reason for using custom WKT for reprojecting to Web Mercator is because I have experienced alignment problems displaying GDAL-produced GeoTIFFs in ArcGIS when reprojecting using "EPSG:3857" as the target SRS (this seems to be related to <a href="http://trac.osgeo.org/gdal/ticket/3962">http://trac.osgeo.org/gdal/ticket/3962</a>). </div>
<div><br></div><div>I generated the WKT below by copying the 102113 definition in $GDAL_DATA/esri_extra.wkt and merging information from an ESRI Web Mercator *.prj file. Using this custom WKT as the target SRS with gdalwarp normally results in rasters that align properly with ArcGIS Web Mercator data (without on-the-fly datum transformation).</div>
<div><br></div><div><br></div><div>----------------------------------------------------------------------------------------------</div><div>Custom WKT for reprojecting to ESRI Web Mercator using GDAL (<a href="http://jgreenlaw.org/gdal/esri-3857.prf">http://jgreenlaw.org/gdal/esri-3857.prf</a>):</div>
<div>----------------------------------------------------------------------------------------------</div><div><br></div><div><div>PROJCS['WGS_1984_Web_Mercator_Auxiliary_Sphere',</div><div> GEOGCS['GCS_WGS_1984',</div>
<div> DATUM['D_WGS_1984',</div><div> SPHEROID['WGS_1984',6378137.0,298.257223563]</div><div> ],</div><div> PRIMEM['Greenwich',0.0],</div><div> UNIT['Degree',0.0174532925199433]</div>
<div> ],</div><div> PROJECTION['Mercator_Auxiliary_Sphere'],</div><div> PARAMETER['False_Easting',0.0],</div><div> PARAMETER['False_Northing',0.0],</div><div> PARAMETER['Central_Meridian',0.0],</div>
<div> PARAMETER['latitude_of_origin',0.0],</div><div> PARAMETER['Standard_Parallel_1',0.0],</div><div> PARAMETER['Auxiliary_Sphere_Type',0.0],</div><div> UNIT['Meter',1.0],</div>
<div> EXTENSION["PROJ4",</div><div> "+proj=merc +a=6378137 +b=6356752.314245179 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0.0 +k=1.0 +units=m +nadgrids=@null +wktext +no_defs"</div><div> ],</div>
<div> AUTHORITY["EPSG","3857"]</div><div>]</div></div><div><br></div></div>