Télécharger cadrcl.eso

Retour à la liste

Numérotation des lignes :

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

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