Gfoss -- Geographic Free and Open Source Software - Italian mailing list

Trovare extent del layer da console di python

_‹ Previous Topic Next Topic _›

Classic

List

Threaded

4 messages Options

Romina Di Paolo

Trovare extent del layer da console di python

Ciao a tutti,

sto provando a creare uno script nella console di qgis che come prima cosa estragga da ciascuna cartella lo shape chiamato "particella.shp" e poi esegua un controllo sull'estensione del layer per verificare se è all'interno della provincia che mi interessa oppure no, nel primo caso vorrei che il layer venisse aggiunto nella toc altrimenti scartato.

Ammesso che mi sia riuscita a spiegare vi allego parte dello script

import sys,os,string

import processing

folder = 'C:\\Dati_geografici\\Particelle catastali\\Catastale'

lista_folder = (os.listdir(folder))

for sub_folder in lista_folder:

new_folder = folder+'\\'+sub_folder

lista_shp = (os.listdir(new_folder))

for shp in lista_shp:

if string.find(shp,'particella.shp') >=0:

MANCA FUNZIONE PER CONTROLLARE L'ESTENSIONE DEL LAYER

qgis.utils.iface.addVectorLayer(shp,shp,"ogr")

cercando tra le API di qgis ho visto che esiste questa classe QgsVectorLayer che tra le varie cose calcola anche l'extent ma non so come richiamarla e soprattutto quale operatore di contronto usare per verificare che il mio layer cada dentro un determinato extent.

Altro problema il caricamento del layer nella toc mediante addVectorLayer mi restituisce questo errore "Il layer non è valido: Il layer D197__particella.shp non è valido e non può essere aggiunto alla mappa" che nonc apisco da cosa dipenda.

Grazie in anticipo!

_______________________________________________
[hidden email]
http://lists.gfoss.it/cgi-bin/mailman/listinfo/gfoss
Questa e' una lista di discussione pubblica aperta a tutti.
I messaggi di questa lista non hanno relazione diretta con le posizioni dell'Associazione GFOSS.it.
666+40 iscritti al 5.6.2014

Luigi Pirelli-2

Re: Trovare extent del layer da console di python

senza scomodare gdal e rimanendo in pyqgis

1) carichi il layer in memoria (senza visualizzarlo)

vl = QgsVectorLayer(shp, "mioshape", "ogr")

2) prendi l'extent del vettore

extent = vl.extent()

come vedrai extent e' un QgsRectangle =>

3) ho bisogno di un QgsGeometry per usarne gli operatori spaziali

geomExtent = QgsGeometry.fromRect(extent)

4) alla stessa cosa per generarti la geometria dell'extent che vuoi contollare

geomExtentDiControllo = QgsGeometry.fromRect(...un qgis rectangle... )

ti consiglio di vederti: http://qgis.org/api/classQgsGeometry.html
per costruirti il QgsREctangle di controllo

5) usi l'operatore di QgsGeometry che meglio ti si addice

isInside = geomExtent.within(geomExtentDiControllo)

oppure

isInside = geomExtent.intersects(geomExtentDiControllo)

if isInside:
QgsMapLayerRegistri.instance().addMapLayer(vl)

ciao Luigi Pirelli

2015-03-06 10:36 GMT+01:00 Romina Di Paolo <[hidden email]>:

> Ciao a tutti,
>
> sto provando a creare uno script nella console di qgis che come prima cosa
> estragga da ciascuna cartella lo shape chiamato "particella.shp" e poi
> esegua un controllo sull'estensione del layer per verificare se è
> all'interno della provincia che mi interessa oppure no, nel primo caso
> vorrei che il layer venisse aggiunto nella toc altrimenti scartato.
> Ammesso che mi sia riuscita a spiegare vi allego parte dello script
>
> import sys,os,string
> import processing
>
>
> folder = 'C:\\Dati_geografici\\Particelle catastali\\Catastale'
> lista_folder = (os.listdir(folder))
>
> for sub_folder in lista_folder:
> new_folder = folder+'\\'+sub_folder
> lista_shp = (os.listdir(new_folder))
>
> for shp in lista_shp:
> if string.find(shp,'particella.shp') >=0:
> MANCA FUNZIONE PER CONTROLLARE L'ESTENSIONE DEL LAYER
> qgis.utils.iface.addVectorLayer(shp,shp,"ogr")
>
>
> cercando tra le API di qgis ho visto che esiste questa classe QgsVectorLayer
> che tra le varie cose calcola anche l'extent ma non so come richiamarla e
> soprattutto quale operatore di contronto usare per verificare che il mio
> layer cada dentro un determinato extent.
> Altro problema il caricamento del layer nella toc mediante addVectorLayer mi
> restituisce questo errore "Il layer non è valido: Il layer
> D197__particella.shp non è valido e non può essere aggiunto alla mappa" che
> nonc apisco da cosa dipenda.
>
> Grazie in anticipo!
>
> _______________________________________________
> [hidden email]
> http://lists.gfoss.it/cgi-bin/mailman/listinfo/gfoss
> Questa e' una lista di discussione pubblica aperta a tutti.
> I messaggi di questa lista non hanno relazione diretta con le posizioni
> dell'Associazione GFOSS.it.
> 666+40 iscritti al 5.6.2014

