C PTORI     SOURCE    PV        20/03/24    21:21:03     10554          
C    CE SOUS-PROGRAMME RAMENE UN TORE SUR SES COORDONNEES PROPRES
C    LES DONNEES DU TORE SONT : LE CENTRE UN POINT DE L'AXE ET UN CENTRE
C    SECONDAIRE
C
      SUBROUTINE PTORI(IOP,FER,XPROJ,NDEB,NUMNP,IPC,IPA,IP1,tcval,isens)
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)
-INC CCREEL
-INC SMCOORD

-INC PPARAM
-INC CCOPTIO
      real*8 tcval(*)
      SEGMENT/FER/(NFI(ITT),MAI(IPP),ITOUR)
      SEGMENT XPROJ(3,IMAX)
*      tcval(1)     2    3     4     5    6
*      SAVE XORIG,YORIG,ZORIG,XAXE,YAXE,ZAXE
*            7   8   9
*      SAVE XPC,YPC,ZPC
*                  10   11
*      SAVE ISENS,GRAY,PRAY
      IF (IOP.EQ.2) GOTO 100
      IMCT=MAI(ITOUR+1)
      INCT=MAI(1)+1
      IMAX=(IMCT**2)+10
      CALL LIRENT(IMAX,0,IRETOU)
      IF (IRETOU.NE.0) IMAX=MAX(1,IMAX)
      NDEB=IMCT+1
      SEGINI XPROJ
      SEGACT MCOORD*mod
C   CENTRE DU TORE
      IREF=IPC*4-3
      XPC=XCOOR(IREF)
      YPC=XCOOR(IREF+1)
      ZPC=XCOOR(IREF+2)
      tcval(7)=xpc
      tcval(8)=ypc
      tcval(9)=zpc
C   POINT DE L'AXE
      IREF=IPA*4-3
      XPA=XCOOR(IREF)
      YPA=XCOOR(IREF+1)
      ZPA=XCOOR(IREF+2)
C   CENTRE SECONDAIRE
      IREF=IP1*4-3
      XP1=XCOOR(IREF)
      YP1=XCOOR(IREF+1)
      ZP1=XCOOR(IREF+2)
      XAXE=XPC-XPA
      YAXE=YPC-YPA
      ZAXE=ZPC-ZPA
      DAXE=SQRT(XAXE**2+YAXE**2+ZAXE**2)
      IF (DAXE.EQ.0.D0) CALL ERREUR(21)
      XAXE=XAXE/DAXE
      YAXE=YAXE/DAXE
      ZAXE=ZAXE/DAXE
      tcval(4)=xaxe
      tcval(5)=yaxe
      tcval(6)=zaxe
C  CALCUL DU GRAND RAYON
      GRAY=SQRT((XPC-XP1)**2+(YPC-YP1)**2+(ZPC-ZP1)**2)
      tcval(10)=gray
      IF (IIMPI.NE.0) WRITE (IOIMP,9901) GRAY
 9901 FORMAT (/,' GRAND RAYON ',G12.5)
      IF (GRAY.EQ.0.D0) CALL ERREUR(21)
C  DEROULONS LE TORE (DANS LES DEUX SENS)
      IREF=4*NFI(IMCT)-3
      XV1=XCOOR(IREF)-XPC
      YV1=XCOOR(IREF+1)-YPC
      ZV1=XCOOR(IREF+2)-ZPC
      XI=XV1*XAXE+YV1*YAXE+ZV1*ZAXE
      XC1=XV1-XI*XAXE
      YC1=YV1-XI*YAXE
      ZC1=ZV1-XI*ZAXE
