1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
|
*DECK SYB4TC
SUBROUTINE SYB4TC (DELTAR,DDELTA,ANGLES,NHMAX,IXRAYO,IS1,NSECT4,
& IFAC,NUMREG,RAYONS,ZZW,ZZE,ZZR,HXRAYO)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Compute the track lengths and interception lengths in a rectangular
* cell.
*
*Copyright:
* Copyright (C) 2002 Ecole Polytechnique de Montreal
* This library is free software; you can redistribute it and/or
* modify it under the terms of the GNU Lesser General Public
* License as published by the Free Software Foundation; either
* version 2.1 of the License, or (at your option) any later version
*
*Author(s): A. Hebert
*
*Parameters: input
* DELTAR used to compute the interception.
* DDELTA undefined.
* ANGLES angular values.
* NHMAX number of interceptions.
* IXRAYO tube indices.
* IS1 index of the first sector.
* NSECT4 number of sectors.
* IFAC undefined.
* NUMREG region indices of the tube sectors.
* RAYONS radius.
* ZZW position of the west interception (left).
* ZZE position of the east interception (right).
*
*Parameters: input/output
* ZZR tracking information.
* HXRAYO preceding/next interception lengths.
*
*-----------------------------------------------------------------------
*
*----
* SUBROUTINE ARGUMENTS
*----
INTEGER NHMAX,IXRAYO(NHMAX),IS1,NSECT4,IFAC,NUMREG(NSECT4,*)
REAL DELTAR,DDELTA,ANGLES(NSECT4),RAYONS(*),ZZW,ZZE,ZZR(*),
& HXRAYO(NHMAX+1)
*----
* LOCAL VARIABLES
*----
*:AE Aire de la courbure Est (a Droite)
*:AW Aire de la courbure West (a Gauche)
*:GE Abcisse du point Est (Trajectoire precedente)
*:GW Abcisse du point West (Trajectoire precedente)
*:HE Abcisse du point Est
*:HW Abcisse du point West (Intersection)
*:XW Abcisse du point West (Region)
*:IHE No de l'abcisse Est
*:IRE No du Rayon Cote Est
C
*:IL No de l'Intersection
*:IREGC No de la Region Courante
*:IREGS No de la Region Suivante
*:ISC No du Secteur Courant
*:IS2 No du Secteur Suivant
*:JRC No de la Couronne Courante
*:JRS No de la Couronne Suivante
*
DDELT2 = DDELTA * DDELTA
DELTA2 = DELTAR * DELTAR
*
AW = 0.
GW = HXRAYO(1)
HW = ZZW
XW = HW
JRC = IXRAYO(1)
ISC = IS1
ISF = ISC
IF (NSECT4 .GT. 1) THEN
IF (IFAC .EQ. 1) THEN
ISF = 1 - ISF
ELSEIF (IFAC .EQ. 2) THEN
ISF = ISF - NSECT4 / 2
ELSEIF (IFAC .EQ. 3) THEN
ISF = 1 + NSECT4 / 2 - ISF
ENDIF
IF (ISF .LE. 0) ISF = ISF + NSECT4
ENDIF
IREGC = NUMREG(ISF, JRC)
IL = 0
*
HXRAYO(1) = HW
*----
* Boucle des Volumes internes
* Debut
*----
DO IHE = 2, NHMAX
*
* No de Couronne Suivante
JRS = IXRAYO(IHE)
*
* Soit : Meme Couronne => Changement de Secteur
* Sinon : Changement de Couronne => No de Rayon
IRE = 0
IF (JRC .EQ. JRS) THEN
IF (ISC .LT. NSECT4) THEN
IS2 = ISC + 1
ELSE
IS2 = 1
ENDIF
ELSE
IS2 = ISC
IRE = MIN(JRC, JRS)
ENDIF
*
* Soit Changement de Secteur
IF (JRC .EQ. JRS) THEN
HE = DELTAR * ANGLES(ISC)
*
* Sinon Changement de Couronne
ELSE
H2 = RAYONS(IRE) * RAYONS(IRE) - DELTA2
IF (H2 .GT. 0) THEN
HE = SQRT(H2)
IF (JRS .EQ. IRE) THEN
HE = - HE
ENDIF
ELSE
HE = 0.
ENDIF
*
ENDIF
*
* Protection contre les longueurs negatives
IF (HE .LT. HW) THEN
HE = HW
ENDIF
*
* +++
* Debut
* Region Suivante
ISF = IS2
IF (NSECT4 .GT. 1) THEN
IF (IFAC .EQ. 1) THEN
ISF = 1 - ISF
ELSEIF (IFAC .EQ. 2) THEN
ISF = ISF - NSECT4 / 2
ELSEIF (IFAC .EQ. 3) THEN
ISF = 1 + NSECT4 / 2 - ISF
ENDIF
IF (ISF .LE. 0) ISF = ISF + NSECT4
ENDIF
IREGS = NUMREG(ISF, JRS)
IF (IREGS .NE. IREGC) THEN
*
* Nouvelles Valeurs
AE = 0.
GE = HXRAYO(IHE)
*
* Moyenne entre l'intervalle precedent et Actuel
ZZH = HE - XW
ZZH = (ZZH + GE - GW) * 0.5
*
* Ajout de la courbure
IF (JRC .NE. JRS) THEN
H2 = HE - GE
H2 = H2 * H2
XCORDE = H2 + DDELT2
IF (XCORDE .GT. 0.) THEN
XCORDE = SQRT(XCORDE)
XUNITE = XCORDE / RAYONS(IRE)
XUNITE = XUNITE / 2.
XALPHA = ASIN(XUNITE)
XUNITE = XALPHA - COS(XALPHA) * XUNITE
AE = XUNITE * RAYONS(IRE) * RAYONS(IRE)
AE = AE / DDELTA
ENDIF
*
IF (JRS .EQ. IRE) THEN
AE = - AE
ENDIF
ENDIF
*
* Longueur Moyenne
* Nouvelle Abcisse
IL = IL + 1
ZZR(IL) = ZZH + AE - AW
C
C Valeurs Precedentes (Region)
AW = AE
GW = GE
XW = HE
ENDIF
* Fin Region Suivante
* +++
*
* Valeurs Precedentes (Intersection)
HXRAYO(IHE) = HE
HW = HE
*
* Suivants
IREGC = IREGS
ISC = IS2
JRC = JRS
*
* - - - - - - - - - - - - -
* Boucle des Volumes internes
* Fin
* - - - - - - - - - - - - -
ENDDO
*
* - - - - - - - - - - - - -
* Dernier Intervalle
* Debut
* - - - - - - - - - - - - -
*
* Intervalle
ZZH = ZZE - XW
*
* Protection contre les Volumes Negatifs
ZZH = MAX(ZZH, 0.)
*
* Moyenne entre l'intervalle precedent et Actuel
ZZH = (ZZH + HXRAYO(NHMAX+1) - GW) * 0.5
*
* Ajout de la Courbure Eventuelle
* Mise a jour de la Nouvelle abcisse
IL = IL + 1
ZZR(IL) = ZZH - AW
HXRAYO(NHMAX+1) = ZZE
* - - - - - - - - - - - - -
* Dernier Intervalle
* Fin
* - - - - - - - - - - - - -
*
RETURN
END
|