Télécharger la plaquette Cast3M


fr en

Télécharger ippid.eso

Retour à la liste

Numérotation des lignes : 

ippid

  1. C IPPID     SOURCE    BP208322  21/01/28    21:15:14     10868          
  2.         SUBROUTINE IPPID(INUA,NBC,MTR1,X1,Y1,N,INDJ,NN,MM,NNMM,XP,XELIM)
  3. C-----------------------------------------------------------------------
  4. C NOM         : IPPID
  5. C DESCRIPTION : Interpolation dans un NUAGE par Ponderation Inverse a la
  6. C               Distance
  7. C LANGAGE     : ESOPE
  8. C AUTEUR      : Francois DI PAOLA
  9. C-----------------------------------------------------------------------
  10. C APPELE PAR  : IPLNU2
  11. C APPELE      :
  12. C-----------------------------------------------------------------------
  13. C ENTREES     : INUA (Objet de type NUAGE)
  14. C               NBC nombe de points (cad de NNMM-uplets) du nuage
  15. C               NNMM nombe de composantes du nuage
  16. C               NN nombre de composantes connues du point cible
  17. C               MM nombre de composantes a calculer du point cible
  18. C               MTR1.IADD contient la correspondance entre les numeros
  19. C               des composantes du point cible et celles du nuage :
  20. C               - La composante connue I (entre 1 et NN) est a la
  21. C                 position IADD(I) dans le nuage
  22. C               - La composante a calculer J (entre 1 et MM) est a la
  23. C                 position IADD(J) dans le nuage
  24. C               X1 tableau contenant les valeurs des NN composantes
  25. C               connues au point cible
  26. C               XP est l'exposant choisit pour la ponderation
  27. C SORTIES     : Y1 tableau contenant les valeurs des MM composantes
  28. C               interpolees au point cible
  29. C-----------------------------------------------------------------------
  30. C VERSION    : v1, 27/05/2016, version initiale
  31. C HISTORIQUE : v1, 27/05/2016, creation
  32. C HISTORIQUE :
  33. C HISTORIQUE :
  34. C-----------------------------------------------------------------------
  35. C Priere de PRENDRE LE TEMPS de completer les commentaires
  36. C en cas de modification de ce sous-programme afin de faciliter
  37. C la maintenance !
  38. C-----------------------------------------------------------------------
  39. C REMARQUES : - L'interpolation est exacte, c'est-a-dire que l'on
  40. C               retrouve les valeurs du nuage si l'on interpole en point
  41. C               du nuage
  42. C             - Pour plus d'information consultez :
  43. C               * Shepard, Donald (1968). « A two-dimensional interpolation
  44. C                 function for irregularly-spaced data ». Proceedings of
  45. C                 the 1968 ACM National Conference: 517–524.
  46. C               https://fr.wikipedia.org/wiki/Pond%C3%A9ration_inverse_%C3%A0_la_distance
  47. C               https://en.wikipedia.org/wiki/Inverse_distance_weighting
  48. C-----------------------------------------------------------------------
  49. C
  50. C
  51.       IMPLICIT INTEGER(I-N)
  52.       IMPLICIT    REAL*8(A-H,O-Z)
  53. -INC CCREEL
  54.  
  55. -INC PPARAM
  56. -INC CCOPTIO
  57. -INC SMNUAGE
  58.       SEGMENT MTR1
  59.         INTEGER IADD(NNMM)
  60.       ENDSEGMENT
  61.       REAL*8  X1,Y1
  62.       DIMENSION X1(N,*)
  63.       DIMENSION Y1(N,*)
  64.       REAL*8  DIS(NBC),W(NBC)
  65. C
  66. C     Critere pour le calcul des poids avec des distance trop grandes
  67.       XGRANDR=XGRAND**(1.D0/ABS(XP))
  68. C
  69. C     Calcul du nombre de points dans le nuage
  70.       MNUAG1=INUA
  71. C
  72. C     Si nuage vide, operation impossible
  73.       IF (NBC.EQ.0) THEN
  74.         CALL ERREUR(680)
  75.         RETURN
  76.       ENDIF
  77. C
  78. C     Calcul preliminaire des distances euclidiennes des points du nuage
  79. C     au point XI (sur les NN composantes connues)
  80.       NZER=0
  81. C     Boucle sur les points du nuage
  82.       DO I=1,NBC
  83.         DI=XZERO
  84. C       Boucle sur les NN composantes connues du champ et calcul de la
  85. C       distance euclidienne du point I
  86.         DO J=1,NN
  87.           X1J=X1(INDJ,J)
  88.           JJ=IADD(J)
  89.           NUAVF1=MNUAG1.NUAPOI(JJ)
  90.           XNJ=NUAVF1.NUAFLO(I)
  91.           DXJ=X1J-XNJ
  92.           DI=DI+(DXJ*DXJ)
  93.         ENDDO
  94.         DI=SQRT(DI)
  95. C       Si le point I du nuage est a une distance nulle, on arrete le
  96. C       calcul des distances car il pese pour la totalite de
  97. C       l'interpolation, ce qui permet d'avoir une interpolation exacte 
  98. c       -> facile : on fait le travail ici !
  99.         IF (ABS(DI).LT.XELIM) THEN
  100. c         boucle sur les MM composantes a traiter
  101.           DO J=1,MM 
  102.             II=IADD(J+NN)
  103.             NUAVF1=MNUAG1.NUAPOI(II)
  104.             Y1(INDJ,J)=NUAVF1.NUAFLO(I)
  105.           ENDDO
  106.           RETURN
  107.         ENDIF
  108. c       si le point n'est pas a une distance nulle on stocke la distance
  109.         DIS(I)=DI
  110.       ENDDO
  111.  
  112. C     Calcul des ponderations associees a chaque point du nuage
  113. C     (inverses a la distance)
  114. C     Si un point du nuage est a une distance nulle, 
  115.       WSOM=XZERO
  116.       DO I=1,NBC
  117. C       Si la distance est trop elevee, erreur division par zero
  118.         IF (DIS(I).GT.XGRANDR) THEN
  119.           CALL ERREUR(835)
  120.           RETURN
  121.         ENDIF
  122. C       On peut diviser plus sereinement
  123.         W(I)=(1.D0/DIS(I))**(XP)
  124.         WSOM=WSOM+W(I)
  125.       ENDDO
  126. C
  127. C     Calcul des valeurs des MM composantes inconnues au point X1
  128.       DO J=1,MM
  129.         FJ=XZERO
  130.         II=IADD(J+NN)
  131.         NUAVF1=MNUAG1.NUAPOI(II)
  132. C       Boucle sur les points du nuage
  133.         DO I=1,NBC
  134.           WI=W(I)
  135.           FI=NUAVF1.NUAFLO(I)
  136.           FJ=FJ+(WI*FI)
  137.         ENDDO
  138.         FJ=FJ/WSOM
  139.         Y1(INDJ,J)=FJ
  140.       ENDDO
  141.       RETURN
  142.       END
  143.  
  144.  
  145.  
  146.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales