Skip to Content

par_covar.py

  1. import fiona
  2. import numpy as np
  3. from shapely.geometry import shape
  4. # extraction des point x,y du fichier shapefile
  5. x, y = zip(*[(shape(point['geometry']).x, shape(point['geometry']).y) for point in fiona.open('points.shp'])
  6. # calcul de la matrice de covariance
  7. cov = np.cov(x,y) # calcul basé sur (n-1)
  8. print cov
  9. [[ 12208.04067225   2314.67141452]
  10.  [  2314.67141452   2615.4249297 ]]
  11. # extraction des valeurs et  vecteurs propres
  12. eigenval,eigenvec = np.linalg.eig(cov)
  13. # rayons de l'ellipse (racine carrée des valeurs proprez
  14. sigmay,sigmax = np.sqrt(eigenval)
  15. print sigmax, sigmay
  16. 45.67393805, 112.85989981
  17. # angle de l'ellipse à partir de la matrice de covariance
  18. theta = np.arctan((2*cov[:,0][1])/(cov[:,0][0] - cov[:,1][1]))/2
  19. # changement d'axe pour le N géographique 
  20. theta = 90 - np.rad2deg(theta)
  21. print theta
  22. 77.1191513516
  23. # angle de l'ellipse à partir des vecteurs propres 
  24. theta = 90 - np.rad2deg(np.arccos(eigenvec[0, 0]))
  25. print theta
  26. 77.1191513516
  27. # aire de l'ellipse
  28. aire = np.pi * sigmax*sigmay
  29. print aire
  30. 16194.143806210201
  31. # excentricité de l'ellipse
  32. exc = np.sqrt(1 - ((min(sigmax,sigmay)**2)/(max(sigmax,sigmay)**2)))
  33. print exc
  34. 0.91445132916786864
  • 1184 lectures


Theme by Dr. Radut.