Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

I have a flow accumulation grid and a shapefile containing catchment areas. Now I want to do some calculations on them so I need to extract data from cells in the raster adjacent to certain catchment areas.

In that process I try to use Extract by Mask on my raster, with one of the catchment polygons as mask (in a loop later on, but now I just need it to work once).

The problem is that the extracted raster is shifted half a cell to the left, compared to both the original raster and the masking polygon. This is messing up the rest of my process, since I try to convert the raster to points and intersect with other polygons later on. And unless the extracted raster is overlapping the original one, I wont get the right polygon.

The original raster and the polygons are aligning perfectly, since the polygons are created from the raster. I have tried to check if there somehow have been some kind of change in the projections, but can't find any. I have also tried in new ArcMap documents, restarting everything several times (this solves the problem surprisingly often) and executing the tools from the python prompt instead of calling my script. I have even tried to execute the Extract by Mask tool from the Toolbox, with same result.

If anyone can help me with a fix, or at least tell me what I am doing wrong, that would be greatly appreciated!

This is the (very simple) code that I use in the python prompt in ArcMap:

inpgs = r'C:\GIS\Catchm_ex.shp'
flowacc = r'C:\GIS\flowacc'
arcpy.MakeFeatureLayer_management(inpgs, 'currentMask', '"HydroID" = 269634')
outRaster = arcpy.sa.ExtractByMask(flowacc, 'currentMask')

I have some sample files here

I use Windows 7, ArcInfo 10.0 and Python 2.6.

EDIT: I have tried Extract by Polygon now, and that seems to work fine. Is it just my machine malfunctioning on ExtractByMask? Only drawback is that I'm having problems using that function in my script, since I can't work out how to get the coordinates of my masking polygon (only covering 9 cells, with a point in center if that helps).

EDIT2: I have now got it working with Extract by polygon by using the code below (modified from Esri's support pages). Unfortunately, this does not solve my problem. It appears as if the output (from this tool as well as Extract by mask) is sometimes shifted half a cell to the left, and sometimes not. Mostly it is, but occasionally it actually is located where it should be.

So now I'm really confused. Is this just happening to me, on my computer? If anyone else have experienced the same, please write a comment so I know that I'm not all alone :) To me this problem seems to be a too big of an issue to be a bug, which leaves me with this being a problem that I have created myself along the way. Any ideas on how? Or how to avoid it?

def ReadPoints(infc):
    try:
        import arcpy, numpy
    except ImportError:
        sys.exit('Unable to import arcpy')

    # Identify the geometry field
    desc = arcpy.Describe(infc)
    shapefieldname = desc.ShapeFieldName
    points = []

    # Create search cursor
    rows = arcpy.SearchCursor(infc)

    # Enter for loop for each feature/row
    for row in rows:
        # Create the geometry object
        feat = row.getValue(shapefieldname)
        partnum = 0

        # Step through each part of the feature
        for part in feat:
            # Step through each vertex in the feature
            for pnt in feat.getPart(partnum):
                if pnt:
                    points.append([int(pnt.X), int(pnt.Y)])

    pt = numpy.array(points)
    pts = []
    for i in range(len(pt)):
        t = pt[i]
        pts.append(arcpy.Point(t[0],t[1]))
    return pts
share|improve this question
add comment

3 Answers

up vote 1 down vote accepted

This may not be helpful now, but this is a known issue with the Extract by Mask tool that was obviously not fixed for ArcGIS 10: http://support.esri.com/en/knowledgebase/techarticles/detail/34448

I don't know the scripting solution . . .

share|improve this answer
add comment

Is it possible to convert the whole image to a points and then select the points based on interesctions between polygon and points layers?. I know this does not give any insight to strange results you are getting but it might be a useful "work around"

share|improve this answer
    
Yes, i could do that. Only problem is that the raster is a bit too big, so it'll take quite some time to compute. I might try that anyway, as opposed to not have it working at all though. Thanks for your input! –  Martin Nov 20 '12 at 12:50
    
Not sure if this will help but you can import symbology, so maybe you could set the symbology for one and then import it to all the rest. Have you had a look at the image Analysis window, it has a lot of features that might help. Also consider creating a mosaic dataset if the images are large. With a mosaic you can add functions that could improve the visual display. hope that helps :) –  dango Nov 20 '12 at 13:38
    
I have not worked much with Image analysis before, but now I might have to give it a try. I'll report back if this can get me anywhere :) –  Martin Nov 20 '12 at 14:05
add comment

In your raster environments, set the snap raster to be the same as your input raster (the one being extracted). This will avoid any shifting.

share|improve this answer
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.