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

[Frank Broniewski]

Ok, this is the script made up yesterday:
http://paste.pocoo.org/show/336317/

and an example script which is using it:
http://paste.pocoo.org/show/336318/

The convert_datasource funtion can convert an input datasource from one 
format to another, and in an optional step do a coordinate 
transformation from the input srs to another srs.

I did try to make the funtion as pythonic as I am able to do ;-), but I 
think there are still many improvements to make. Escpecially the 
transformation part should be factored out into another function I 
think. I will see what I can do over the weekend and post another, 
improved version next week.


Frank


Am 10.02.2011 20:24, schrieb 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
>


-- 
Frank BRONIEWSKI

METRICO s.à r.l.
géomètres
technologies d'information géographique
rue des Romains 36
L-5433 NIEDERDONVEN

tél.: +352 26 74 94 - 28
fax.: +352 26 74 94 99
http://www.metrico.lu