Télécharger la plaquette Cast3M


fr en

Télécharger deplac.eso

Retour à la liste

Numérotation des lignes : 

deplac

  1. C DEPLAC    SOURCE    SP204843  24/03/15    21:15:03     11871          
  2.  
  3. C  CE SOUS-PROGRAMME A POUR PROPOS DE DEPLACER L'ENSEMBLE DES POINTS
  4. C  CONTENUS DANS UN OBJET
  5. C  09/2003 : Modifications dans le cas IDIM=1.
  6.  
  7.       SUBROUTINE DEPLAC
  8.  
  9.       IMPLICIT INTEGER(I-N)
  10.       IMPLICIT REAL*8 (A-H,O-Z)
  11.  
  12. -INC CCREEL
  13.  
  14. -INC PPARAM
  15. -INC CCOPTIO
  16. -INC SMCOORD
  17. -INC SMELEME
  18. -INC SMCHPOI
  19.  
  20.       CHARACTER*4 MCLE(11)
  21.       CHARACTER*4 MCL2(3)
  22.       CHARACTER*4 MCL3(7)
  23.       CHARACTER*4 MCL4(2)
  24.       CHARACTER*4 MCL5(2)
  25.       DIMENSION XU(3),XV(3),XW(3),XP(3)
  26.       REAL*8 XXX,RAP
  27.       SEGMENT ICPR(nbpts)
  28.  
  29.       DATA NCLE / 11 /
  30.       DATA MCLE / 'PLUS','MOIN','TOUR','HOMO','AFFI','SYME','PROJ',
  31.      .            'COOR','MILI','BARS','DEDU' /
  32.       DATA MCL2 / 'POIN','DROI','PLAN' /
  33.       DATA MCL3 / 'PLAN','SPHE','CYLI','CONI','TORI','DROI','CERC' /
  34.       DATA MCL4 / 'DIRE','CONI' /
  35.       DATA MCL5 / 'CYLI','CART' /
  36.  
  37. C  Lecture de l'option de la directive DEPLACER
  38.       CALL LIRMOT(MCLE,NCLE,ICLE,1)
  39.       IF (IERR.NE.0) RETURN
  40. C  Cas IDIM=1, seules les options suivantes sont disponibles :
  41. C   'PLUS','MOIN','HOMO','SYME','MILI','DEDU'
  42.       IF (IDIM.EQ.1) THEN
  43.         IF ((ICLE.NE.1).AND.(ICLE.NE.2).AND.(ICLE.NE.4).AND.
  44.      .      (ICLE.NE.6).AND.(ICLE.NE.9).AND.(ICLE.NE.11)) THEN
  45.           MOTERR(1:4)=MCLE(ICLE)
  46.           INTERR(1)=IDIM
  47.           CALL ERREUR(971)
  48.           RETURN
  49.         ENDIF
  50.       ENDIF
  51.  
  52.       idimp1=IDIM+1
  53.       ITYPLU=2
  54.  
  55. C  Option BARS : Elements de BARSOUM
  56. C ---------------
  57. C  Deplacement des noeuds milieu autour de la pointe de fissure
  58.       IF (ICLE.EQ.10) THEN
  59.         CALL BARSOU
  60.         RETURN
  61.       ENDIF
  62.  
  63. C  Option DEDU :
  64. C ---------------
  65. C  Transformation affine definie par au plus (IDIM+1) points
  66.       IF (ICLE.EQ.11) THEN
  67.         CALL DEDU(1)
  68.         RETURN
  69.       ENDIF
  70.  
  71.       CALL LIROBJ('MAILLAGE',MELEME,0,IRETOU)
  72.       ITYPSA=MELEME
  73.       IF (IRETOU.NE.1) THEN
  74.         ITYPLU=1
  75.         CALL LIROBJ('POINT',IP,1,IRETOU)
  76.         IF (IERR.NE.0) RETURN
  77.         ITYPSA=IP
  78.         CALL CRELEM(IP)
  79.         MELEME=IP
  80.       ENDIF
  81.       SEGACT MCOORD*MOD
  82.  
  83. C  Option MILI :
  84. C ---------------
  85.       IF (ICLE.EQ.9) GOTO 1900
  86.  
  87. C  Remplissage du tableau ICPR des points a deplacer
  88.       SEGINI ICPR
  89. C  Est-ce bien utile apres un SEGINI ?
  90. c       DO i=1,ICPR(/1)
  91. c         ICPR(i)=0
  92. c       ENDDO
  93. c     petite numerotation locale dans ICPR
  94.       NICPR=0
  95.       SEGACT MELEME
  96.       IPT1=MELEME
  97.       DO io=1,MAX(1,LISOUS(/1))
  98.         IF (LISOUS(/1).NE.0) THEN
  99.           IPT1=LISOUS(io)
  100.           SEGACT IPT1
  101.         ENDIF
  102.         DO i=1,IPT1.NUM(/1)
  103.           DO j=1,IPT1.NUM(/2)
  104.             I1=IPT1.NUM(i,j)
  105.             if(ICPR(I1).eq.0) then
  106.               NICPR=NICPR+1
  107.               ICPR(I1)=NICPR
  108.             endif
  109.           ENDDO
  110.         ENDDO
  111.         IF (LISOUS(/1).NE.0) SEGDES IPT1
  112.       ENDDO
  113.       SEGDES MELEME
  114.  
  115. C  Lecture des donnees supplementaires (selon ICLE) et action
  116. C   ICLE=1 : PLUS    1 POINT (VECTEUR)
  117. C   ICLE=2 : MOIN    1 POINT (VECTEUR)
  118. C   ICLE=3 : TOUR    1 NOMBRE (ANGLE) 1 OU 2 POINTS (AXE 2 OU 3 D)
  119. C   ICLE=4 : HOMO    1 NOMBRE (RAPPORT) 1 POINT (CENTRE)
  120. C   ICLE=5 : AFFI    1 NOMBRE (RAPPORT) 1 POINT (INVARIANT)
  121. C                                       1 POINT(VECTEUR)
  122. C   ICLE=6 : SYME    1 MOT (POIN DROI OU PLAN) 1,2 OU 3 POINTS
  123. C   ICLE=7 : PROJ    1 MOT (PLAN SPHE CYLI CONI TORI DROI CERC)
  124. C                          DES DONNEES  EN RAPPORT
  125. C   ICLE=8 : COOR    1 MOT (CART CYLI)
  126.       GOTO (100,200,300,400,500,600,700,1800),ICLE
  127.  
  128. C  Option PLUS :
  129. C ---------------
  130.  100  XSENS=1.D0
  131.       GOTO 201
  132. C  Option MOIN :
  133. C ---------------
  134.  200  XSENS=-1.D0
  135.  201  CALL LIROBJ('POINT',IP,0,IRETOU)
  136.       IF (IERR.NE.0) GOTO 10000
  137.       IF (IRETOU.EQ.0) THEN
  138.         CALL LIROBJ('CHPOINT ',IPCH,1,IRETOU)
  139.         IF (IERR.NE.0) GOTO 10000
  140.         CALL DEPCHP(ICPR,IPCH,XSENS)
  141.         GOTO 1000
  142.       ENDIF
  143.       IREF=idimp1*(IP-1)
  144.       DO j=1,IDIM
  145.         XV(j)=XCOOR(IREF+j)*XSENS
  146.       ENDDO
  147.       DO i=1,ICPR(/1)
  148.         IF (ICPR(i).NE.0) THEN
  149.           IREF=idimp1*(i-1)
  150.           DO j=1,IDIM
  151.             XCOOR(IREF+j)=XCOOR(IREF+j)+XV(j)
  152.           ENDDO
  153.         ENDIF
  154.       ENDDO
  155.       GOTO 1000
  156.  
  157. C  Option TOUR :
  158. C ---------------
  159. c     lecture de l'angle
  160.  300  CALL LIRREE(ANG,0,IRETOU)
  161.       IF (IERR.NE.0) GOTO 10000
  162.       IPCH1=0
  163. C  -- bp (nouvelle option 06/2013) : TOUR CHPOint --
  164.       IF (IRETOU.EQ.0) THEN
  165.         CALL LIROBJ('CHPOINT ',IPCH1,1,IRET1)
  166.         IF (IERR.NE.0) GOTO 10000
  167. c       creation d'un MPOVAL coherent avec ICPR
  168. c       avec comme composantes COS et SIN
  169.         N=NICPR
  170.         NC=2
  171.         SEGINI,MPOVAL
  172.         DO I=1,N
  173.           VPOCHA(I,1)=1.D0
  174. c           VPOCHA(I,2)=0.D0 =sin(0)
  175.         ENDDO
  176. c       remplissage du MPOVAL depuis le MCHPO1
  177.         MCHPO1=IPCH1
  178.         SEGACT,MCHPO1
  179.         DO IP=1,MCHPO1.IPCHP(/1)
  180.           MSOUP1=MCHPO1.IPCHP(IP)
  181.           SEGACT,MSOUP1
  182.           IPT1=MSOUP1.IGEOC
  183.           SEGACT,IPT1
  184.           MPOVA1=MSOUP1.IPOVAL
  185.           SEGACT,MPOVA1
  186.           N1  = MPOVA1.VPOCHA(/1)
  187.           NC1 = MPOVA1.VPOCHA(/2)
  188.           IF (NC1.ne.1) THEN
  189.             call ERREUR(180)
  190.             GOTO 10000
  191.           ENDIF
  192.           DO 303 J1=1,N1
  193.             I  = ICPR(IPT1.NUM(1,J1))
  194.             IF(I.eq.0) GOTO 303
  195.             ANG=MPOVA1.VPOCHA(J1,1)/180.D0*XPI
  196.             VPOCHA(I,1)=COS(ANG)
  197.             VPOCHA(I,2)=SIN(ANG)
  198.  303      CONTINUE
  199.           SEGDES,MPOVA1,IPT1
  200.           SEGDES,MSOUP1
  201.         ENDDO
  202.         SEGDES,MCHPO1
  203.       ELSE
  204.         ANG=ANG/180.D0*XPI
  205.         CO=COS(ANG)
  206.         SI=SIN(ANG)
  207.       ENDIF
  208. c     lecture des points centre et axe + fabrication repere de rotation
  209.       CALL LIROBJ('POINT',IP1,1,IRETOU)
  210.       IF (IDIM.EQ.3) CALL LIROBJ('POINT',IP2,1,IRETOU)
  211.       IF (IERR.NE.0) GOTO 10000
  212.       IREF=idimp1*(IP1-1)
  213.       XP(1)=XCOOR(IREF+1)
  214.       XP(2)=XCOOR(IREF+2)
  215.       IF (IDIM.LT.3) THEN
  216.         XP(3)=0.D0
  217.         XW(1)=0.D0
  218.         XW(2)=0.D0
  219.         XW(3)=1.D0
  220.       ELSE
  221.         XP(3)=XCOOR(IREF+3)
  222.         IREF=idimp1*(IP2-1)
  223.         XW(1)=XCOOR(IREF+1)-XP(1)
  224.         XW(2)=XCOOR(IREF+2)-XP(2)
  225.         XW(3)=XCOOR(IREF+3)-XP(3)
  226.         XN=SQRT(XW(1)**2+XW(2)**2+XW(3)**2)
  227.         IF (XN.EQ.0.D0) CALL ERREUR(21)
  228.         IF (IERR.NE.0) GOTO 10000
  229.         XW(1)=XW(1)/XN
  230.         XW(2)=XW(2)/XN
  231.         XW(3)=XW(3)/XN
  232.       ENDIF
  233.       XU(1)=-XW(2)
  234.       XU(2)=XW(1)
  235.       XN=SQRT(XU(1)**2+XU(2)**2)
  236.       IF (XN.GE.0.1D0) THEN
  237.         XU(3)=0.D0
  238.         XU(1)=XU(1)/XN
  239.         XU(2)=XU(2)/XN
  240.       ELSE
  241.         XU(1)=0.D0
  242.         XU(2)=-XW(3)
  243.         XU(3)=XW(2)
  244.         XN=SQRT(XU(2)**2+XU(3)**2)
  245.         XU(2)=XU(2)/XN
  246.         XU(3)=XU(3)/XN
  247.       ENDIF
  248.       XV(1)=XW(2)*XU(3)-XW(3)*XU(2)
  249.       XV(2)=XW(3)*XU(1)-XW(1)*XU(3)
  250.       XV(3)=XW(1)*XU(2)-XW(2)*XU(1)
  251. c     boucle sur les noeuds + rotation
  252.       DO i=1,ICPR(/1)
  253.         IF (ICPR(i).NE.0) THEN
  254.           IREF=idimp1*(i-1)
  255.           XD=XCOOR(IREF+1)-XP(1)
  256.           YD=XCOOR(IREF+2)-XP(2)
  257.           ZD=0.D0
  258.           IF (IDIM.EQ.3) ZD=XCOOR(IREF+3)-XP(3)
  259.           XE=XD*XU(1)+YD*XU(2)+ZD*XU(3)
  260.           YE=XD*XV(1)+YD*XV(2)+ZD*XV(3)
  261.           ZE=XD*XW(1)+YD*XW(2)+ZD*XW(3)
  262.           IF(IPCH1.ne.0) THEN
  263.             CO = VPOCHA(ICPR(i),1)
  264.             SI = VPOCHA(ICPR(i),2)
  265.           ENDIF
  266.           XD=XE*CO-YE*SI
  267.           YD=XE*SI+YE*CO
  268.           ZD=ZE
  269.           XCOOR(IREF+1)=XD*XU(1)+YD*XV(1)+ZD*XW(1)+XP(1)
  270.           XCOOR(IREF+2)=XD*XU(2)+YD*XV(2)+ZD*XW(2)+XP(2)
  271.           IF (IDIM.EQ.3) XCOOR(IREF+3)=XD*XU(3)+YD*XV(3)+ZD*XW(3)+XP(3)
  272.         ENDIF
  273.       ENDDO
  274.       IF(IPCH1.ne.0) SEGSUP,MPOVAL
  275.       GOTO 1000
  276.  
  277. C  Option HOMO :
  278. C ---------------
  279.  400  CALL LIRREE(RAP,1,IRETOU)
  280.       CALL LIROBJ('POINT',IP,1,IRETOU)
  281.       IF (RAP.EQ.0.D0) CALL ERREUR(36)
  282.       IF (IERR.NE.0) GOTO 10000
  283.       IREF=idimp1*(IP-1)
  284.       DO j=1,IDIM
  285.         XP(j)=XCOOR(IREF+j)
  286.       ENDDO
  287.       DO i=1,ICPR(/1)
  288.         IF (ICPR(i).NE.0) THEN
  289.           IREF=idimp1*(i-1)
  290.           DO j=1,IDIM
  291.             XCOOR(IREF+j)=XP(j)+RAP*(XCOOR(IREF+j)-XP(j))
  292.           ENDDO
  293.         ENDIF
  294.       ENDDO
  295.       GOTO 1000
  296.  
  297. C  Option AFFI :
  298. C ---------------
  299.  500  CALL LIRREE(RAP,1,IRETOU)
  300.       CALL LIROBJ('POINT',IPC,1,IRETOU)
  301.       CALL LIROBJ('POINT',IPV,1,IRETOU)
  302.       IF (RAP.EQ.0.D0) CALL ERREUR(36)
  303.       IF (IERR.NE.0) GOTO 10000
  304.       RAP=RAP-1.
  305.       IREF=idimp1*(IPC-1)
  306.       XP(1)=XCOOR(IREF+1)
  307.       XP(2)=XCOOR(IREF+2)
  308.       XP(3)=0.D0
  309.       IF (IDIM.EQ.3) XP(3)=XCOOR(IREF+3)
  310.       IREF=idimp1*(IPV-1)
  311.       XV(1)=XCOOR(IREF+1)-XP(1)
  312.       XV(2)=XCOOR(IREF+2)-XP(2)
  313.       XV(3)=0.D0
  314.       IF (IDIM.EQ.3) XV(3)=XCOOR(IREF+3)-XP(3)
  315.       XN=XV(1)**2+XV(2)**2+XV(3)**2
  316.       IF (XN.EQ.0.D0) CALL ERREUR(21)
  317.       IF (IERR.NE.0) GOTO 10000
  318.       XN=SQRT(XN)
  319.       XV(1)=XV(1)/XN
  320.       XV(2)=XV(2)/XN
  321.       XV(3)=XV(3)/XN
  322.       DO i=1,ICPR(/1)
  323.         IF (ICPR(i).NE.0) THEN
  324.           IREF=idimp1*(i-1)
  325.           XU(1)=XCOOR(IREF+1)-XP(1)
  326.           XU(2)=XCOOR(IREF+2)-XP(2)
  327.           XU(3)=0.D0
  328.           IF (IDIM.EQ.3) XU(3)=XCOOR(IREF+3)-XP(3)
  329.           SCA=(XU(1)*XV(1)+XU(2)*XV(2)+XU(3)*XV(3))*RAP
  330.           DO j=1,IDIM
  331.             XCOOR(IREF+j)=XP(j)+XU(j)+SCA*XV(j)
  332.           ENDDO
  333.         ENDIF
  334.       ENDDO
  335.       GOTO 1000
  336.  
  337. C  Option SYME :
  338. C ---------------
  339.  600  CALL LIRMOT(MCL2,3,JCLE,1)
  340.       IF (IERR.NE.0) GOTO 10000
  341.       IF ((IDIM.EQ.1).AND.(JCLE.NE.1)) THEN
  342.         MOTERR(1:4)=MCL2(JCLE)
  343.         INTERR(1)=IDIM
  344.         CALL ERREUR(971)
  345.         GOTO 10000
  346.       ENDIF
  347.       GOTO (610,620,630),JCLE
  348.  
  349. C  Option SYME POIN :
  350. C --------------------
  351.  610  CALL LIROBJ('POINT',IP,1,IRETOU)
  352.       IF (IERR.NE.0) GOTO 10000
  353.       IREF=idimp1*(IP-1)
  354.       DO j=1,IDIM
  355.         XP(j)=2*XCOOR(IREF+j)
  356.       ENDDO
  357.       DO i=1,ICPR(/1)
  358.         IF (ICPR(i).NE.0) THEN
  359.           IREF=idimp1*(i-1)
  360.           DO j=1,IDIM
  361.             XCOOR(IREF+j)=XP(j)-XCOOR(IREF+j)
  362.           ENDDO
  363.         ENDIF
  364.       ENDDO
  365.       GOTO 1000
  366.  
  367. C  Option SYME DROI :
  368. C --------------------
  369.  620  CALL LIROBJ('POINT',IP1,1,IRETOU)
  370.       CALL LIROBJ('POINT',IP2,1,IRETOU)
  371.       IF (IERR.NE.0) GOTO 10000
  372.       IREF=idimp1*(IP1-1)
  373.       XP(1)=XCOOR(IREF+1)
  374.       XP(2)=XCOOR(IREF+2)
  375.       XP(3)=0.D0
  376.       IF (IDIM.EQ.3) XP(3)=XCOOR(IREF+3)
  377.       IREF=idimp1*(IP2-1)
  378.       XU(1)=XCOOR(IREF+1)-XP(1)
  379.       XU(2)=XCOOR(IREF+2)-XP(2)
  380.       XU(3)=0.D0
  381.       IF (IDIM.EQ.3) XU(3)=XCOOR(IREF+3)-XP(3)
  382.       XN=XU(1)**2+XU(2)**2+XU(3)**2
  383.       IF (XN.EQ.0.D0) CALL ERREUR(21)
  384.       IF (IERR.NE.0) GOTO 10000
  385.       XN=SQRT(XN)
  386.       XU(1)=XU(1)/XN
  387.       XU(2)=XU(2)/XN
  388.       XU(3)=XU(3)/XN
  389.       DO i=1,ICPR(/1)
  390.         IF (ICPR(i).NE.0) THEN
  391.           IREF=idimp1*(i-1)
  392.           XV(1)=XCOOR(IREF+1)-XP(1)
  393.           XV(2)=XCOOR(IREF+2)-XP(2)
  394.           XV(3)=0.D0
  395.           IF (IDIM.GE.3) XV(3)=XCOOR(IREF+3)-XP(3)
  396.           SCA=2*(XU(1)*XV(1)+XU(2)*XV(2)+XU(3)*XV(3))
  397.           DO j=1,IDIM
  398.             XCOOR(IREF+j)=XP(j)+SCA*XU(j)-XV(j)
  399.           ENDDO
  400.         ENDIF
  401.       ENDDO
  402.       GOTO 1000
  403.  
  404. C  Option SYME PLAN :
  405. C --------------------
  406.  630  IF (IDIM.EQ.2) CALL ERREUR(21)
  407.       IF (IERR.NE.0) GOTO 10000
  408.       CALL LIROBJ('POINT',IP1,1,IRETOU)
  409.       CALL LIROBJ('POINT',IP2,1,IRETOU)
  410.       CALL LIROBJ('POINT',IP3,1,IRETOU)
  411.       IF (IERR.NE.0) GOTO 10000
  412.       IREF=idimp1*(IP1-1)
  413.       XP(1)=XCOOR(IREF+1)
  414.       XP(2)=XCOOR(IREF+2)
  415.       XP(3)=XCOOR(IREF+3)
  416.       IREF=idimp1*(IP2-1)
  417.       XU(1)=XCOOR(IREF+1)-XP(1)
  418.       XU(2)=XCOOR(IREF+2)-XP(2)
  419.       XU(3)=XCOOR(IREF+3)-XP(3)
  420.       IREF=idimp1*(IP3-1)
  421.       XV(1)=XCOOR(IREF+1)-XP(1)
  422.       XV(2)=XCOOR(IREF+2)-XP(2)
  423.       XV(3)=XCOOR(IREF+3)-XP(3)
  424.       XW(1)=XU(2)*XV(3)-XU(3)*XV(2)
  425.       XW(2)=XU(3)*XV(1)-XU(1)*XV(3)
  426.       XW(3)=XU(1)*XV(2)-XU(2)*XV(1)
  427.       XN=XW(1)**2+XW(2)**2+XW(3)**2
  428.       IF (XN.EQ.0.D0) CALL ERREUR(21)
  429.       IF (IERR.NE.0) GOTO 10000
  430.       XN=SQRT(XN)
  431.       XW(1)=XW(1)/XN
  432.       XW(2)=XW(2)/XN
  433.       XW(3)=XW(3)/XN
  434.       DO i=1,ICPR(/1)
  435.         IF (ICPR(i).NE.0) THEN
  436.           IREF=idimp1*(i-1)
  437.           XU(1)=XCOOR(IREF+1)-XP(1)
  438.           XU(2)=XCOOR(IREF+2)-XP(2)
  439.           XU(3)=XCOOR(IREF+3)-XP(3)
  440.           SCA=2*(XU(1)*XW(1)+XU(2)*XW(2)+XU(3)*XW(3))
  441.           XCOOR(IREF+1)=XP(1)+XU(1)-SCA*XW(1)
  442.           XCOOR(IREF+2)=XP(2)+XU(2)-SCA*XW(2)
  443.           XCOOR(IREF+3)=XP(3)+XU(3)-SCA*XW(3)
  444.         ENDIF
  445.       ENDDO
  446.       GOTO 1000
  447.  
  448. C  Option SYME PROJ :
  449. C --------------------
  450.  700  CALL LIRMOT(MCL4,2,ITYYP,0)
  451. C  Lecture du vecteur de projection ou du centre (projection conique)
  452.       CALL LIROBJ('POINT',IP1,1,IRETOU)
  453.       IF (IERR.NE.0) GOTO 10000
  454.       IREF=idimp1*(IP1-1)
  455.       XVECS=XCOOR(IREF+1)
  456.       YVECS=XCOOR(IREF+2)
  457.       ZVECS=0.D0
  458.       IF (IDIM.GE.3) ZVECS=XCOOR(IREF+3)
  459.       ICARR=1
  460.       IF (ITYYP.NE.0) ICARR=0
  461.       CALL LIRMOT(MCL3,7,JCLE,ICARR)
  462.       IF (IERR.NE.0) GOTO 10000
  463.       IF (JCLE.EQ.0) THEN
  464.         JCLE=ITYYP+2
  465.         ITYYP=0
  466.       ENDIF
  467.       IF ((IDIM.EQ.2).AND.(JCLE.LE.5)) CALL ERREUR(34)
  468.       IF ((IDIM.EQ.3).AND.(JCLE.GE.6)) CALL ERREUR(34)
  469.       IF (IERR.NE.0) GOTO 10000
  470. C  Projection cylindrique normalisation du vecteur
  471.       IF (ITYYP.NE.2) THEN
  472.         SVECS=SQRT(XVECS*XVECS+YVECS*YVECS+ZVECS*ZVECS)
  473.         IF (SVECS.EQ.0.D0) CALL ERREUR(21)
  474.         IF (IERR.NE.0) GOTO 10000
  475.         XVEC=XVECS/SVECS
  476.         YVEC=YVECS/SVECS
  477.         ZVEC=ZVECS/SVECS
  478.       ENDIF
  479.       GOTO (710,720,730,740,750,760,770),JCLE
  480.  
  481. C  Option PROJ PLAN :
  482. C --------------------
  483.  710  CALL LIROBJ('POINT',IP1,1,IRETOU)
  484.       CALL LIROBJ('POINT',IP2,1,IRETOU)
  485.       CALL LIROBJ('POINT',IP3,1,IRETOU)
  486.       IF (IERR.NE.0) GOTO 10000
  487.       IREF=idimp1*(IP1-1)
  488.       XPT1=XCOOR(IREF+1)
  489.       YPT1=XCOOR(IREF+2)
  490.       ZPT1=XCOOR(IREF+3)
  491.       IREF=idimp1*(IP2-1)
  492.       XV2=XCOOR(IREF+1)-XPT1
  493.       YV2=XCOOR(IREF+2)-YPT1
  494.       ZV2=XCOOR(IREF+3)-ZPT1
  495.       IREF=idimp1*(IP3-1)
  496.       XV3=XCOOR(IREF+1)-XPT1
  497.       YV3=XCOOR(IREF+2)-YPT1
  498.       ZV3=XCOOR(IREF+3)-ZPT1
  499.       XV1=YV2*ZV3-ZV2*YV3
  500.       YV1=ZV2*XV3-XV2*ZV3
  501.       ZV1=XV2*YV3-YV2*XV3
  502.       SV1=SQRT(XV1*XV1+YV1*YV1+ZV1*ZV1)
  503.       IF (SV1.EQ.0.D0) CALL ERREUR(21)
  504.       IF (IERR.NE.0) GOTO 10000
  505.       XV1=XV1/SV1
  506.       YV1=YV1/SV1
  507.       ZV1=ZV1/SV1
  508.       GOTO 780
  509.  
  510. C  Option PROJ SPHE :
  511. C --------------------
  512.  720  CALL LIROBJ('POINT',IPCEN,1,IRETOU)
  513.       CALL LIROBJ('POINT',IP1,1,IRETOU)
  514.       IF (IERR.NE.0) GOTO 10000
  515.       IREF=(IPCEN-1)*idimp1
  516.       XPT1=XCOOR(IREF+1)
  517.       YPT1=XCOOR(IREF+2)
  518.       ZPT1=XCOOR(IREF+3)
  519.       IREF=(IP1-1)*idimp1
  520.       XV1=XCOOR(IREF+1)-XPT1
  521.       YV1=XCOOR(IREF+2)-YPT1
  522.       ZV1=XCOOR(IREF+3)-ZPT1
  523.       ANGLE=SQRT(XV1*XV1+YV1*YV1+ZV1*ZV1)
  524.       IF (ANGLE.EQ.0.D0) CALL ERREUR(21)
  525.       IF (IERR.NE.0) GOTO 10000
  526.       GOTO 780
  527.  
  528. C  Option PROJ CYLI :
  529. C --------------------
  530.  730  CALL LIROBJ('POINT',IPOIN1,1,IRETOU)
  531.       CALL LIROBJ('POINT',IPOIN2,1,IRETOU)
  532.       CALL LIROBJ('POINT',IP1,1,IRETOU)
  533.       IF (IERR.NE.0) GOTO 10000
  534.       IREF=(IPOIN1-1)*idimp1
  535.       XPT1=XCOOR(IREF+1)
  536.       YPT1=XCOOR(IREF+2)
  537.       ZPT1=XCOOR(IREF+3)
  538.       IREF=(IPOIN2-1)*idimp1
  539.       XV1=XCOOR(IREF+1)-XPT1
  540.       YV1=XCOOR(IREF+2)-YPT1
  541.       ZV1=XCOOR(IREF+3)-ZPT1
  542.       S=SQRT(XV1*XV1+YV1*YV1+ZV1*ZV1)
  543.       IF (S.EQ.0.D0) CALL ERREUR(21)
  544.       IF (IERR.NE.0) GOTO 10000
  545.       XV1=XV1/S
  546.       YV1=YV1/S
  547.       ZV1=ZV1/S
  548.       IREF=(IP1-1)*idimp1
  549.       XV2=XCOOR(IREF+1)-XPT1
  550.       YV2=XCOOR(IREF+2)-YPT1
  551.       ZV2=XCOOR(IREF+3)-ZPT1
  552.       XV3=YV1*ZV2-ZV1*YV2
  553.       YV3=ZV1*XV2-XV1*ZV2
  554.       ZV3=XV1*YV2-YV1*XV2
  555.       ANGLE=SQRT(XV3*XV3+YV3*YV3+ZV3*ZV3)
  556.       IF (ANGLE.EQ.0) CALL ERREUR(21)
  557.       IF (IERR.NE.0) RETURN
  558.       GOTO 780
  559.  
  560. C  Option PROJ CONI :
  561. C --------------------
  562.  740  CALL LIROBJ('POINT',IPT1,1,IRETOU)
  563.       CALL LIROBJ('POINT',IP1,1,IRETOU)
  564.       CALL LIROBJ('POINT',IP2,1,IRETOU)
  565.       IF (IERR.NE.0) GOTO 10000
  566.       IREF=(IPT1-1)*idimp1
  567.       XPT1=XCOOR(IREF+1)
  568.       YPT1=XCOOR(IREF+2)
  569.       ZPT1=XCOOR(IREF+3)
  570.       IREF=(IP1-1)*idimp1
  571.       XV1=XCOOR(IREF+1)-XPT1
  572.       YV1=XCOOR(IREF+2)-YPT1
  573.       ZV1=XCOOR(IREF+3)-ZPT1
  574.       SV1=SQRT(XV1*XV1+YV1*YV1+ZV1*ZV1)
  575.       IF (SV1.EQ.0.D0) CALL ERREUR(21)
  576.       IF (IERR.NE.0) GOTO 10000
  577.       XV1=XV1/SV1
  578.       YV1=YV1/SV1
  579.       ZV1=ZV1/SV1
  580.       IREF=(IP2-1)*idimp1
  581.       XV2=XCOOR(IREF+1)-XPT1
  582.       YV2=XCOOR(IREF+2)-YPT1
  583.       ZV2=XCOOR(IREF+3)-ZPT1
  584.       SV2=SQRT(XV2*XV2+YV2*YV2+ZV2*ZV2)
  585.       IF (SV2.EQ.0.D0) CALL ERREUR(21)
  586.       IF (IERR.NE.0) GOTO 10000
  587.       XV2=XV2/SV2
  588.       YV2=YV2/SV2
  589.       ZV2=ZV2/SV2
  590.       ANGLE=(XV1*XV2+YV1*YV2+ZV1*ZV2)**2
  591.       GOTO 780
  592.  
  593. C  Option PROJ TORI
  594. C ------------------
  595.  750  CALL LIROBJ('POINT',IPT1,1,IRETOU)
  596.       CALL LIROBJ('POINT',IPAX,1,IRETOU)
  597.       CALL LIROBJ('POINT',IPCS,1,IRETOU)
  598.       CALL LIROBJ('POINT',IP1 ,1,IRETOU)
  599.       IF (IERR.NE.0) GOTO 10000
  600.       IREF=(IPT1-1)*idimp1
  601.       XPT1=XCOOR(IREF+1)
  602.       YPT1=XCOOR(IREF+2)
  603.       ZPT1=XCOOR(IREF+3)
  604.       IREF=(IPAX-1)*idimp1
  605.       XV1=XCOOR(IREF+1)-XPT1
  606.       YV1=XCOOR(IREF+2)-YPT1
  607.       ZV1=XCOOR(IREF+3)-ZPT1
  608.       SV1=XV1*XV1+YV1*YV1+ZV1*ZV1
  609.       IF (SV1.EQ.0.D0) CALL ERREUR(21)
  610.       IF (IERR.NE.0) GOTO 10000
  611.       SV1=SQRT(SV1)
  612.       XV1=XV1/SV1
  613.       YV1=YV1/SV1
  614.       ZV1=ZV1/SV1
  615.       IREF=(IPCS-1)*idimp1
  616.       XV2=XCOOR(IREF+1)-XPT1
  617.       YV2=XCOOR(IREF+2)-YPT1
  618.       ZV2=XCOOR(IREF+3)-ZPT1
  619.       SCA=XV2*XV1+YV2*YV1+ZV2*ZV1
  620.       XPT1=XPT1+SCA*XV1
  621.       YPT1=YPT1+SCA*YV1
  622.       ZPT1=ZPT1+SCA*ZV1
  623.       XV2=XV2-SCA*XV1
  624.       YV2=YV2-SCA*YV1
  625.       ZV2=ZV2-SCA*ZV1
  626.       GR2=XV2*XV2+YV2*YV2+ZV2*ZV2
  627.       IF (GR2.EQ.0.D0) CALL ERREUR(21)
  628.       IF (IERR.NE.0) GOTO 10000
  629.       IREF=(IP1-1)*idimp1
  630.       XV2=XCOOR(IREF+1)-XPT1
  631.       YV2=XCOOR(IREF+2)-YPT1
  632.       ZV2=XCOOR(IREF+3)-ZPT1
  633.       SCA=XV2*XV1+YV2*YV1+ZV2*ZV1
  634.       XC=XV2-SCA*XV1
  635.       YC=YV2-SCA*YV1
  636.       ZC=ZV2-SCA*ZV1
  637.       SC=XC*XC+YC*YC+ZC*ZC
  638.       IF (SC.EQ.0.D0) CALL ERREUR(21)
  639.       IF (IERR.NE.0) GOTO 10000
  640.       RAP=SQRT(GR2/SC)
  641.       XC=XC*RAP
  642.       YC=YC*RAP
  643.       ZC=ZC*RAP
  644.       PR2=(XV2-XC)**2+(YV2-YC)**2+(ZV2-ZC)**2
  645.       IF (PR2.EQ.0.D0) CALL ERREUR(21)
  646.       IF (IERR.NE.0) GOTO 10000
  647.       XV2=PR2
  648.       YV2=GR2
  649.       GOTO 780
  650.  
  651. C  Option PROJ DROI :
  652. C --------------------
  653.  760  CALL LIROBJ('POINT',IP1,1,IRETOU)
  654.       CALL LIROBJ('POINT',IP2,1,IRETOU)
  655.       IF (IERR.NE.0) GOTO 10000
  656.       IREF=(IP1-1)*idimp1
  657.       XPT1=XCOOR(IREF+1)
  658.       YPT1=XCOOR(IREF+2)
  659.       ZPT1=0.D0
  660.       IREF=(IP2-1)*idimp1
  661.       XV2=XCOOR(IREF+1)-XPT1
  662.       YV2=XCOOR(IREF+2)-YPT1
  663.       ZV2=0.D0
  664.       XV3=0.D0
  665.       YV3=0.D0
  666.       ZV3=1.D0
  667.       XV1=YV2*ZV3-ZV2*YV3
  668.       YV1=ZV2*XV3-XV2*ZV3
  669.       ZV1=XV2*YV3-YV2*XV3
  670.       SV1=SQRT(XV1*XV1+YV1*YV1+ZV1*ZV1)
  671.       IF (SV1.EQ.0.D0) CALL ERREUR(21)
  672.       IF (IERR.NE.0) GOTO 10000
  673.       XV1=XV1/SV1
  674.       YV1=YV1/SV1
  675.       ZV1=ZV1/SV1
  676.       GOTO 780
  677.  
  678. C  Option PROJ CERC :
  679. C --------------------
  680.  770  CALL LIROBJ('POINT',IPCEN,1,IRETOU)
  681.       CALL LIROBJ('POINT',IP1,1,IRETOU)
  682.       IF (IERR.NE.0) GOTO 10000
  683.       IREF=(IPCEN-1)*idimp1
  684.       XPT1=XCOOR(IREF+1)
  685.       YPT1=XCOOR(IREF+2)
  686.       ZPT1=0.D0
  687.       IREF=(IP1-1)*idimp1
  688.       XV1=XCOOR(IREF+1)-XPT1
  689.       YV1=XCOOR(IREF+2)-YPT1
  690.       ZV1=0.D0
  691.       ANGLE=SQRT(XV1*XV1+YV1*YV1+ZV1*ZV1)
  692.       IF (ANGLE.EQ.0.D0) CALL ERREUR(21)
  693.       IF (IERR.NE.0) GOTO 10000
  694.       GOTO 780
  695. C  Traitement commun - Options PROJ
  696.  780  DO i=1,ICPR(/1)
  697.         IF (ICPR(i).NE.0) THEN
  698.           IREF=(i-1)*idimp1
  699.           XPOIN=XCOOR(IREF+1)
  700.           YPOIN=XCOOR(IREF+2)
  701.           ZPOIN=0.D0
  702.           IF (IDIM.GE.3) ZPOIN=XCOOR(IREF+3)
  703.           IF (ITYYP.EQ.2) THEN
  704.             XVEC=XVECS-XPOIN
  705.             YVEC=YVECS-YPOIN
  706.             ZVEC=ZVECS-ZPOIN
  707.             SVEC=SQRT(XVEC*XVEC+YVEC*YVEC+ZVEC*ZVEC)
  708.             IF (SVEC.EQ.0.D0) CALL ERREUR(21)
  709.             IF (IERR.NE.0) RETURN
  710.             XVEC=XVEC/SVEC
  711.             YVEC=YVEC/SVEC
  712.             ZVEC=ZVEC/SVEC
  713.           ENDIF
  714.           GOTO (810,820,830,840,850,860,870),JCLE
  715. C  Option PROJ PLAN
  716.  810      XV2=XPOIN-XPT1
  717.           YV2=YPOIN-YPT1
  718.           ZV2=ZPOIN-ZPT1
  719.           DENUM=XV2*XV1+YV2*YV1+ZV2*ZV1
  720.           DENOM=XVEC*XV1+YVEC*YV1+ZVEC*ZV1
  721.           IF (DENOM.EQ.0.D0) CALL ERREUR(21)
  722.           IF (IERR.NE.0) GOTO 10000
  723.           RAP=-DENUM/DENOM
  724.           XPOIN=XPOIN+RAP*XVEC
  725.           YPOIN=YPOIN+RAP*YVEC
  726.           ZPOIN=ZPOIN+RAP*ZVEC
  727.           GOTO 880
  728. C  Option PROJ SPHE
  729.  820      XV1=XPT1-XPOIN
  730.           YV1=YPT1-YPOIN
  731.           ZV1=ZPT1-ZPOIN
  732.           SCA=XVEC*XV1+YVEC*YV1+ZVEC*ZV1
  733.           XV1=XVEC*SCA
  734.           YV1=YVEC*SCA
  735.           ZV1=ZVEC*SCA
  736.           S2=(XPOIN+XV1-XPT1)**2+(YPOIN+YV1-YPT1)**2+(ZPOIN+ZV1-ZPT1)**2
  737.           IF (S2.GT.(ANGLE*ANGLE)) CALL ERREUR(21)
  738.           IF (IERR.NE.0) GOTO 10000
  739.           C=SQRT(ANGLE*ANGLE-S2)
  740.           IF (SCA.LT.0.D0) C=-C
  741.           XPOIN=XPOIN+XV1-C*XVEC
  742.           YPOIN=YPOIN+YV1-C*YVEC
  743.           ZPOIN=ZPOIN+ZV1-C*ZVEC
  744.           GOTO 880
  745. C  Option PROJ CYLI
  746.  830      XV2=YVEC*ZV1-ZVEC*YV1
  747.           YV2=ZVEC*XV1-XVEC*ZV1
  748.           ZV2=XVEC*YV1-YVEC*XV1
  749.           S2V2=XV2*XV2+YV2*YV2+ZV2*ZV2
  750.           IF (S2V2.EQ.0.D0) CALL ERREUR(21)
  751.           IF (IERR.NE.0) GOTO 10000
  752.           C2=(XVEC*XV1+YVEC*YV1+ZVEC*ZV1)**2
  753.           IF (C2.EQ.1.D0) CALL ERREUR(21)
  754.           IF (IERR.NE.0) GOTO 10000
  755.           XV3=XPT1-XPOIN
  756.           YV3=YPT1-YPOIN
  757.           ZV3=ZPT1-ZPOIN
  758.           XV4=YV3*ZV1-ZV3*YV1
  759.           YV4=ZV3*XV1-XV3*ZV1
  760.           ZV4=XV3*YV1-YV3*XV1
  761.           DNUM=(XV4*XVEC+YV4*YVEC+ZV4*ZVEC)**2
  762.           DIS2=DNUM/S2V2
  763.           IF (DIS2.GT.ANGLE**2) CALL ERREUR(21)
  764.           IF (IERR.NE.0) GOTO 10000
  765.           DMU=SQRT((ANGLE**2-DIS2)/(1.D0-C2))
  766.           DNUM=XV2*XV4+YV2*YV4+ZV2*ZV4
  767.           DLA=DNUM/S2V2
  768.           DMU=SIGN(DMU,DLA)
  769.           XPOIN=XPOIN+XVEC*(DLA-DMU)
  770.           YPOIN=YPOIN+YVEC*(DLA-DMU)
  771.           ZPOIN=ZPOIN+ZVEC*(DLA-DMU)
  772.           GOTO 880
  773. C  Option PROJ CONI
  774.  840      XV2=XPOIN-XPT1
  775.           YV2=YPOIN-YPT1
  776.           ZV2=ZPOIN-ZPT1
  777.           VECV1=XVEC*XV1+YVEC*YV1+ZVEC*ZV1
  778.           VEC2=XVEC*XVEC+YVEC*YVEC+ZVEC*ZVEC
  779.           V2V1=XV2*XV1+YV2*YV1+ZV2*ZV1
  780.           VECV2=XVEC*XV2+YVEC*YV2+ZVEC*ZV2
  781.           V22=XV2*XV2+YV2*YV2+ZV2*ZV2
  782.           A=VECV1*VECV1-ANGLE*VEC2
  783.           B=2*(V2V1*VECV1-ANGLE*VECV2)
  784.           C=V2V1*V2V1-ANGLE*V22
  785.           DELTA=B*B-4*A*C
  786.           IF (DELTA.LT.0.D0) CALL ERREUR(21)
  787.           IF (IERR.NE.0) GOTO 10000
  788.           DEL=SQRT(DELTA)
  789.           X1=(-B+DEL)/(2*A)
  790.           X2=(-B-DEL)/(2*A)
  791.           X=X2
  792.           IF (ABS(X1).LT.ABS(X2)) X=X1
  793.           XPOIN=XPOIN+X*XVEC
  794.           YPOIN=YPOIN+X*YVEC
  795.           ZPOIN=ZPOIN+X*ZVEC
  796.           GOTO 880
  797. C  Option PROJ TORI
  798.  850      PR2=XV2
  799.           GR2=YV2
  800.           XOP=XPOIN-XPT1
  801.           YOP=YPOIN-YPT1
  802.           ZOP=ZPOIN-ZPT1
  803.           OPV=XOP*XVEC+YOP*YVEC+ZOP*ZVEC
  804.           R2=XOP*XOP+YOP*YOP+ZOP*ZOP-GR2-PR2
  805.           VA=XVEC*XV1+YVEC*YV1+ZVEC*ZV1
  806.           QR2VA2=4*GR2*VA*VA
  807.           OPA=XOP*XV1+YOP*YV1+ZOP*ZV1
  808.           HR2PV=8*GR2*OPA*VA
  809.           R=4*GR2*OPA*OPA-4*PR2*GR2
  810.           RLMD=0.D0
  811. C  Resolution de l'equation du 4eme degre par iteration
  812.           DO iter=1,25
  813.             EXP1=RLMD*(RLMD+2*OPV)+R2
  814.             FDLM=EXP1*EXP1+QR2VA2*RLMD*RLMD+HR2PV*RLMD+R
  815.             FPDLM=4*EXP1*(RLMD+OPV)+QR2VA2*2*RLMD+HR2PV
  816.             IF (FPDLM.EQ.0.D0) CALL ERREUR(40)
  817.             IF (IERR.NE.0) GOTO 10000
  818.             CORR=FDLM/FPDLM
  819.             RLMD=RLMD-CORR
  820.             IF ((RLMD.NE.0.D0).AND.(ABS(CORR/RLMD).LT.1E-5)) GOTO 851
  821.           ENDDO
  822.           CALL ERREUR(40)
  823.           GOTO 10000
  824.  851      XPOIN=XPOIN+RLMD*XVEC
  825.           YPOIN=YPOIN+RLMD*YVEC
  826.           ZPOIN=ZPOIN+RLMD*ZVEC
  827.           GOTO 880
  828. C  Option PROJ DROI
  829.  860      XV2=XPOIN-XPT1
  830.           YV2=YPOIN-YPT1
  831.           ZV2=0.D0
  832.           DENUM=XV2*XV1+YV2*YV1+ZV2*ZV1
  833.           DENOM=XVEC*XV1+YVEC*YV1+ZVEC*ZV1
  834.           IF (DENOM.EQ.0.D0) CALL ERREUR(21)
  835.           IF (IERR.NE.0) GOTO 10000
  836.           RAP=-DENUM/DENOM
  837.           XPOIN=XPOIN+RAP*XVEC
  838.           YPOIN=YPOIN+RAP*YVEC
  839.           ZPOIN=ZPOIN+RAP*ZVEC
  840.           GOTO 880
  841. C  Option PROJ CERC
  842.  870      XV1=XPT1-XPOIN
  843.           YV1=YPT1-YPOIN
  844.           ZV1=0.D0
  845.           SCA=XVEC*XV1+YVEC*YV1+ZVEC*ZV1
  846.           XV1=XVEC*SCA
  847.           YV1=YVEC*SCA
  848.           ZV1=ZVEC*SCA
  849.           S2=(XPOIN+XV1-XPT1)**2+(YPOIN+YV1-YPT1)**2+(ZPOIN+ZV1-ZPT1)**2
  850.           IF (S2.GT.ANGLE**2) CALL ERREUR(21)
  851.           IF (IERR.NE.0) GOTO 10000
  852.           C=SQRT(ANGLE*ANGLE-S2)
  853.           IF (SCA.LT.0.D0) C=-C
  854.           XPOIN=XPOIN+XV1-C*XVEC
  855.           YPOIN=YPOIN+YV1-C*YVEC
  856.           ZPOIN=ZPOIN+ZV1-C*ZVEC
  857.           GOTO 880
  858. C  Traitement commun des options
  859.  880      XCOOR(IREF+1)=XPOIN
  860.           XCOOR(IREF+2)=YPOIN
  861.           IF (IDIM.EQ.3) XCOOR(IREF+3)=ZPOIN
  862.         ENDIF
  863.       ENDDO
  864.       GOTO 1000
  865.  
  866. C  Option COOR : Changement de systeme de coordonnees
  867. C ---------------
  868.  1800 CALL LIRMOT(MCL5,2,ICYLI,1)
  869.       IF (IERR.NE.0) GOTO 10000
  870.       IF (ICYLI.EQ.1) THEN
  871.         CALL LIROBJ('POINT   ',IP1,1,IRETOU)
  872.         IF (IERR.NE.0)  GO TO 10000
  873.         CALL LIROBJ('POINT   ',IP2,1,IRETOU)
  874.         IF (IERR.NE.0)  GO TO 10000
  875.         IF (IDIM.EQ.3) CALL  LIROBJ('POINT   ',IP3,1,IRETOU)
  876.         IF (IERR.NE.0)  GO TO 10000
  877.       ENDIF
  878. C  Transformation des coordonnees cartesiennes en cylindriques
  879.       IF (ICYLI.EQ.1) THEN
  880.         IREF=(IP1-1)*idimp1
  881.         XP11=XCOOR(IREF+1)
  882.         XP12=XCOOR(IREF+2)
  883.         XP13=0.D0
  884.         IF (IDIM.EQ.3) XP13=XCOOR(IREF+3)
  885.         IREF=(IP2-1)*idimp1
  886.         XP21=XCOOR(IREF+1)-XP11
  887.         XP22=XCOOR(IREF+2)-XP12
  888.         XP23=0.D0
  889.         IF(IDIM.EQ.3) XP23 = XCOOR(IREF+3)-XP13
  890.         XL=SQRT(XP21*XP21+XP22*XP22+XP23*XP23)
  891.         IF (XL.EQ.0.D0) THEN
  892.           CALL ERREUR(21)
  893.           GOTO 10000
  894.         ENDIF
  895.         XP21=XP21/XL
  896.         XP22=XP22/XL
  897.         YP21=-XP22
  898.         YP22=XP21
  899.         YP23=0.D0
  900.         IF (IDIM.EQ.3) THEN
  901.           XP23=XP23/XL
  902.           XP31=XP21
  903.           XP32=XP22
  904.           XP33=XP23
  905.           IREF=(IP3-1)*idimp1
  906.           XP21=XCOOR(IREF+1)-XP11
  907.           XP22=XCOOR(IREF+2)-XP12
  908.           XP23=XCOOR(IREF+3)-XP13
  909.           XXX=XP21*XP31+XP22*XP32+XP23*XP33
  910.           XP21=XP21-XXX*XP31
  911.           XP22=XP22-XXX*XP32
  912.           XP23=XP23-XXX*XP33
  913.           XL=SQRT(XP21*XP21+XP22*XP22+XP23*XP23)
  914.           IF (XL.EQ.0.D0) THEN
  915.             CALL ERREUR(21)
  916.             GOTO 10000
  917.           ENDIF
  918.           XP21=XP21/XL
  919.           XP22=XP22/XL
  920.           XP23=XP23/XL
  921.           YP21=XP32*XP23-XP33*XP22
  922.           YP22=XP33*XP21-XP31*XP23
  923.           YP23=XP31*XP22-XP32*XP21
  924.         ENDIF
  925.         SCA=180.D0/XPI
  926.         DO i=1,ICPR(/1)
  927.           IF (ICPR(i).NE.0) THEN
  928.             IREF=(i-1)*idimp1
  929.             XD=XCOOR(IREF+1)-XP11
  930.             YD=XCOOR(IREF+2)-XP12
  931.             ZD=0.D0
  932.             IF (IDIM.EQ.3) THEN
  933.               ZD=XCOOR(IREF+3)-XP13
  934.               ZF=XD*XP31+YD*XP32+ZD*XP33
  935.               XD=XD-ZF*XP31
  936.               YD=YD-ZF*XP32
  937.               ZD=ZD-ZF*XP33
  938.             ENDIF
  939.             XA=XD*XP21+YD*XP22+ZD*XP23
  940.             YA=XD*YP21+YD*YP22+ZD*YP23
  941.             RHO=SQRT(XA*XA+YA*YA)
  942.             TETA=ATAN2(YA,XA)
  943.             XCOOR(IREF+1)=RHO
  944.             XCOOR(IREF+2)=TETA*SCA
  945.             IF (IDIM.EQ.3) XCOOR(IREF+3)=ZF
  946.           ENDIF
  947.         ENDDO
  948. C  Transformation des coordonnees cylindriques en cartesiennes
  949. C  L'axe X est conserve, ainsi que l'axe Z.
  950.       ELSE IF (ICYLI.EQ.2) THEN
  951.         SCA=XPI/180.D0
  952.         DO i=1,ICPR(/1)
  953.           IF (ICPR(I).NE.0) THEN
  954.             IREF=(i-1)*idimp1
  955.             XD=XCOOR(IREF+1)
  956.             YD=XCOOR(IREF+2)*SCA
  957.             XCOOR(IREF+1)=XD*COS(YD)
  958.             XCOOR(IREF+2)=XD*SIN(YD)
  959.           ENDIF
  960.         ENDDO
  961.       ENDIF
  962.       GOTO 1000
  963.  
  964. C  Fin normale de l'operateur - Ecriture des objets (ICLE = 1 a 8)
  965. C -----------------------------------------------------------------
  966.  1000 SEGSUP ICPR
  967.       IF (ITYPLU.EQ.1) CALL ECROBJ('POINT',ITYPSA)
  968.       IF (ITYPLU.EQ.2) CALL ECROBJ('MAILLAGE',ITYPSA)
  969.       RETURN
  970.  
  971. C  Option MILI :
  972. C ---------------
  973.  1900 IF (IDIM.GE.2) THEN
  974.         CALL ECROBJ('MAILLAGE',MELEME)
  975.         CALL CHANLG
  976.         CALL LIROBJ('MAILLAGE',MELEME,1,IRETOU)
  977.         IF (IERR.NE.0) RETURN
  978.       ENDIF
  979.       SEGACT MELEME
  980.       IPT1=MELEME
  981.       INSO=LISOUS(/1)
  982.       DO j=1,MAX(1,INSO)
  983.         IF (INSO.NE.0) THEN
  984.           IPT1=LISOUS(j)
  985.           SEGACT IPT1
  986.         ENDIF
  987.         IF (IPT1.ITYPEL.EQ.3) THEN
  988.           DO i=1,IPT1.NUM(/2)
  989.             IRef1=idimp1*(IPT1.NUM(1,i)-1)
  990.             IRef2=idimp1*(IPT1.NUM(3,i)-1)
  991.             IRef3=idimp1*(IPT1.NUM(2,i)-1)
  992.             X1=XCOOR(IRef1+1)
  993.             X2=XCOOR(IRef2+1)
  994.             X3=XCOOR(IRef3+1)
  995.             X2X1=X2-X1
  996.             X3M=(X2+X1)*0.5-X3
  997.             XNUM=X2X1*X3M
  998.             XDEN=X2X1*X2X1
  999.             IF (IDIM.GE.2) THEN
  1000.               Y1=XCOOR(IRef1+2)
  1001.               Y2=XCOOR(IRef2+2)
  1002.               Y3=XCOOR(IRef3+2)
  1003.               Y2Y1=Y2-Y1
  1004.               Y3M=(Y2+Y1)*0.5-Y3
  1005.               XNUM=XNUM+Y2Y1*Y3M
  1006.               XDEN=XDEN+Y2Y1*Y2Y1
  1007.               IF (IDIM.GE.3) THEN
  1008.                 Z1=XCOOR(IRef1+3)
  1009.                 Z2=XCOOR(IRef2+3)
  1010.                 Z3=XCOOR(IRef3+3)
  1011.                 Z2Z1=Z2-Z1
  1012.                 Z3M=(Z2+Z1)*0.5-Z3
  1013.                 XNUM=XNUM+Z2Z1*Z3M
  1014.                 XDEN=XDEN+Z2Z1*Z2Z1
  1015.               ENDIF
  1016.             ENDIF
  1017.             IF (XDEN.LT.1.D-20)  THEN
  1018.                 wrITE (6,*) 'noeud double ' , IPT1.NUM(1,i)
  1019.                 xnum = 0.D0
  1020.                 xden = 1.D0
  1021. c              RETURN
  1022.             endif
  1023.             XLAMBD=XNUM/XDEN
  1024.             XCOOR(IRef3+1)=X3+XLAMBD*X2X1
  1025.             IF (IDIM.GE.2) THEN
  1026.               XCOOR(IRef3+2)=Y3+XLAMBD*Y2Y1
  1027.               IF (IDIM.GE.3) THEN
  1028.                 XCOOR(IRef3+3)=Z3+XLAMBD*Z2Z1
  1029.               ENDIF
  1030.             ENDIF
  1031.           ENDDO
  1032.         ENDIF
  1033.         SEGDES IPT1
  1034.       ENDDO
  1035.       SEGDES MELEME
  1036.       RETURN
  1037.  
  1038. C  Traitement des erreurs (options ICLE = 1 a 8)
  1039. C -----------------------------------------------
  1040. 10000 SEGSUP ICPR
  1041.       RETURN
  1042.  
  1043.       END
  1044.  
  1045.  
  1046.  
  1047.  
  1048.  
  1049.  
  1050.  
  1051.  
  1052.  
  1053.  
  1054.  
  1055.  
  1056.  

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