zonal_stats.org 13.3 KB

Construction d’un data frame avec les stats zonales

Init

Petit bout de python vide pour initialiser la session. Ca permet d’avoir la version de python utilisée

import geopandas as gp
import pandas as pd

Simplification du fichier TERLAB

  • Eliminer les géométries qui ne contiennent pas d’information de rendement.

Génération des fichiers de stats zonales

Le code est ici /work/OT/theia/oso/jordi/src/maassp/jobs/zonal_stats.pbs

Comme l’image pour les stats zonales contient 720 bandes, le fichier vecteur contient > 720 * 4 attributs, ce qui dépasse le nb de colonnes dans un fichier sqlite. Le shapefile n’a pas ce problème, mais est limité à 4GB, ce qui est trop petit.

J’ai ensuite séparé la série en 2 demi-séries pour pouvoir générer 2 shapefiles, mais on a ensuite le pb de faire le lien entre les polygones des 2 shapefiles générés, parce qu’on perd l’identifiant.

On pourrait commencer par ajouter un identifiant unique à chaque polygone du fichier TERLAB de façon à pouvoir faire le lien entre les différents fichiers de stats locales (reflectances, DEM, masques) générés à partir du même fichier TERLAB.

export ${iota2dir}=
export terlabdir='/work/OT/theia/oso/shapes/TERLAB/'
export terlab31=${terlabdir}'SURFACES-2017-PARCELLES-GRAPHIQUES-CONSTATEES_031_20180210'
cp ${terlab31}.shp ${terlabdir}'terlab31.shp'
cp ${terlab31}.dbf ${terlabdir}'terlab31.dbf'
cp ${terlab31}.shx ${terlabdir}'terlab31.shx'
cp ${terlab31}.prj ${terlabdir}'terlab31.prj'
python ${iota2dir}/scripts/VectorTools/AddFieldID.py ${terlabdir}'terlab31.shp'

Import des fichiers dans des DF

Le fichier TERLAB qui contient le dept 31 est:

terlabdir = '/work/OT/theia/oso/shapes/TERLAB/'
terlab31 = terlabdir+'terlab31.shp'
df_terlab = gp.read_file(terlab31)
print(df_terlab.head())

Les fichiers shape avec les stats sont

zs_dir='/work/OT/theia/oso/jordi/maassp_data'
zs_rf1=zs_dir+'/zonal_stats_refl_1.shp'
zs_rf2=zs_dir+'/zonal_stats_refl_2.shp'
zs_dem=zs_dir+'/zonal_stats_dem.shp'
zs_mask=zs_dir+'/zonal_stats_mask.shp'

Lire les fichiers dans un DF pandas

df_rf1 = gp.read_file(zs_rf1)
df_rf2 = gp.read_file(zs_rf2)
df_dem = gp.read_file(zs_dem).to_crs
df_mask = gp.read_file(zs_mask)
print(len(df_rf1))
print(len(df_rf2))
print(len(df_dem))
print(len(df_mask))

Il paraît que le fichier de DEM a perdu quelques polygones?

Vérification des projections :

print(df_dem.crs)
print(df_mask.crs)
print(df_rf1.crs)
print(df_rf1.crs)

Reprojection du DEM

df_dem = df_dem.to_crs(df_mask.crs)
print(df_dem.crs)

Fusion dans un seul DF

for c in df_terlab.columns:
    print(c)

for c in df_mask.columns:
    print(c)

print(df_dem.geometry.head())
print(df_mask.geometry.head())
print(df_rf1.geometry.head())
print(df_rf2.geometry.head())

On fait un spatial join.

df_dem_mask = gp.sjoin(df_mask, df_mask, op='intersects')
print(df_dem_mask.head())
print(len(df_dem_mask))

On a beaucoup trop d’éléments dans le nouveau DF, ce qui veut dire qu’il y a trop d’intersections. Si on le fait entre les 2 fichiers de réflectances, qui sont issus de la même donnée :

df_dem_mask = gp.sjoin(df_rf2, df_rf1, op='intersects')
print(df_dem_mask.head())
print(len(df_dem_mask))

Essayons de merger sur plusieurs colonnes

df_tmp = df_mask
df_dem_mask = df_tmp.merge(df_rf1, on=['CODE_CULTU', 'SURF_ADM'])
print(df_dem_mask.head())
print(len(df_dem_mask))

Templates


import matplotlib.pyplot as plt
# plt.figure(figsize=(10,5))
# df31.hist(column='RENDNORME')
# plt.savefig(matplot_lib_filename)
matplot_lib_filename