[GIS] Arcpy/REST: how to save raster dataset from URL

arcgis-rest-apiarcgis-serverarcpydownloadraster

I have a script which downloads features from an ArcGIS MapServer via the REST API. Getting vector features is relatively trivial since the features can be retrieved as JSON. However, when retrieving a map as an image, I'm getting stuck at the point where I have retrieved a response.

The relevant portion of my code, which gives me a valid image URL, is this:

# Export the map to JSON.
urlGetImage = '%s/export?bbox=%s&bboxSR=%s&format=png&f=json' % (
    theBasicUrl, bbox, bboxSR)
urlGetImageResponse = urllib2.urlopen(urlGetImage)
jsonImage = json.loads(urlGetImageResponse.read())

# Read the image URL from the JSON and save the image.
urllib.urlretrieve(jsonImage['href'], os.path.join(tempFolder, 'test.png'))

Note the parameter where I retrieve the response as JSON. This is not essential. I don't care what I retrieve it as, so long as I can retrieve the generated raster somehow.

A sample response is this:

{"width": 400, "href": "https://[SERVER NAME]:[PORT]/arcgis/rest/directories/arcgisoutput/imagery/DEM_MapServer/_ags_map26932bd6681a4546a951125dd7c60e88.png", "scale": 20326865.193372156, "extent": {"xmin": 2996004.971208999, "ymin": -3232131.5223958646, "ymax": -1080867.3202355732, "xmax": 5147269.17336929, "spatialReference": {"wkid": 102100, "latestWkid": 3857}}, "height": 400}

This gives me a path to download an image from, as well as the relevant spatial reference WKID. I then save it to an image.

My question is this: how do I go from here (an image on the drive) to a (referenced) raster dataset? Or have I already gone wrong in getting to this point?

EDIT:

To clarify, I need to work with a MapService, not an ImageService, as that's how the rasters will be published. When I export the result above as a KMZ (rather than JSON), and then load it in Google Earth, it appears in the correct location. So in theory, getting a raster image from the map server would seem possible. However, trying to convert this into a raster layer via arcpy.KMLToLayer_conversion gives:

  "ERROR 000401: No features were found for processing."

… and even if I extract the image information in the KML, download the image, and replace the path in the KML (which by default points to the generated image on the server), it still complains with this error.

Best Answer

Are you working with an ImageService? From what I can tell is that you are working with a map service and if that is the case, the map export option will only return a static image of the map service or individual layer.

Looking at what gets returned in your example, we can see there is a png being returned. This is not a true raster but rather just an image of it. Therefore, it will not be georectified and it will not have the pixel values of the raster (they will just be RGB values).

These appear to be services on your own servers? If so, you can just publish the raster by itself as an image service and then you can get the image back as a tif, which will have the spatial reference information as well as the pixel values.

EDIT

After further discussion, we determined that you just needed to just georeference the output map image. I have written a small function based on your code above to do this and it worked for me:

import urllib, urllib2, json, os, arcpy
arcpy.env.overwriteOutput = True

def exportMap(url, bbox, bboxSR, out_png):

    # get bounding box properties (assumes is comma delimited)
    xmin, ymin, xmax, ymax = map(float, bbox.split(','))
    size = ','.join(map(str, [int(abs(xmin-xmax)), int(abs(ymin-ymax))]))

    # make requests         
    urlGetImage = '%s/export?bbox=%s&bboxSR=%s&size=%s&format=png&f=json' %(url, bbox, bboxSR, size)
    urlGetImageResponse = urllib2.urlopen(urlGetImage)
    jsonImage = json.loads(urlGetImageResponse.read())

    # Read the image URL from the JSON and save the image.
    tmp_png = os.path.join(arcpy.env.scratchFolder, 'temp_png.png')
    urllib.urlretrieve(jsonImage['href'], tmp_png)

    # read the raster extent with arcpy and determine shift params
    ext = arcpy.Describe(tmp_png).extent
    ext_json = jsonImage['extent']
    x_shift = ext_json['xmin'] + ext.XMin
    y_shift = ext_json['ymax'] + ext.YMax #might have to adjust these?

    # run shift tool
    arcpy.env.outputCoordinateSystem = ext_json['spatialReference']['wkid']
    arcpy.management.Shift(tmp_png, out_png, x_shift, y_shift)
    arcpy.management.Delete(tmp_png)
    return out_png

if __name__ == '__main__':

    # test it
    url = 'http://arcgis.dnr.state.mn.us/arcgis/rest/services/elevation/mn_dem_1m/MapServer'
    bbox = '230720.3137,4954797.5862,232077.6289,4955591.3378'
    sr = 26915
    out_png = r'C:\TEMP\Map_export.png'
    exportMap(url, bbox, sr, out_png)
    print 'done'

note: if the map service is tiled, it will not align exactly because it is pulling out the actual tiles for the output image.

Related Question