Romina Di Paolo

Re: Trovare extent del layer da console di python

Grazie Luigi per l'esempio fornito,

Quindi da quello che hai riportato nell'esempio da una parte ho la mia lista di shp (extent = vl.extent) che deve essere confrontata con un extent di controllo ovvero geomExtentDiControllo. Il mio problema è che non so come accedere alle proprietà geometriche di quest'ultimo.

Ho provato così

folder = 'S:\\Dati_geografici_RP\\I.6_Particelle catastali\\Catastale'
lista_folder = (os.listdir(folder))
geomExtentDiControllo = QgsGeometry.fromRect(QgsRectangle (xmin =442882.70,ymin = 5015200.28,xmax =487902.38,ymax =5087429.92))

for sub_folder in lista_folder:
new_folder = folder+'\\'+sub_folder
lista_shp = (os.listdir(new_folder))

for shp in lista_shp:
if string.find(shp,'PARTICELLA.shp') >=0:
vl = QgsVectorLayer(shp, "mioshape", "ogr")
extent = vl.extent()
geomExtent = QgsGeometry.fromRect(extent)
isInside = geomExtent.within(geomExtentDiControllo)
if isInside:
QgsMapLayerRegistri.instance().addMapLayer(vl)

In grassetto la parte che mi è ancora lacunosa infatti il programma finisce senza restituire nulla!

Sto usando come riferimento questo perchè sono ancora alle prime armi come si può capire!

http://www.padido.eu/gfoss/qgis/workshop/workshop/esplorare/python_in_qgis_tutorial2.html

2015-03-06 11:10 GMT+01:00 Luigi Pirelli <[hidden email]>:

Luigi Pirelli-2

Re: Trovare extent del layer da console di python

non avendo i dati sottomano, su due piedi farei il debug della nonna
1) controllerei che le x e y di geomExtentDiControllo siano lon e lat
2) controllerei che gli shape siano nella stessa proiezione con cui
hai espresso il bounding box di geomExtentDiControllo
3) controllerei che gli extent degli shape (fanne un print) per
verificare a mano se within debba darti true o false
4) verificherei anche con l'operatore intersects di QgsGeometry

ciao Luigi Pirelli

2015-03-09 12:44 GMT+01:00 Romina Di Paolo <[hidden email]>:

>
> Grazie Luigi per l'esempio fornito,
> Quindi da quello che hai riportato nell'esempio da una parte ho la mia lista
> di shp (extent = vl.extent) che deve essere confrontata con un extent di
> controllo ovvero geomExtentDiControllo. Il mio problema è che non so come
> accedere alle proprietà geometriche di quest'ultimo.
>
> Ho provato così
>
> folder = 'S:\\Dati_geografici_RP\\I.6_Particelle catastali\\Catastale'
> lista_folder = (os.listdir(folder))
> geomExtentDiControllo = QgsGeometry.fromRect(QgsRectangle (xmin
> =442882.70,ymin = 5015200.28,xmax =487902.38,ymax =5087429.92))
>
> for sub_folder in lista_folder:
> new_folder = folder+'\\'+sub_folder
> lista_shp = (os.listdir(new_folder))
>
> for shp in lista_shp:
> if string.find(shp,'PARTICELLA.shp') >=0:
> vl = QgsVectorLayer(shp, "mioshape", "ogr")
> extent = vl.extent()
> geomExtent = QgsGeometry.fromRect(extent)
> isInside = geomExtent.within(geomExtentDiControllo)
> if isInside:
> QgsMapLayerRegistri.instance().addMapLayer(vl)
>
> In grassetto la parte che mi è ancora lacunosa infatti il programma finisce
> senza restituire nulla!
>
> Sto usando come riferimento questo perchè sono ancora alle prime armi come
> si può capire!
> http://www.padido.eu/gfoss/qgis/workshop/workshop/esplorare/python_in_qgis_tutorial2.html
>
>
>
>
> 2015-03-06 11:10 GMT+01:00 Luigi Pirelli <[hidden email]>:
>
> senza scomodare gdal e rimanendo in pyqgis
>
> 1) carichi il layer in memoria (senza visualizzarlo)
>
> vl = QgsVectorLayer(shp, "mioshape", "ogr")
>
> 2) prendi l'extent del vettore
>
> extent = vl.extent()
>
> come vedrai extent e' un QgsRectangle =>
>
> 3) ho bisogno di un QgsGeometry per usarne gli operatori spaziali
>
> geomExtent = QgsGeometry.fromRect(extent)
>
> 4) alla stessa cosa per generarti la geometria dell'extent che vuoi
> contollare
>
> geomExtentDiControllo = QgsGeometry.fromRect(...un qgis rectangle... )
>
> ti consiglio di vederti: http://qgis.org/api/classQgsGeometry.html
> per costruirti il QgsREctangle di controllo
>
> 5) usi l'operatore di QgsGeometry che meglio ti si addice
>
> isInside = geomExtent.within(geomExtentDiControllo)
>
> oppure
>
> isInside = geomExtent.intersects(geomExtentDiControllo)
>
> if isInside:
> QgsMapLayerRegistri.instance().addMapLayer(vl)
>
> ciao Luigi Pirelli