[gdal-dev] [python-ogr] Transform datasource problems

[Even Rouault]

Le jeudi 10 février 2011 09:09:45, Frank Broniewski a écrit :
> Even,
> 
> thank you for your verbose and helpful answer!
> 
> Am 09.02.2011 20:24, schrieb Even Rouault:
> > Several points :
> > 
> > 1) The implementation of Layer.GetFeature() is only efficient for just a
> > few drivers (shapefile). Other drivers will need to sequentially read
> > from the first feature ... So use GetNextFeature()
> 
> Are these informations somewhere documented? I had a look at the docs
> here [1],[2], but didn't find anything related. What about iterating
> over layers in a datasource, or fields in a feature? Are there similar
> concerns?

Well, the doc for OGRLayer::GetFeature() at 
http://gdal.org/ogr/classOGRLayer.html#cb7625383f161e5a04aeea2173dce411 
mentions "Use OGRLayer::TestCapability(OLCRandomRead) to establish if this 
layer supports efficient random access reading via GetFeature(); however, the 
call should always work if the feature exists as a fallback implementation 
just scans all the features in the layer looking for the desired feature."

DataSource.GetLayer() and Feature.GetFieldAsXXXX() are efficents since the 
datasource maintains a layer array and the feature maintains a field array

> 
> > 2) The real bug is that you don't rewrite the feature to its layer. So
> > you need to add a target_layer.SetFeature(target_feature)
> 
> Ah, Ok, so the link/object relation between the layer and the feature is
> lost, once the feature is fetched?

Exactly. Or to be more precise, any update done to the feature is only "in 
memory" until the feature is commited back to the layer with 
Layer.SetFeature().

> Since most things are done by
> reference in Python, this was not obvious to me. Thanks for pointing
> this out!

The GDAL Python bindings are known not to be very python'ish. They tend to 
make the C/C++ API accessible to Python in the most straighforward way. The 
main reason is that we use swig as a facility to have a common API between 
Python, Perl, Java and C#.

You should have a look at http://trac.osgeo.org/gdal/wiki/GdalOgrInPython and 
particularly the "Gotchas" section.

> 
> > 3) But what you are doing is really what ogr2ogr.py does, although your
> > approach will be less efficient as you first copy the layer and then
> > rewrite each of its feature.
> 
> Yes, I know. And I had a deep look at the ogr2ogr.py script before
> writing to the list, but frankly this is too much for me. It looks more
> like a c/c++ script than anything else to me - no functions, no classes,
> nothing from the python zen ..

Yes this is "by design" as explained at the beginning of the file. Basically 
you can open ogr2ogr.cpp and ogr2ogr.py side by side and see how the C/C++ API 
maps to Python.

> I really got lost in the middle of the
> script :-) but that's just a cause of my insufficient programming skills.

Yes, ogr2ogr code has become more complicated over the time, with the addition 
of new options. Hoewever, what you want to do should be doable in ~ 30 lines I 
think.

Or, just blindly use ogr2ogr.py and call its main() with the appropriate 
arguments ;-)

> 
> > For best performance, you should rather :
> > a) create a target layer from scratch
> > b) copy the feature definition from the source layer to the target layer
> > c) iterate over the source features, for each source feature, create a
> > target feature (Feature.SetFrom(other_feature)). reproject the geometry,
> > assign it to the target feature, and Layer.CreateFeature(target_feature)
> > 
> > Basically a stripped version of ogr2ogr.py
> 
> I will make a new version and post it here on the list, so everyone can
> have a look and maybe learn something from it
> 
> 
> Frank
> 
> 
> [1] http://www.gdal.org/ogr/ogr_formats.html
> [2] http://www.gdal.org/ogr/classOGRLayer.html