Télécharger la plaquette Cast3M


fr en

Télécharger ptori.eso

Retour à la liste

Numérotation des lignes : 

ptori

  1. C PTORI     SOURCE    PV        20/03/24    21:21:03     10554          
  2. C    CE SOUS-PROGRAMME RAMENE UN TORE SUR SES COORDONNEES PROPRES
  3. C    LES DONNEES DU TORE SONT : LE CENTRE UN POINT DE L'AXE ET UN CENTRE
  4. C    SECONDAIRE
  5. C
  6.       SUBROUTINE PTORI(IOP,FER,XPROJ,NDEB,NUMNP,IPC,IPA,IP1,tcval,isens)
  7.       IMPLICIT INTEGER(I-N)
  8.       IMPLICIT REAL*8 (A-H,O-Z)
  9. -INC CCREEL
  10. -INC SMCOORD
  11.  
  12. -INC PPARAM
  13. -INC CCOPTIO
  14.       real*8 tcval(*)
  15.       SEGMENT/FER/(NFI(ITT),MAI(IPP),ITOUR)
  16.       SEGMENT XPROJ(3,IMAX)
  17. *      tcval(1)     2    3     4     5    6
  18. *      SAVE XORIG,YORIG,ZORIG,XAXE,YAXE,ZAXE
  19. *            7   8   9
  20. *      SAVE XPC,YPC,ZPC
  21. *                  10   11
  22. *      SAVE ISENS,GRAY,PRAY
  23.       IF (IOP.EQ.2) GOTO 100
  24.       IMCT=MAI(ITOUR+1)
  25.       INCT=MAI(1)+1
  26.       IMAX=(IMCT**2)+10
  27.       CALL LIRENT(IMAX,0,IRETOU)
  28.       IF (IRETOU.NE.0) IMAX=MAX(1,IMAX)
  29.       NDEB=IMCT+1
  30.       SEGINI XPROJ
  31.       SEGACT MCOORD*mod
  32. C   CENTRE DU TORE
  33.       IREF=IPC*4-3
  34.       XPC=XCOOR(IREF)
  35.       YPC=XCOOR(IREF+1)
  36.       ZPC=XCOOR(IREF+2)
  37.       tcval(7)=xpc
  38.       tcval(8)=ypc
  39.       tcval(9)=zpc
  40. C   POINT DE L'AXE
  41.       IREF=IPA*4-3
  42.       XPA=XCOOR(IREF)
  43.       YPA=XCOOR(IREF+1)
  44.       ZPA=XCOOR(IREF+2)
  45. C   CENTRE SECONDAIRE
  46.       IREF=IP1*4-3
  47.       XP1=XCOOR(IREF)
  48.       YP1=XCOOR(IREF+1)
  49.       ZP1=XCOOR(IREF+2)
  50.       XAXE=XPC-XPA
  51.       YAXE=YPC-YPA
  52.       ZAXE=ZPC-ZPA
  53.       DAXE=SQRT(XAXE**2+YAXE**2+ZAXE**2)
  54.       IF (DAXE.EQ.0.D0) CALL ERREUR(21)
  55.       XAXE=XAXE/DAXE
  56.       YAXE=YAXE/DAXE
  57.       ZAXE=ZAXE/DAXE
  58.       tcval(4)=xaxe
  59.       tcval(5)=yaxe
  60.       tcval(6)=zaxe
  61. C  CALCUL DU GRAND RAYON
  62.       GRAY=SQRT((XPC-XP1)**2+(YPC-YP1)**2+(ZPC-ZP1)**2)
  63.       tcval(10)=gray
  64.       IF (IIMPI.NE.0) WRITE (IOIMP,9901) GRAY
  65.  9901 FORMAT (/,' GRAND RAYON ',G12.5)
  66.       IF (GRAY.EQ.0.D0) CALL ERREUR(21)
  67. C  DEROULONS LE TORE (DANS LES DEUX SENS)
  68.       IREF=4*NFI(IMCT)-3
  69.       XV1=XCOOR(IREF)-XPC
  70.       YV1=XCOOR(IREF+1)-YPC
  71.       ZV1=XCOOR(IREF+2)-ZPC
  72.       XI=XV1*XAXE+YV1*YAXE+ZV1*ZAXE
  73.       XC1=XV1-XI*XAXE
  74.       YC1=YV1-XI*YAXE
  75.       ZC1=ZV1-XI*ZAXE
  76. C  CENTRE SECONDAIRE
  77.       SC1=SQRT(XC1**2+YC1**2+ZC1**2)
  78.       U1=XC1/SC1
  79.       V1=YC1/SC1
  80.       W1=ZC1/SC1
  81.       XD1=XPC+U1*GRAY
  82.       YD1=YPC+V1*GRAY
  83.       ZD1=ZPC+W1*GRAY
  84.       XORIG=U1
  85.       YORIG=V1
  86.       ZORIG=W1
  87.       tcval(1)=xorig
  88.       tcval(2)=yorig
  89.       tcval(3)=zorig
  90.       U2=YAXE*W1-ZAXE*V1
  91.       V2=ZAXE*U1-XAXE*W1
  92.       W2=XAXE*V1-YAXE*U1
  93.       XVOR1=XV1+XPC-XD1
  94.       YVOR1=YV1+YPC-YD1
  95.       ZVOR1=ZV1+ZPC-ZD1
  96.       COT2=ATAN2(U2*(YAXE*ZVOR1-ZAXE*YVOR1)+V2*
  97.      #     (ZAXE*XVOR1-XAXE*ZVOR1)+W2*(XAXE*YVOR1-YAXE*XVOR1),
  98.      #     XAXE*XVOR1+YAXE*YVOR1+ZAXE*ZVOR1)
  99.       RAYON=0.D0
  100.       COT1=0.D0
  101.       DO 1 I=INCT,max(IMCT,mai(itour+2))
  102.          II=NFI(I)
  103.          IREF=4*II-3
  104.          XV2=XCOOR(IREF)-XPC
  105.          YV2=XCOOR(IREF+1)-YPC
  106.          ZV2=XCOOR(IREF+2)-ZPC
  107.          XI=XV2*XAXE+YV2*YAXE+ZV2*ZAXE
  108.          XC2=XV2-XI*XAXE
  109.          YC2=YV2-XI*YAXE
  110.          ZC2=ZV2-XI*ZAXE
  111.          ANG=ATAN2(XAXE*(YC1*ZC2-YC2*ZC1)+YAXE*(ZC1*XC2-ZC2*XC1)+
  112.      #        ZAXE*(XC1*YC2-XC2*YC1),XC1*XC2+YC1*YC2+ZC1*ZC2)
  113.          COT1=COT1+ANG
  114.          XPROJ(1,I)=COT1*GRAY
  115. C  CENTRE SECONDAIRE
  116.          SC2=SQRT(XC2**2+YC2**2+ZC2**2)
  117.          if (sc2.eq.0.D0) call erreur(21)
  118.          U1=XC2/SC2
  119.          V1=YC2/SC2
  120.          W1=ZC2/SC2
  121.          XD2=XPC+U1*GRAY
  122.          YD2=YPC+V1*GRAY
  123.          ZD2=ZPC+W1*GRAY
  124.          U2=YAXE*W1-ZAXE*V1
  125.          V2=ZAXE*U1-XAXE*W1
  126.          W2=XAXE*V1-YAXE*U1
  127.          XVOR2=XV2+XPC-XD2
  128.          YVOR2=YV2+YPC-YD2
  129.          ZVOR2=ZV2+ZPC-ZD2
  130.          ANG=ATAN2(U2*(YAXE*ZVOR2-ZAXE*YVOR2)+V2*
  131.      #        (ZAXE*XVOR2-XAXE*ZVOR2)+W2*(XAXE*YVOR2-YAXE*XVOR2),
  132.      #        XAXE*XVOR2+YAXE*YVOR2+ZAXE*ZVOR2)
  133.          ADIF=ANG-COT2
  134.          ADIF=ADIF-INT(ADIF/(2*XPI))*2*XPI
  135.          IF (ADIF.GE.XPI) ADIF=ADIF-2*XPI
  136.          IF (ADIF.LE.-XPI) ADIF=ADIF+2*XPI
  137.          COT2=COT2+ADIF
  138.          XPROJ(2,I)=COT2
  139.          XPROJ(3,I)=XCOOR(IREF+3)
  140.          if (i.le.imct) RAYON=RAYON+XVOR2**2+YVOR2**2+ZVOR2**2
  141.          XC1=XC2
  142.          YC1=YC2
  143.          ZC1=ZC2
  144.  1    CONTINUE
  145.       RAYON2=RAYON/(IMCT-INCT+1)
  146.       PRAY=SQRT(RAYON2)
  147.       tcval(11)=pray
  148.       IF (IIMPI.NE.0) WRITE (IOIMP,9902) PRAY
  149.  9902 FORMAT (' PETIT RAYON ',G12.5)
  150.       DO 2 I=INCT,max(IMCT,mai(itour+2))
  151.          XPROJ(2,I)=XPROJ(2,I)*PRAY
  152.          II=NFI(I)
  153.          NFI(I)=I
  154.          IREF=4*II-3
  155.          XV=XCOOR(IREF)-XPC
  156.          YV=XCOOR(IREF+1)-YPC
  157.          ZV=XCOOR(IREF+2)-ZPC
  158.          XI=XV*XAXE+YV*YAXE+ZV*ZAXE
  159.          XC=XV-XI*XAXE
  160.          YC=YV-XI*YAXE
  161.          ZC=ZV-XI*ZAXE
  162. C  CENTRE SECONDAIRE
  163.          SC=SQRT(XC**2+YC**2+ZC**2)
  164.          U1=XC/SC
  165.          V1=YC/SC
  166.          W1=ZC/SC
  167.          XC=XPC+U1*GRAY
  168.          YC=YPC+V1*GRAY
  169.          ZC=ZPC+W1*GRAY
  170.          XVOR=XV+XPC-XC
  171.          YVOR=YV+YPC-YC
  172.          ZVOR=ZV+ZPC-ZC
  173.          RAY2=XVOR**2+YVOR**2+ZVOR**2
  174.          RAP=RAY2/RAYON2
  175.          IF (i.gt.imct.or.(RAP.LE.1.02D0.AND.RAP.GE.0.98D0)) GOTO 2
  176.          IF (IIMPI.NE.0) WRITE (IOIMP,9903) I,RAP
  177.  9903    FORMAT(' POINT ',I6,' ERREUR RELATIVE SUR LE PETIT RAYON ',G12
  178.      $        .5)
  179.          CALL ERREUR(21)
  180.          RETURN
  181.  2    CONTINUE
  182. C  IL FAUT TOURNER DANS LE BON SENS
  183.       SURF=0.D0
  184.       DO 3 IT=1,ITOUR
  185.          II1=MAI(IT-1+1)+1
  186.          II2=MAI(IT+1)
  187.          XV1=XPROJ(1,II2)
  188.          YV1=XPROJ(2,II2)
  189.          DO 31 I=II1,II2
  190.             XV2=XPROJ(1,I)
  191.             YV2=XPROJ(2,I)
  192.             IF (XPROJ(3,I).EQ.0.D0) XPROJ(3,I)=SQRT((XV2-XV1)**2+(YV2
  193.      $           -YV1)**2)
  194.             SURF=SURF+XV1*YV2-XV2*YV1
  195.             XV1=XV2
  196.             YV1=YV2
  197.  31      CONTINUE
  198.  3    CONTINUE
  199.       ISENS=1
  200.       IF (SURF.GT.0.D0) GOTO 5
  201.       ISENS=-1
  202.       DO 4 I=INCT,max(IMCT,mai(itour+2))
  203.          XPROJ(1,I)=-XPROJ(1,I)
  204.  4    CONTINUE
  205.  5    CONTINUE
  206.       RETURN
  207. C     TRANSFORMATION INVERSE
  208.  100  CONTINUE
  209.       xorig=tcval(1)
  210.       yorig=tcval(2)
  211.       zorig=tcval(3)
  212.       xaxe=tcval(4)
  213.       yaxe=tcval(5)
  214.       zaxe=tcval(6)
  215.       xpc=tcval(7)
  216.       ypc=tcval(8)
  217.       zpc=tcval(9)
  218.       gray= tcval(10)
  219.       pray=tcval(11)
  220.       SEGACT MCOORD*mod
  221.       XREP1=XORIG
  222.       YREP1=YORIG
  223.       ZREP1=ZORIG
  224.       XREP2=YAXE*ZREP1-ZAXE*YREP1
  225.       YREP2=ZAXE*XREP1-XAXE*ZREP1
  226.       ZREP2=XAXE*YREP1-YAXE*XREP1
  227.       IF (NDEB.GT.NUMNP) GOTO 102
  228.       NBPTA=nbpts
  229.       NBPTS=NBPTA+NUMNP-NDEB+1
  230.       SEGADJ MCOORD
  231.       DO 101 I=NDEB,NUMNP
  232.          ANG=XPROJ(1,I)*ISENS/GRAY
  233.          SI=SIN(ANG)
  234.          CO=COS(ANG)
  235.          XC1=(XREP1*CO+XREP2*SI)*GRAY
  236.          YC1=(YREP1*CO+YREP2*SI)*GRAY
  237.          ZC1=(ZREP1*CO+ZREP2*SI)*GRAY
  238.          SC1=SQRT(XC1**2+YC1**2+ZC1**2)
  239.          ANG=XPROJ(2,I)/PRAY
  240.          SI=SIN(ANG)
  241.          CO=COS(ANG)
  242.          XV1=(XAXE*CO+XC1*SI/SC1)*PRAY
  243.          YV1=(YAXE*CO+YC1*SI/SC1)*PRAY
  244.          ZV1=(ZAXE*CO+ZC1*SI/SC1)*PRAY
  245.          XCOOR(NBPTA*(IDIM+1)+1)=XPC+XC1+XV1
  246.          XCOOR(NBPTA*(IDIM+1)+2)=YPC+YC1+YV1
  247.          XCOOR(NBPTA*(IDIM+1)+3)=ZPC+ZC1+ZV1
  248.          XCOOR((NBPTA+1)*(IDIM+1))=XPROJ(3,I)
  249.          NBPTA=NBPTA+1
  250.  101  CONTINUE
  251.  102  CONTINUE
  252.       SEGSUP XPROJ
  253.       RETURN
  254.       END
  255.  
  256.  
  257.  
  258.  
  259.  
  260.  
  261.  
  262.  

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