C  CENTRE SECONDAIRE
      SC1=SQRT(XC1**2+YC1**2+ZC1**2)
      U1=XC1/SC1
      V1=YC1/SC1
      W1=ZC1/SC1
      XD1=XPC+U1*GRAY
      YD1=YPC+V1*GRAY
      ZD1=ZPC+W1*GRAY
      XORIG=U1
      YORIG=V1
      ZORIG=W1
      tcval(1)=xorig
      tcval(2)=yorig
      tcval(3)=zorig
      U2=YAXE*W1-ZAXE*V1
      V2=ZAXE*U1-XAXE*W1
      W2=XAXE*V1-YAXE*U1
      XVOR1=XV1+XPC-XD1
      YVOR1=YV1+YPC-YD1
      ZVOR1=ZV1+ZPC-ZD1
      COT2=ATAN2(U2*(YAXE*ZVOR1-ZAXE*YVOR1)+V2*
     #     (ZAXE*XVOR1-XAXE*ZVOR1)+W2*(XAXE*YVOR1-YAXE*XVOR1),
     #     XAXE*XVOR1+YAXE*YVOR1+ZAXE*ZVOR1)
      RAYON=0.D0
      COT1=0.D0
      DO 1 I=INCT,max(IMCT,mai(itour+2))
         II=NFI(I)
         IREF=4*II-3
         XV2=XCOOR(IREF)-XPC
         YV2=XCOOR(IREF+1)-YPC
         ZV2=XCOOR(IREF+2)-ZPC
         XI=XV2*XAXE+YV2*YAXE+ZV2*ZAXE
         XC2=XV2-XI*XAXE
         YC2=YV2-XI*YAXE
         ZC2=ZV2-XI*ZAXE
         ANG=ATAN2(XAXE*(YC1*ZC2-YC2*ZC1)+YAXE*(ZC1*XC2-ZC2*XC1)+
     #        ZAXE*(XC1*YC2-XC2*YC1),XC1*XC2+YC1*YC2+ZC1*ZC2)
         COT1=COT1+ANG
         XPROJ(1,I)=COT1*GRAY
C  CENTRE SECONDAIRE
         SC2=SQRT(XC2**2+YC2**2+ZC2**2)
         if (sc2.eq.0.D0) call erreur(21)
         U1=XC2/SC2
         V1=YC2/SC2
         W1=ZC2/SC2
         XD2=XPC+U1*GRAY
         YD2=YPC+V1*GRAY
         ZD2=ZPC+W1*GRAY
         U2=YAXE*W1-ZAXE*V1
         V2=ZAXE*U1-XAXE*W1
         W2=XAXE*V1-YAXE*U1
         XVOR2=XV2+XPC-XD2
         YVOR2=YV2+YPC-YD2
         ZVOR2=ZV2+ZPC-ZD2
         ANG=ATAN2(U2*(YAXE*ZVOR2-ZAXE*YVOR2)+V2*
     #        (ZAXE*XVOR2-XAXE*ZVOR2)+W2*(XAXE*YVOR2-YAXE*XVOR2),
     #        XAXE*XVOR2+YAXE*YVOR2+ZAXE*ZVOR2)
         ADIF=ANG-COT2
         ADIF=ADIF-INT(ADIF/(2*XPI))*2*XPI
         IF (ADIF.GE.XPI) ADIF=ADIF-2*XPI
         IF (ADIF.LE.-XPI) ADIF=ADIF+2*XPI
         COT2=COT2+ADIF
         XPROJ(2,I)=COT2
         XPROJ(3,I)=XCOOR(IREF+3)
         if (i.le.imct) RAYON=RAYON+XVOR2**2+YVOR2**2+ZVOR2**2
         XC1=XC2
         YC1=YC2
         ZC1=ZC2
 1    CONTINUE
      RAYON2=RAYON/(IMCT-INCT+1)
      PRAY=SQRT(RAYON2)
      tcval(11)=pray
      IF (IIMPI.NE.0) WRITE (IOIMP,9902) PRAY
 9902 FORMAT (' PETIT RAYON ',G12.5)
      DO 2 I=INCT,max(IMCT,mai(itour+2))
         XPROJ(2,I)=XPROJ(2,I)*PRAY
         II=NFI(I)
         NFI(I)=I
         IREF=4*II-3
         XV=XCOOR(IREF)-XPC
         YV=XCOOR(IREF+1)-YPC
         ZV=XCOOR(IREF+2)-ZPC
         XI=XV*XAXE+YV*YAXE+ZV*ZAXE
         XC=XV-XI*XAXE
         YC=YV-XI*YAXE
         ZC=ZV-XI*ZAXE
C  CENTRE SECONDAIRE
         SC=SQRT(XC**2+YC**2+ZC**2)
         U1=XC/SC
         V1=YC/SC
         W1=ZC/SC
         XC=XPC+U1*GRAY
         YC=YPC+V1*GRAY
         ZC=ZPC+W1*GRAY
         XVOR=XV+XPC-XC
         YVOR=YV+YPC-YC
         ZVOR=ZV+ZPC-ZC
         RAY2=XVOR**2+YVOR**2+ZVOR**2
         RAP=RAY2/RAYON2
         IF (i.gt.imct.or.(RAP.LE.1.02D0.AND.RAP.GE.0.98D0)) GOTO 2
         IF (IIMPI.NE.0) WRITE (IOIMP,9903) I,RAP
 9903    FORMAT(' POINT ',I6,' ERREUR RELATIVE SUR LE PETIT RAYON ',G12
     $        .5)
         CALL ERREUR(21)
         RETURN
 2    CONTINUE
C  IL FAUT TOURNER DANS LE BON SENS
      SURF=0.D0
      DO 3 IT=1,ITOUR
         II1=MAI(IT-1+1)+1
         II2=MAI(IT+1)
         XV1=XPROJ(1,II2)
         YV1=XPROJ(2,II2)
         DO 31 I=II1,II2
            XV2=XPROJ(1,I)
            YV2=XPROJ(2,I)
            IF (XPROJ(3,I).EQ.0.D0) XPROJ(3,I)=SQRT((XV2-XV1)**2+(YV2
     $           -YV1)**2)
            SURF=SURF+XV1*YV2-XV2*YV1
            XV1=XV2
            YV1=YV2
 31      CONTINUE
 3    CONTINUE
      ISENS=1
      IF (SURF.GT.0.D0) GOTO 5
      ISENS=-1
      DO 4 I=INCT,max(IMCT,mai(itour+2))
         XPROJ(1,I)=-XPROJ(1,I)
 4    CONTINUE
 5    CONTINUE
      RETURN
C     TRANSFORMATION INVERSE
 100  CONTINUE
      xorig=tcval(1)
      yorig=tcval(2)
      zorig=tcval(3)
      xaxe=tcval(4)
      yaxe=tcval(5)
      zaxe=tcval(6)
      xpc=tcval(7)
      ypc=tcval(8)
      zpc=tcval(9)
      gray= tcval(10)
      pray=tcval(11)
      SEGACT MCOORD*mod
      XREP1=XORIG
      YREP1=YORIG
      ZREP1=ZORIG
      XREP2=YAXE*ZREP1-ZAXE*YREP1
      YREP2=ZAXE*XREP1-XAXE*ZREP1
      ZREP2=XAXE*YREP1-YAXE*XREP1
      IF (NDEB.GT.NUMNP) GOTO 102
      NBPTA=nbpts
      NBPTS=NBPTA+NUMNP-NDEB+1
      SEGADJ MCOORD
      DO 101 I=NDEB,NUMNP
         ANG=XPROJ(1,I)*ISENS/GRAY
         SI=SIN(ANG)
         CO=COS(ANG)
         XC1=(XREP1*CO+XREP2*SI)*GRAY
         YC1=(YREP1*CO+YREP2*SI)*GRAY
         ZC1=(ZREP1*CO+ZREP2*SI)*GRAY
         SC1=SQRT(XC1**2+YC1**2+ZC1**2)
         ANG=XPROJ(2,I)/PRAY
         SI=SIN(ANG)
         CO=COS(ANG)
         XV1=(XAXE*CO+XC1*SI/SC1)*PRAY
         YV1=(YAXE*CO+YC1*SI/SC1)*PRAY
         ZV1=(ZAXE*CO+ZC1*SI/SC1)*PRAY
         XCOOR(NBPTA*(IDIM+1)+1)=XPC+XC1+XV1
         XCOOR(NBPTA*(IDIM+1)+2)=YPC+YC1+YV1
         XCOOR(NBPTA*(IDIM+1)+3)=ZPC+ZC1+ZV1
         XCOOR((NBPTA+1)*(IDIM+1))=XPROJ(3,I)
         NBPTA=NBPTA+1
 101  CONTINUE
 102  CONTINUE
      SEGSUP XPROJ
      RETURN
      END