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 am using Python 2.6.5 (32bit) with Numpy 1.3 and Gdal 1.9.1 installed on Windows 7 64bit. I am trying to read an 800 MB Imagine (.img) raster dataset into a Numpy array to do some raster algebra, but as soon as I run the following code, Python.exe crashes.

from osgeo import gdal

g = gdal.Open(r'path\to\dataset', gdal.GA_Readonly)
b = g.GetRasterBand(1)
data = b.ReadAsArray()

Python.exe crashes on the b.ReadAsArray() call. I did some Google searching and found articles dated from Gdal 1.6 that mentioned this issue with Windows 7 64bit, but they also mentioned that it had been fixed in the latest development versions at that time.

Has anyone else had this issue? Any workarounds?

UPDATE:

I decided to debug the code in PyDev to try to determine where it is failing. From what I can tell (still not getting any error messages) it is failing on line 22 of gdal_array.py.

_mod = imp.load_module('_gdal_array', fp, pathname, description)

When I step into the line of code above, it brings me into the init.py module of numpy. When I reach the end of the numpy.__init__.py module, it steps out back to the above line of code. Then, when I hit the step into button, which should bring me to the next line within gdal_array.py, the script just terminates with no error messages or anything.

UPDATE #2:

I uninstalled GDAL 1.9.1 and installed GDAL 1.6.1 from the Python Cheeseshop and OSGeo's Windows Binaries. Still had the same issue.

share|improve this question
    
I to was having this issue. Are you using the gdal python bindings from Tamas at gis.internals? If so, move your additions to your PATH to the front. Some other library was causing some issue for me. –  Jay Laura Dec 21 '12 at 2:23
    
I believe I did download from his internals. I will try to adjust my path when I get in the office tomorrow. Thanks for the advice. –  Brian Dec 21 '12 at 3:04
1  
If that does not work, I recently switched to using these packages - lfd.uci.edu/~gohlke/pythonlibs –  Jay Laura Dec 21 '12 at 3:07
    
I tried moving things around in my System Path (PATH environment variable in Windows) with no luck. I also uninstalled my version of GDAL and installed the version of GDAL on the link you provided and still had the same issue. –  Brian Dec 21 '12 at 15:54
    
Hmmm....the version that I linked to was only the bindings, so you still need to have GDAL core from Tamas. If the other calls worked though, that is working fine. Three things to try that are long shots (in order of 'longness'. 1) Update your version of Numpy. 2) gdal_translate to gtiff and try the code on that image. 3) append ReadAsArray() with .astype(numpy.float32). Is the image public? I can test it on my machine. Can you post the stack trace if none of those work? –  Jay Laura Dec 21 '12 at 16:50
show 9 more comments

3 Answers 3

up vote 4 down vote accepted

As suspected by many commenters, this was an issue with my install. Apparently, i was not paying close enough attention when installing GDAL and the Python Bindings.

I installed GDAL Core and plugins (dll's) from gisinternals.com, but somehow I didn't think to install the Python Bindings from there as well. The Python Bindings I installed were from a different site (can't remember which one at this point).

When I reinstalled GDAL and Python Bindings all from gisinternals.com, I was able to successfully ReadAsArray.

Thank you to all who commented and answered and I apologize for my ignorance.

share|improve this answer
add comment

It's possible it's a memory issue. When you use ReadAsArray it brings the data into memory, and while 800mb isn't massive, it's not tiny either. Have you tried reading the array in chunks?

data = b.ReadAsArray(x_offset, y_offset, x_size, y_size)

You should be able to loop through the array and process it chunk at a time, though depending on the processing you're doing you should probably look to reading in areas with overlap to avoid edge effects.

share|improve this answer
    
I did try using chunks. I tried data = b.ReadAsArray(0,0, 500, 500) with the same result. –  Brian Dec 21 '12 at 0:13
    
Hmm. I'm assuming you've tried other image formats? Also, was there a specific error message? –  om_henners Dec 21 '12 at 0:32
    
I haven't tried any other formats yet. there was no error message, just a popup that said "python.exe has stopped working". –  Brian Dec 21 '12 at 3:02
    
I converted the .img file to a GeoTIFF this morning and tried again. No luck. –  Brian Dec 21 '12 at 15:52
add comment

Sorry I'm late to this party, but your basic problem is that 32-bit Python can't store very large rasters in memory. You can read your large raster into memory in bite-size chunks, but then you're fairly limited in terms of what you can process efficiently without extremely inefficient/frequent read/writes to disk.

What I do instead (which sacrifices some efficiency due to disk read/writes) is to call (via EXE-wrapping) the 64-bit gisinternals.com version of whatever gdal method you need. Be careful with using python's subprocess module in a loop (i.e. you may want/need to call subprocess sequentially) as you may inadvertently generate too many open threads for your windows box and get ominous system warnings. You sacrifice a little bit in the way to disk read/writes with this gdal approach, but your processing efficiency only drops (i.e. relative to a lightening fast in-memory computation, if your box/library could support it) by about a factor or 10.

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.