Tell me more ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.
up vote 3 down vote favorite

In terms of processing speed what's the fastest way of creating and then adding values to a raster?

The values will originally come from ASCII files which can be about 500MB in size - ideally I would like to create an Esri GRID but that is not essential.

I have Python (GDAL) and/or ArcOjects (9.3.1) in my arsenal

I cannot use the ASCII to Raster tool due to the format of the values in the ASCII file . I'm looking for a solution of beaing able to write values to the raster directly.

As an example I may have an XYZ asccii as like:

    @File_Version: 3
@Coordinate_Type_is: 3
@Export_Type_is: 1
@Number_of_Projects 1
@Project_Type_Name: 3 cstone25 
@Project_Unit_is: meters
#File_Version____________-> 3
#Project_Name____________-> cstone25
#Project_Type____________-> 3
#Project_Unit____________-> meters
#OpenWorks_Project_______-> EXPLORATION_UTM31
#Master_Project_______-> 
#Coordinate_type_________-> 3
#Number_of_points_in_hzd_-> 1
zpxc1_ufa_nra_env_sele+20+100_Max_Peak_Amp
#Horizon_internal_id_____->_120
#Horizon_extremes_are____->   7245798912.00000
#Horizon_onset_is__Minimum_____-> 1
#Horizon_type_is_Time______-> 1
#Horizon_color_is________-> 0 255 0
 6148.0 17392.0     529940.17623    6274994.18251 394123616.00000
 6148.0 17396.0     529940.83286    6275019.17388 326897280.00000
 6148.0 17400.0     529941.48949    6275044.16526 567459712.00000
 6148.0 17404.0     529942.14613    6275069.15663 429482368.00000
 6148.0 17408.0     529942.80276    6275094.14801 133182280.00000
 6148.0 17412.0     529943.45939    6275119.13938 32126824.00000
 6148.0 17416.0     529944.11603    6275144.13076 17940496.00000
 6148.0 17420.0     529944.77266    6275169.12213 320429248.00000
 6148.0 17424.0     529945.42929    6275194.11351 453848192.00000
 6148.0 17428.0     529946.08592    6275219.10488 336978496.00000
 6148.0 17432.0     529946.74256    6275244.09626 290245120.00000
 6148.0 17436.0     529947.39919    6275269.08763 323590976.00000
 6148.0 17440.0     529948.05582    6275294.07901 249979616.00000
 6148.0 17444.0     529948.71246    6275319.07038 221433360.00000
 6148.0 17448.0     529949.36909    6275344.06176 244277504.00000

or I might have a header than values:

!
!     FILE NAME :  Rob_Grid                
!     FORMATTED FILE CREATION DATE:  APR 25 2011
!     FORMATTED FILE CREATION TIME:  17:11
!
@Rob_Grid HEADER                         ,      GRID,         5
      15,  0.1000000E+31,            ,      7,      1
     245,    232,   262410.0    ,   273960.0    ,   6426990.    ,   6439190.    
    200.0000    ,   0.000000    ,   0.000000    
@
  0.1000000E+31  0.1000000E+31  0.1000000E+31  0.1000000E+31  0.1000000E+31
  0.1000000E+31  0.1000000E+31  0.1000000E+31  0.1000000E+31  0.1000000E+31
  0.1000000E+31  0.1000000E+31  0.1000000E+31   1023.936       1021.771    
   1020.476       1019.986       1019.994       1020.656       1022.158    
   1024.154       1026.193       1027.803       1028.656       1028.808    
   1028.548       1028.080       1027.858       1028.684       1030.780

NOTE: this question is about how to create and populate the raster values the fastest

Thanks Rob

share|improve this question
1  
Two more pieces of information would help, Rob: (1) exactly what is the logical formatting of your ASCII file and (2) does it represent a regular grid of values or (potentially) irregularly spaced points? – whuber♦ May 9 '11 at 15:02
Hi Whuber - Regular spacing so the cell size is known and the XY extents are also known (usually by the header having origin X&Y and then number of rows and columns). I deal with two formats of ASCII files either XYZ or an XY corner & rows/columns/cell size and then all of the values tab or space seperated – Rob Clark May 9 '11 at 15:19
1  
A small sample of what the ascii files actually look like would help. A verbal description is not enough to account for the variable ways it could be interpreted. – matt wilkie♦ May 9 '11 at 16:30
I agree with Matt, it would be fairly simple to program/insert a header into an existing file to produce a file of useable format once the format is known. – Dan Patterson May 10 '11 at 1:26
Can you copy and paste the header in the body of the post? I assume you can't just open and view the file as is in catalog? – Jakub May 10 '11 at 3:13
show 1 more comment

2 Answers

up vote 1 down vote accepted

I'm not sure if it's the fastest, but you might consider using IRasterEdit.

The sample here could be modified to write out to a different format (instead of TIFF), by providing appropriate format to IRasterWorkspace2.CreateRasterDataSet. I think also you might want to loop through instead of writing the entire raster with one call to IRasterEdit.Write, as done in the sample.

I think this can also be done using IPixelBlockCursor, but the remarks give me the impression it's been deprecated: "A legacy interface in ArcGIS 9.2".

share|improve this answer
Thanks Kirk - I'm using IRasterEdit now – Rob Clark May 11 '11 at 15:04
up vote 3 down vote

I built a Python script with GDAL/NumPy that converts some custom ASCII format (similar to yours) into a GeoTIFF for ArcGIS or other GIS software. Everything you need to know to do this is provided in this GDAL tutorial.

Basically, the script (which I can't post) reads the data into NumPy data structures, and does the following:

  1. Determine the GeoTransform tuple from the X and Y origins and grid spacings. I use something like (xdim.min() - dx/2, dx, 0.0, ydim.min() - dy/2, 0.0, dy) (both dx and dy are positive, which is why I use min())
  2. Structure your raster heights into a NumPy array with shape (len(ydim), len(xdim)), then use WriteArray for each band
  3. If you have incomplete raster coverages, be sure to set SetNoDataValue in both the raster and NumPy array

Lastly, as special note for ArcGIS: when you view your data, it may appear black, so you should Calculate Statistics. I've added this step in my Python script too.

share|improve this answer
If performance is your only concern, you can write against the C or C++ APIs directly and shave off a little extra time, but I suspect that the Python bindings will be sufficient to beat out anything Arc*. – scw May 25 '11 at 5:09
Re performance: The only bottleneck I experienced with my Python conversion script is reading the ASCII data, which I try to reduce by using NumPy's loadtxt or fromfile built-ins. From a 10MB ASCII file: reading NumPy data 7.672 sec; writing GDAL raster 0.203 sec; running ArcGIS's gp.CalculateStatistics_management 10.703 sec. So the GDAL Python API is really fast. Writing a C/C++ program will ultimately fix the ASCII reading bottleneck. – Mike Toews May 25 '11 at 5:48

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.