Le portail francophone de la géomatique

pourquoi ne pas avoir choisi une méthode d'interpolation classique, par exemple de distance inverse (inverse distance weighted)?

Bonjour,

Réponse courte, sur ArcGIS: http://resources.arcgis.com/fr/help/mai … 0077000000
L'outil propose deux méthodes d'interpolation: linéaire (de base, qui correspondrait sauf erreur à la réponse de Renaud) ou "voisins naturels", qui aura l'intérêt de donner un poids aux surfaces de vos  "triangles 3D". Explication ESRI ici: http://help.arcgis.com/fr/arcgisdesktop … 003r000000

Mais voir aussi la réponse de Baobazz! En repartant de votre semis de points, vous pouvez créer directement un MNT raster, avec pour ça plusieurs méthodes d'interpolation possibles. Ce serait j'imagine l'option la plus naturelle si vous travaillez avec du libre (GRASS GDAL, R...). Du coup, IDW est la plus classique, mais le champ des possibles est large (krigeage...) et il n'y a pas forcément UNE méthode mathématique qui s'impose aux autres comme la bonne. Voir ici pour une discussion similaire: http://gis.stackexchange.com/questions/ … -csv-files.

Bon courage!

Pour ceux qui ont des souvenirs de géométrie vectorielle:
Point donné I, triangle ABC
On calcule les produits vectoriels suivants (dans le plan, c'est à dire x1y2-x2y1)
V1 = IB^IC
V2 = IC^IA
V3 = IA^IB
D = V1+V2+V3
Si D est nul, le triangle est dégénéré (en fait, D est le double de la surface du triangle), sinon :
ZI = (V1*ZA+V2*ZB+V3*ZC)/D



En prime, le code Pascal de la fonction correspondante dans le logiciel Cumulus que j'ai développé pour le traitement de scans laser :

Function TDelaunay.InterpoleZ(var P : TPointDelaunay) : boolean;
var
   t : integer;
   A,B,C,D : extended;
begin
   T := QuelTriangle(P);
   if T < 0 then begin
      InterpoleZ := false;
      P.Z := 0
   end
   else with triangle[T],ListePoint do begin
      InterpoleZ :=  true;
      A := PVect(P,Sommet[TrSommet[1]],Sommet[TrSommet[2]]);
      B := PVect(P,Sommet[TrSommet[2]],Sommet[TrSommet[0]]);
      C := PVect(P,Sommet[TrSommet[0]],Sommet[TrSommet[1]]);
      D := A+B+C;
      if D <> 0 then
         P.Z := (A*Sommet[TrSommet[0]].Z+
                        B*Sommet[TrSommet[1]].Z+
                        C*Sommet[TrSommet[2]].Z)/D
      else P.Z := -1;
   end;
end;

Dernière modification par Yves Egels (Fri 14 June 2013 14:16)

Mais si, c'est indiqué au début de ma réponse, ce sont des produits vectoriels en dimension 2.
https://fr.wikipedia.org/wiki/Produit_vectoriel
J'aurai dû écrire en Latex, mais c'est trop long...
en développant complètement :
V1= IB^IC se calcule le la façon suivante :

V1 = (XB-XI)(YC-YI) - (XC-XI)(YB-YI) idem mutatis mutandis pour V2 et V3