1 : C CADRCL SOURCE PV 06/04/16 21:15:12 5405 2 : C TRACE DE DEFORMES 3 : C CALCUL DU CADRE GLOBAL 4 : C PROJECTION SUCCESSIVE DE CHAQUE DEFORME (SELON ICLE) 5 : C 1995 option DIRE P.PEGON JRC-ISPRA 6 : C 7 : C PP option DIRE 8 : SUBROUTINE CADRCL(KABCOR,LABCO2,IOEIL,XPROJ, 9 : # ICLE,XMINT,YMINT,XMAXT,YMAXT,zmint,zmaxt,cgrav,diloc,LDIRE, 10 : > axez) 11 : IMPLICIT INTEGER(I-N) 12 : REAL*8 XO,XG,XP,XN,SN,XV,SV,UI,UJ 13 : DIMENSION XO(3),XP(3),XN(3),XG(3),XV(3),UI(3),UJ(3),axez(3) 14 : dimension cgrav(*) 15 : C+PP option DIRE 16 : dimension diloc(*) 17 : C+PP 18 : COMMON /CCADRC/XN,XG,XO,UI,UJ 19 : C ATTENTION 20 : SEGMENT KABCOR(0) 21 : SEGMENT KABCO2(2,0) 22 : SEGMENT LABCO2(3,0) 23 : SEGMENT ICOR2(0) 24 : SEGMENT KXPRO2(NVEC) 25 : SEGMENT XCORD(IDIM,ITE) 26 : SEGMENT XCOR2(IDIM,ITE) 27 : -INC CCOPTIO 28 : -INC SMCOORD 29 : SEGMENT XPROJ(3,ITE) 30 : SEGMENT XPRO2(3,ITE) 31 : C+PP option DIRE 32 : REAL*8 SSN 33 : LOGICAL LDIRE 34 : C+PP 35 : C ICLE = 0 RECHERCHE DES MIN MAX SUR TOUTES LES DEFORMES 36 : C CALCUL DE LA PROJECTION 37 : IF (ICLE.NE.0) GOTO 1000 38 : XMINT=1E30 39 : YMINT=XMINT 40 : ZMINT=XMINT 41 : XMAXT=-1E30 42 : YMAXT=XMAXT 43 : ZMAXT=XMAXT 44 : NCORD=KABCOR(/1) 45 : IF (IDIM.EQ.3) GOTO 50 46 : DO 10 ICORD=1,NCORD 47 : XCORD=KABCOR(ICORD) 48 : ITE=XCORD(/2) 49 : DO 19 J=1,ITE 50 : XMINT=MIN(XMINT,XCORD(1,J)) 51 : XMAXT=MAX(XMAXT,XCORD(1,J)) 52 : YMINT=MIN(YMINT,XCORD(2,J)) 53 : YMAXT=MAX(YMAXT,XCORD(2,J)) 54 : ZMINT=0 55 : ZMAXT=0 56 : 19 CONTINUE 57 : 10 CONTINUE 58 : IF (LABCO2.EQ.0) GOTO 60 59 : DO 61 ICORD=1,NCORD 60 : IF (LABCO2(3,ICORD).EQ.0) GOTO 61 61 : KABCO2=LABCO2(1,ICORD) 62 : XCORD=KABCOR(ICORD) 63 : ITE=XCORD(/2) 64 : NVEC=KABCO2(/2) 65 : IF (NVEC.EQ.0) GOTO 61 66 : DO 7 IVEC=1,NVEC 67 : XCOR2=KABCO2(1,IVEC) 68 : ICOR2=KABCO2(2,IVEC) 69 : ITE=XCOR2(/2) 70 : DO 18 J=1,ITE 71 : IF (ICOR2(J).EQ.0) GOTO 18 72 : UX=XCOR2(1,J)-XCORD(1,J) 73 : UY=XCOR2(2,J)-XCORD(2,J) 74 : U1=XCOR2(1,J)-UX/3-UY/5 75 : V1=XCOR2(2,J)-UY/3+UX/5 76 : XMINT=MIN(XMINT,U1) 77 : XMAXT=MAX(XMAXT,U1) 78 : YMINT=MIN(YMINT,V1) 79 : YMAXT=MAX(YMAXT,V1) 80 : U1=XCOR2(1,J)-UX/3+UY/5 81 : V1=XCOR2(2,J)-UY/3-UX/5 82 : XMINT=MIN(XMINT,U1) 83 : XMAXT=MAX(XMAXT,U1) 84 : YMINT=MIN(YMINT,V1) 85 : YMAXT=MAX(YMAXT,V1) 86 : XMINT=MIN(XMINT,XCOR2(1,J)) 87 : XMAXT=MAX(XMAXT,XCOR2(1,J)) 88 : YMINT=MIN(YMINT,XCOR2(2,J)) 89 : YMAXT=MAX(YMAXT,XCOR2(2,J)) 90 : 18 CONTINUE 91 : 7 CONTINUE 92 : 61 CONTINUE 93 : 60 CONTINUE 94 : RETURN 95 : 50 CONTINUE 96 : SEGACT MCOORD 97 : XMINT=1E30 98 : XMAXT=-XMINT 99 : YMINT=XMINT 100 : YMAXT=-YMINT 101 : ZMINT=XMINT 102 : ZMAXT=-YMINT 103 : IREF=(IOEIL-1)*4 104 : XO(1)=XCOOR(IREF+1) 105 : XO(2)=XCOOR(IREF+2) 106 : XO(3)=XCOOR(IREF+3) 107 : C POINT MOYEN 108 : DO 1 I=1,3 109 : XG(I)=0.D0 110 : 1 CONTINUE 111 : NPTOT=0 112 : DO 100 ICORD=1,NCORD 113 : XCORD=KABCOR(ICORD) 114 : ITE=XCORD(/2) 115 : NPTOT=NPTOT+ITE 116 : DO 2 I=1,ITE 117 : DO 3 J=1,3 118 : XG(J)=XG(J)+XCORD(J,I) 119 : 3 CONTINUE 120 : 2 CONTINUE 121 : 100 CONTINUE 122 : C+PP option DIRE 123 : IF (LDIRE)THEN 124 : DO J=1,3 125 : XG(J)=cgrav(J) 126 : XN(J)=XO(J)-XG(J) 127 : ENDDO 128 : ELSE 129 : C+PP 130 : DO 4 J=1,3 131 : XG(J)=XG(J)/NPTOT 132 : XN(J)=XO(J)-XG(J) 133 : cgrav(j)=xg(j) 134 : 4 CONTINUE 135 : C+PP option DIRE 136 : ENDIF 137 : C+PP 138 : C NORMALE AU PLAN 139 : SN=SQRT(XN(1)**2+XN(2)**2+XN(3)**2) 139 : IF (SN.EQ.0.) CALL ERREUR(21) 140 : IF (IERR.NE.0) RETURN 141 : C+PP option DIRE 142 : SSN=SN 143 : C+PP 144 : DO 5 J=1,3 145 : XN(J)=XN(J)/SN 146 : 5 CONTINUE 147 : C REPERE LOCAL SUR LE PLAN 148 : C+PP option DIRE 149 : IF (LDIRE)THEN 150 : DO J=1,3 151 : UJ(J)=diloc(J) 152 : ENDDO 153 : ELSE 154 : C+PP 155 : UJ(1)=0.D0 156 : UJ(2)=0.D0 157 : UJ(3)=1.D0 158 : C+PP option DIRE 159 : ENDIF 160 : C+PP 161 : 21 CONTINUE 162 : SV=UJ(1)*XN(1)+UJ(2)*XN(2)+UJ(3)*XN(3) 163 : DO 20 J=1,3 164 : UJ(J)=UJ(J)-SV*XN(J) 165 : 20 CONTINUE 166 : SV=UJ(1)**2+UJ(2)**2+UJ(3)**2 167 : IF (ABS(SV).LT.0.1) THEN 168 : UJ(1)=0.D0 169 : UJ(2)=1.D0 170 : UJ(3)=1.D0 171 : GOTO 21 172 : ENDIF 173 : SV=SQRT(SV) 174 : UJ(1)=UJ(1)/SV 175 : UJ(2)=UJ(2)/SV 176 : UJ(3)=UJ(3)/SV 177 : UI(1)=UJ(2)*XN(3)-UJ(3)*XN(2) 178 : UI(2)=UJ(3)*XN(1)-UJ(1)*XN(3) 179 : UI(3)=UJ(1)*XN(2)-UJ(2)*XN(1) 180 : C PROJECTION CONIQUE SUR LE PLAN 181 : DO 170 ICORD=1,NCORD 182 : XCORD=KABCOR(ICORD) 183 : ITE=XCORD(/2) 184 : DO 12 I=1,ITE 185 : DO 13 J=1,3 186 : XP(J)=XCORD(J,I) 187 : XV(J)=XP(J)-XO(J) 188 : 13 CONTINUE 189 : SV=XV(1)*XN(1)+XV(2)*XN(2)+XV(3)*XN(3) 190 : C+PP option DIRE 191 : IF ((LDIRE.AND.-SV.GE.SSN).OR.(.NOT.LDIRE))THEN 192 : C+PP 193 : SN=(XP(1)-XG(1))*XN(1)+(XP(2)-XG(2))*XN(2)+(XP(3)-XG(3))*XN(3) 194 : DO 14 J=1,3 195 : XP(J)=XP(J)-(SN/SV)*XV(J)-XG(J) 196 : 14 CONTINUE 197 : XPRO=XP(1)*UI(1)+XP(2)*UI(2)+XP(3)*UI(3) 198 : YPRO=XP(1)*UJ(1)+XP(2)*UJ(2)+XP(3)*UJ(3) 199 : ZPRO=-SN 200 : XMINT=MIN(XMINT,XPRO) 201 : XMAXT=MAX(XMAXT,XPRO) 202 : YMINT=MIN(YMINT,YPRO) 203 : YMAXT=MAX(YMAXT,YPRO) 204 : ZMINT=MIN(ZMINT,ZPRO) 205 : ZMAXT=MAX(ZMAXT,ZPRO) 206 : C+PP option DIRE 207 : ENDIF 208 : C+PP 209 : 12 CONTINUE 210 : 170 CONTINUE 211 : IF (LABCO2.EQ.0) GOTO 11 212 : DO 171 ICORD=1,NCORD 213 : IF (LABCO2(3,ICORD).EQ.0) GOTO 171 214 : KABCO2=LABCO2(1,ICORD) 215 : XCORD=KABCOR(ICORD) 216 : ITE=XCORD(/2) 217 : NVEC=KABCO2(/2) 218 : IF (NVEC.EQ.0) GOTO 171 219 : DO 6 IVEC=1,NVEC 220 : XCOR2=KABCO2(1,IVEC) 221 : ICOR2=KABCO2(2,IVEC) 222 : DO 17 I=1,ITE 223 : IF (ICOR2(I).EQ.0) GOTO 17 224 : C+PP a faire meme ss DIRE car sinon que signifie XPRO,YPRO? 225 : DO J=1,3 226 : XP(J)=XCORD(J,I) 227 : XV(J)=XP(J)-XO(J) 228 : ENDDO 229 : SV=XV(1)*XN(1)+XV(2)*XN(2)+XV(3)*XN(3) 230 : C+PP option DIRE 231 : IF (LDIRE.AND.-SV.LT.SSN) THEN 232 : ICOR2(I)=0 233 : GOTO 17 234 : ENDIF 235 : C+PP a faire meme ss DIRE car sinon que signifie XPRO,YPRO? 236 : SN=(XP(1)-XG(1))*XN(1)+(XP(2)-XG(2))*XN(2)+(XP(3)-XG(3))*XN(3) 237 : DO J=1,3 238 : XP(J)=XP(J)-(SN/SV)*XV(J)-XG(J) 239 : ENDDO 240 : XPRO=XP(1)*UI(1)+XP(2)*UI(2)+XP(3)*UI(3) 241 : YPRO=XP(1)*UJ(1)+XP(2)*UJ(2)+XP(3)*UJ(3) 242 : C+PP 243 : DO 15 J=1,3 244 : XP(J)=XCOR2(J,I) 245 : XV(J)=XP(J)-XO(J) 246 : 15 CONTINUE 247 : SV=XV(1)*XN(1)+XV(2)*XN(2)+XV(3)*XN(3) 248 : C+PP option DIRE 249 : IF (LDIRE.AND.-SV.LT.SSN) THEN 250 : ICOR2(I)=0 251 : GOTO 17 252 : ENDIF 253 : C+PP 254 : SN=(XP(1)-XG(1))*XN(1)+(XP(2)-XG(2))*XN(2)+(XP(3)-XG(3))*XN(3) 255 : DO 16 J=1,3 256 : XP(J)=XP(J)-(SN/SV)*XV(J)-XG(J) 257 : 16 CONTINUE 258 : XPROO=XP(1)*UI(1)+XP(2)*UI(2)+XP(3)*UI(3) 259 : YPROO=XP(1)*UJ(1)+XP(2)*UJ(2)+XP(3)*UJ(3) 260 : XMINT=MIN(XMINT,XPROO) 261 : XMAXT=MAX(XMAXT,XPROO) 262 : YMINT=MIN(YMINT,YPROO) 263 : YMAXT=MAX(YMAXT,YPROO) 264 : UX=XPROO-XPRO 265 : UY=YPROO-YPRO 266 : U1=XPROO-UX/3-UY/5 267 : V1=YPROO-UY/3+UX/5 268 : XMINT=MIN(XMINT,U1) 269 : XMAXT=MAX(XMAXT,V1) 270 : U1=XPROO-UX/3+UY/5 271 : V1=YPROO-UY/3-UX/5 272 : XMINT=MIN(XMINT,U1) 273 : XMAXT=MAX(XMAXT,V1) 274 : 17 CONTINUE 275 : 6 CONTINUE 276 : 171 CONTINUE 277 : 11 CONTINUE 278 : * axez pour tourner correctement avec opengl 279 : axez(1)=0 280 : axez(2)=uj(3) 281 : axez(3)=sqrt(1-uj(3)**2) 282 : if (xn(3).lt.0) axez(3)=-axez(3) 283 : * write (6,*) ' axez ',axez(1),axez(2),axez(3) 284 : RETURN 285 : 1000 CONTINUE 286 : C ON REMPLIT XPROJ POUR LA DEFORME CONCERNEE 287 : XCORD=KABCOR(ICLE) 288 : ITE=XCORD(/2) 289 : IF (IDIM.EQ.2) GOTO 1100 290 : C+PP option DIRE 291 : SSN=SQRT(XN(1)**2+XN(2)**2+XN(3)**2) 292 : SSN=0.95*SSN 293 : C+PP 294 : C PROJECTION CONIQUE SUR LE PLAN 295 : DO 612 I=1,ITE 296 : DO 613 J=1,3 297 : XP(J)=XCORD(J,I) 298 : XV(J)=XP(J)-XO(J) 299 : 613 CONTINUE 300 : SV=XV(1)*XN(1)+XV(2)*XN(2)+XV(3)*XN(3) 301 : SN=(XP(1)-XG(1))*XN(1)+(XP(2)-XG(2))*XN(2)+(XP(3)-XG(3))*XN(3) 302 : XPROJ(3,I)=-SN 303 : DO 614 J=1,3 304 : XP(J)=XP(J)-(SN/SV)*XV(J)-XG(J) 305 : 614 CONTINUE 306 : C+PP option DIRE 307 : IF ((LDIRE.AND.-SV.GE.SSN).OR.(.NOT.LDIRE))THEN 308 : C+PP 309 : XPROJ(1,I)=XP(1)*UI(1)+XP(2)*UI(2)+XP(3)*UI(3) 310 : XPROJ(2,I)=XP(1)*UJ(1)+XP(2)*UJ(2)+XP(3)*UJ(3) 311 : C+PP option DIRE 312 : ELSE 313 : XPROJ(1,I)=(XMINT+XMAXT)/2 314 : XPROJ(2,I)=(YMINT+YMAXT)/2 315 : ENDIF 316 : C+PP 317 : 612 CONTINUE 318 : IF (LABCO2.EQ.0) GOTO 618 319 : IF (LABCO2(3,ICLE).EQ.0) GOTO 618 320 : KABCO2=LABCO2(1,ICLE) 321 : NVEC=KABCO2(/2) 322 : SEGINI KXPRO2 323 : LABCO2(2,ICLE)=KXPRO2 324 : IF (NVEC.EQ.0) GOTO 618 325 : DO 8 IVEC=1,NVEC 326 : XCOR2=KABCO2(1,IVEC) 327 : ICOR2=KABCO2(2,IVEC) 328 : SEGINI XPRO2 329 : KXPRO2(IVEC)=XPRO2 330 : DO 617 I=1,ITE 331 : IF (ICOR2(I).EQ.0) GOTO 617 332 : DO 615 J=1,3 333 : XP(J)=XCOR2(J,I) 334 : XV(J)=XP(J)-XO(J) 335 : 615 CONTINUE 336 : SV=XV(1)*XN(1)+XV(2)*XN(2)+XV(3)*XN(3) 337 : SN=(XP(1)-XG(1))*XN(1)+(XP(2)-XG(2))*XN(2)+(XP(3)-XG(3))*XN(3) 338 : XPRO2(3,I)=-SN 339 : DO 616 J=1,3 340 : XP(J)=XP(J)-(SN/SV)*XV(J)-XG(J) 341 : 616 CONTINUE 342 : XPRO2(1,I)=XP(1)*UI(1)+XP(2)*UI(2)+XP(3)*UI(3) 343 : XPRO2(2,I)=XP(1)*UJ(1)+XP(2)*UJ(2)+XP(3)*UJ(3) 344 : 617 CONTINUE 345 : SEGSUP XCOR2 346 : 8 CONTINUE 347 : 618 CONTINUE 348 : SEGSUP XCORD 349 : RETURN 350 : 1100 CONTINUE 351 : DO 1110 I=1,IDIM 352 : DO 1110 J=1,ITE 353 : XPROJ(I,J)=XCORD(I,J) 354 : 1110 CONTINUE 355 : IF (LABCO2.EQ.0) GOTO 1111 356 : IF (LABCO2(3,ICLE).EQ.0) GOTO 1111 357 : KABCO2=LABCO2(1,ICLE) 358 : NVEC=KABCO2(/2) 359 : SEGINI KXPRO2 360 : LABCO2(2,ICLE)=KXPRO2 361 : IF (NVEC.EQ.0) GOTO 1111 362 : DO 9 IVEC=1,NVEC 363 : XCOR2=KABCO2(1,IVEC) 364 : ICOR2=KABCO2(2,IVEC) 365 : SEGINI XPRO2 366 : KXPRO2(IVEC)=XPRO2 367 : DO 1113 J=1,ITE 368 : IF (ICOR2(J).EQ.0) GOTO 1113 369 : DO 1112 I=1,IDIM 370 : XPRO2(I,J)=XCOR2(I,J) 371 : 1112 CONTINUE 372 : 1113 CONTINUE 373 : SEGSUP XCOR2 374 : 9 CONTINUE 375 : 1111 CONTINUE 376 : SEGSUP XCORD 377 : * axez pour tourner correctement avec opengl 378 : axez(1)=0 379 : axez(2)=uj(3) 380 : axez(3)=sqrt(1-uj(3)**2) 381 : if (xn(3).lt.0) axez(3)=-axez(3) 382 : * write (6,*) ' axez ',axez(1),axez(2),axez(3) 383 : 384 : RETURN 385 : END 386 : 387 : 388 : 389 : 390 :
© Cast3M 2003 - All rights reserved.
Disclaimer