[gdal-dev] Gdal-Grid lidar.

Nicolas Cadieux nicolas.cadieux at archeotec.ca
Fri Jun 5 19:00:40 PDT 2015


Thanks Even,

My programing skills are very limited (to mapbasic and thankfully  python).   Restructuring the code is beyond me but I will try to change the variables and recompile.  What to you recommend the new value be? 

I have never done recompiling  in anything but MapBasic so any link toward a tutorial to do this on Windows would be helpful.

I will also look at the rasterize.  This could be a good solution as I am already using the fill no data for anything larger than my 3 meter radius.  

Talking about null data, is it possible to write the no data value in gdal_grid?  I know I can modify it to -32768 for example, but when I open the layer, that no data value become a true value and I need to use another module (translate or warp i think) to explicitly write the no data value in the raster.

The learning curve will be steep but it will be fun!

Thanks for the help.


Nicolas Cadieux M.Sc.
Les Entreprises Archéotec inc. 
8548, rue Saint-Denis Montréal H2P 2H2
Téléphone: 514.381.5112  Fax: 514.381.4995
www.archeotec.ca

On Jun 5, 2015 7:13 PM, Even Rouault <even.rouault at spatialys.com> wrote:
>
> Le samedi 06 juin 2015 00:18:45, Nicolas Cadieux a écrit : 
> > Hi, 
> > I have been using gdal_grid to interpolate LiDAR data for some time now. 
> > This would be a typical command: 
> > 
> > gdal_grid -l P0344_sol2004_XYZI_MTM8_CGVD28_150529 -txe 348733.0 401135.0 
> > -tye 5568587.0 5773796.0 -a 
> > nearest:radius1=3.0:radius2=3.0:angle=0.0:nodata=0.0 -outsize 52402.0 
> > 205209.0 -of GTiff 
> > H:/Nicolas/pytemp/OutPut\P0344_sol2004_XYZI_MTM8_CGVD28_150529.vrt 
> > H:/Nicolas/pytemp/OutPut\P0344_sol2004_XYZI_MTM8_CGVD28_150529.tif --confi 
> > g GDAL_NUM_THREADS ALL_CPUS 
> > 
> > As you can see, this is a huge raster and the csv. file (Pointed with a 
> > .vrt) contains over a billion points.  I am trying to do this in one shot 
> > to avoid boundary effects.  (I know I can split this .csv file into 
> > smaller bits and make this more efficient).   The computer I am using has 
> > 12 cores, 64BG of memory (smaller than the csv file) and a 1TB pcie solid 
> > state drive.  It should manage the task but it has been running for 5 days 
> > with less than 10% done).  I have access to a super computer with gdal if 
> > all else fails.    
> > 
> > This is my question.  The "--config GDAL_NUM_THREADS ALL_CPUS" or "--config 
> > GDAL_NUM_THREADS 12" does not seem to work with the nearest neighbor 
> > algorithm. (It seem to work with the -a count command)  If I use this on a 
> > smaller data set, I see no difference in speed or cpu usage between using 
> > the "ALL_CPUS" or not using this switch.  Is this a bug or simply a 
> > limitation of the nearest neighbor algorithm? 
>
> Nicolas, 
>
> The threading mechanism is generic for all algorithms . Add "--debug on" and 
> you should see traces like "GDAL_GRID: Using 12 threads". Now, a reason for it 
> not being efficient is that the work buffer used by gdalgrid (16 MB) is too small 
> w.r.t to nearest speed. Hum, actually looking more closely at the code, I see 
> that the quad tree to index the points will be rebuilt for each work buffer. 
> This is likely the main cause of the inefficiency and why you don't see 
> measurable differences between non threaded or threaded computations (since the 
> building of the quadtree is not multithreaded). If you compiled from source, 
> you could try changing the "const GUInt32 nDesiredBufferSize = 16*1024*1024;" 
> line in apps/gdal_grid.cpp to a much higher value. A more proper solution 
> would be to avoid rebuilding the quad tree for each work buffer. This is mainly 
> code restructuring. 
>
> An alternative to gdal_grid would be to use gdal_rasterize + gdal_fillnodata. 
> This is usually much much faster, so you could try to check if it gives the 
> results you expect. 
>
> You could also try scipy : 
> http://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.interpolate.griddata.html#scipy.interpolate.griddata 
> Not sure how it behaves with that number of points. 
>
> Even 
>
> > I don't want to move this 
> > to a super computer if it will not multi-thread properly.  
> > 
> > Thanks, 
> > Nicolas Cadieux 
> > 
> > (Gdal1.11.2 downloaded with OSGEO4W 64bit install with QGIS 2.8 on windows 
> > 7_64)  I use GDAL from the command line not from QGIS. 
> > 
> > 
> > Nicolas Cadieux M.Sc. 
> > Les Entreprises Archéotec inc. 
> > 8548, rue Saint-Denis Montréal H2P 2H2 
> > Téléphone: 514.381.5112  Fax: 514.381.4995 
> > www.archeotec.ca 
> > _______________________________________________ 
> > gdal-dev mailing list 
> > gdal-dev at lists.osgeo.org 
> > http://lists.osgeo.org/mailman/listinfo/gdal-dev 
>
> -- 
> Spatialys - Geospatial professional services 
> http://www.spatialys.com 


More information about the gdal-dev mailing list