1 |
PROGRAM DIREZIO |
2 |
C |
3 |
C-O---------------------------------------------------------------------O |
4 |
C-O -O |
5 |
C-O PURPOSE AND METHODS : RICOSTRUIRE EVENTO PER EVENTO LA -O |
6 |
C-O DIREZIONE D'ARRIVO DELLA PARTICELLA (PUNTAMENTO) -O |
7 |
C-O -O |
8 |
C-O INPUTS : DATA FILE, OUTPUT FILE, CONTAINEMENT OR NOT -O |
9 |
C-O OUTPUTS : RZ FILE FOR PAW -O |
10 |
C-O CONTROLS: -O |
11 |
C-O -O |
12 |
C-O -O |
13 |
C-O CREATED FEB-2000 BY EMILIANO MOCCHIUTTI -O |
14 |
C-O -O |
15 |
C-O---------------------------------------------------------------------O |
16 |
C |
17 |
C PARAMETER (NUMEV=100000, |
18 |
PARAMETER (NVAR=1, TR=5.) !5. |
19 |
PARAMETER (NPLA=22,NCHA=96,LENSEV=NPLA*NCHA) |
20 |
PARAMETER (SOG=60.) |
21 |
PARAMETER (NWPAWC=3000000) |
22 |
PARAMETER (NWORD=100000) |
23 |
PARAMETER (NPAR=2, NPAR2=2, NPL=44) |
24 |
PARAMETER (PI=3.14159265358979324) |
25 |
PARAMETER (CALIB=0.0001059994) |
26 |
C |
27 |
INTEGER LREC |
28 |
PARAMETER (LREC=6144) |
29 |
C |
30 |
INTEGER EVCON, EVCOINC, EVCONMIS |
31 |
INTEGER TRIG |
32 |
INTEGER POSTRIX, POSTRIY |
33 |
INTEGER NPU, FLA, FLA2 |
34 |
INTEGER NDEP |
35 |
INTEGER EVENTO |
36 |
INTEGER NLK, NLK1, NLK2 |
37 |
INTEGER F |
38 |
INTEGER FINPAR1, FINPAR |
39 |
INTEGER IFIN2, NPLL |
40 |
INTEGER NSTR |
41 |
INTEGER NINF, NSUP |
42 |
INTEGER NLL |
43 |
INTEGER NPLOTX, NPLOTY, NPLOTX2, NPLOTY2, NNPLO |
44 |
C |
45 |
REAL TRASX, TRASY |
46 |
REAL POSIX, POSIY, POSIXM, POSIYM |
47 |
REAL QQX, QQY, QQZ, QQXM, QQYM, QQZM |
48 |
REAL VALORI(NWORD) |
49 |
REAL ZXY, PX, PY, THET, PH |
50 |
REAL X(NPLA), X2(NPLA), Y, Y2, EY, RN |
51 |
REAL XV(NPLA), YV, EYV(NPLA), PMI, EXV(NPLA) |
52 |
REAL EX |
53 |
REAL YN |
54 |
REAL ZZ |
55 |
REAL YM, EYM |
56 |
REAL THX, THY, THXM ,THYM |
57 |
REAL PPY, PPX |
58 |
REAL VENET, ENET, ERENET |
59 |
REAL ER1, ER2, ERENETM, ENETM |
60 |
REAL ENM, ERENM |
61 |
REAL TMIS, PMIS |
62 |
REAL ZETA, ZETAM, ZETA2, ZETA3 |
63 |
REAL PARZEN, PARZ, PARZE, ERZEN, ER3, ER4 |
64 |
REAL LIMLA, LIMLA2, LIMLA3, LIMLA4, LIMLA5 |
65 |
REAL FFLA |
66 |
REAL PMETX, PMETY |
67 |
REAL ENMEZX, ENMEZY |
68 |
REAL DELTX, DELTY |
69 |
REAL PO3S, FFLA2 |
70 |
REAL XEN(NPL), YEN(NPL), EYEN(NPL) |
71 |
REAL XEN2(NPL), YEN2(NPL), EYEN2(NPL) |
72 |
REAL POSMAX, PAR1 |
73 |
REAL QPLX, QPLY |
74 |
REAL X0ATT, FFLACO3, TOTTO2 |
75 |
REAL SCARTA, BANDI |
76 |
REAL VERT, ENCEN |
77 |
REAL ESTMA |
78 |
REAL NUMEV |
79 |
C |
80 |
DOUBLE PRECISION THETA, CPHI, SPHI, PHI |
81 |
DOUBLE PRECISION POSIXMD, POSIYMD, PMISD, TMISD, PMID |
82 |
DOUBLE PRECISION POSIXD, POSIYD |
83 |
DOUBLE PRECISION MX, QX, MY, QY, QXP, QYP |
84 |
DOUBLE PRECISION THETA2, POSIX2, POSIY2 |
85 |
DOUBLE PRECISION MMX, MMY, TTHX, TTHY |
86 |
C |
87 |
CHARACTER*6 CHOPT |
88 |
CHARACTER*20 FILENAME,NAME |
89 |
CHARACTER*80 CHTITLE |
90 |
CHARACTER*1 ALFA, BETA |
91 |
C |
92 |
EXTERNAL FUN |
93 |
EXTERNAL ENERGIF |
94 |
C |
95 |
COMMON/ENSTRIP/ZXY(2,LENSEV), YN(NCHA) |
96 |
C |
97 |
COMMON/VETFIT/YV(NPLA) |
98 |
C |
99 |
COMMON/DATI/THET, PH, TMIS, PMIS, THX, THY, THXM, THYM, |
100 |
& FLA2, FLA, QQX, QQY, QQZ, QQXM, QQYM, QQZM, |
101 |
& ZZ(2,NPLA), TRASX(NPLA), TRASY(NPLA), EVENTO, |
102 |
& Y(NPLA), EX(NPLA), Y2(NPLA), EY(NPLA), |
103 |
& VENET, PARZ(2,NPLA), PARZEN, FFLA, |
104 |
& ENMEZX(2,NPLA),ENMEZY(2,NPLA), |
105 |
& DELTX(NPLA), DELTY(NPLA), |
106 |
& FFLA2, POSMAX, NSTR(4,NPLA), PAR1, |
107 |
& QPLX(6,NPLA),QPLY(6,NPLA),VERT,ENCEN(12,NPLA), |
108 |
& NNPLO(2,NPLA),ESTMA(3) |
109 |
C |
110 |
C COMMON/DATI2/TOTTO2,SCARTA |
111 |
C |
112 |
COMMON/HCFITD/PA(NPAR),SIGPAR(NPAR),CHI2,PMIN(NPAR),PMAX(NPAR), |
113 |
& STEP(NPAR), |
114 |
& PAR(NPAR2),SIGPA(NPAR2),CHI,PMI2(NPAR2),PMA(NPAR2), |
115 |
& STE(NPAR2) |
116 |
COMMON/QUEST/IQUEST(100) |
117 |
COMMON/PAWC/HMEMOR(NWPAWC) |
118 |
C |
119 |
C ################################################################# |
120 |
C |
121 |
CALL HLIMIT(NWPAWC) |
122 |
CALL HDELET(0) |
123 |
C |
124 |
PRINT *,'FILENAME TO SAVE?' |
125 |
READ (*,5050)NAME |
126 |
PRINT *,'NUMERO DI EVENTI?' |
127 |
READ (*,6006),NUMEV |
128 |
IQUEST(10)=300000 |
129 |
CALL HROPEN(13,'SIMULATION', NAME,'N',LREC,ISTAT) |
130 |
IF (ISTAT.NE.0) GO TO 1000 |
131 |
C |
132 |
CALL HBNT(1,'DATI',' ') |
133 |
CALL HBSET('BSIZE',LREC,IERR) |
134 |
CALL HBNAME(1,'DATI',THET, |
135 |
& 'THET:R,PH:R,TMIS:R,PMIS:R,THX:R,THY:R,THXM:R,THYM:R,FLA2:I, |
136 |
& FLA:I,QQX:R,QQY:R,QQZ:R,QQXM:R,QQYM:R,QQZM:R,ZZ(2,22):R, |
137 |
& TRASX(22):R,TRASY(22):R,EVENTO:I,X(22):R,EX(22):R,Y(22):R, |
138 |
& EY(22):R,VENET:R,PARZ(2,22):R,PARZEN:R,FFLA:R,ENMEZX(2,22):R, |
139 |
& ENMEZY(2,22):R,DELTX(22):R,DELTY(22):R,FFLA2:R,POSMAX:R, |
140 |
& NSTR(4,22):I,PAR1:R,QPLX(6,22):R,QPLY(6,22):R,VERT:R, |
141 |
& ENCEN(12,22):R,NNPLO(2,22):I,ESTMA(3):R') |
142 |
C |
143 |
C CALL HBNT(2,'DATI2',' ') |
144 |
C CALL HBSET('BSIZE',LREC,IERR) |
145 |
C CALL HBNAME(2,'DATI2',TOTTO2,'TOTTO2:R,SCARTA:R') |
146 |
C |
147 |
PRINT *,'CONTAINED OR NOT ? (C/N)' |
148 |
C$ READ (*,6001),ALFA |
149 |
ALFA = 'C' |
150 |
C |
151 |
PRINT *,'PROFILO LONGITUDINALE? (Y/N)' |
152 |
C$ READ (*,6001),BETA |
153 |
BETA = 'Y' |
154 |
C |
155 |
IF (ALFA.EQ.'C'.OR.ALFA.EQ.'c') THEN |
156 |
PRINT *,'DA CHE PIANO VUOI CHE ENTRINO LE PARTICELLE?' |
157 |
C$ READ (*,6005),LIMLA2 |
158 |
LIMLA2 = 0. |
159 |
PRINT *,'PER QUANTI PIANI?' |
160 |
C$ READ (*,6005),LIMLA5 |
161 |
LIMLA5 = 8. |
162 |
PRINT *,'DA CHE PIANO VUOI CHE ESCANO?' |
163 |
C$ READ (*,6005),LIMLA4 |
164 |
LIMLA4 = 12. |
165 |
PRINT *,'LIMITAZIONE SULLA PROIEZIONE SUL PIANO 22: ' |
166 |
C$ READ (*,6006),LIMLA3 |
167 |
LIMLA3 = 999. |
168 |
LIMLA4 = 22. - LIMLA4 |
169 |
LIMLA = 12.1+LIMLA4/TAN((18.964-ABS(LIMLA2+LIMLA4)*.862)/24.2) |
170 |
PRINT *,'LIMLA CALCOLATO = ',LIMLA |
171 |
IF (LIMLA.GT.LIMLA3) LIMLA = LIMLA3 |
172 |
PRINT *,'LIMLA = ',LIMLA |
173 |
ELSE |
174 |
LIMLA = 10000. |
175 |
LIMLA2 = 0. |
176 |
LIMLA5 = 11. |
177 |
LIMLA4 = 11. |
178 |
LIMLA3 = 999. |
179 |
ENDIF |
180 |
C |
181 |
TRIG = 0 |
182 |
EVCONMIS = 0 |
183 |
EVCON = 0 |
184 |
EVCOINC = 0 |
185 |
EVENTO = 0 |
186 |
ENETM = 0. |
187 |
ERENETM =0. |
188 |
ER1 = 0. |
189 |
ER2 = 0. |
190 |
ER3 = 0. |
191 |
ER4 = 0. |
192 |
TOTTO2 = 0. |
193 |
SCARTA = 0. |
194 |
C |
195 |
CALL VZERO(YN,NCHA) |
196 |
CALL VZERO(X,NPLA) |
197 |
CALL VZERO(X2,NPLA) |
198 |
C |
199 |
C CONVERSIONE STRIP-CM |
200 |
C |
201 |
DO NN = 1, NCHA |
202 |
YN(NN) = (NN - 1.) * .244 +.218 !.2448 + .185 |
203 |
IF (NN.GT.32) YN(NN) = YN(NN) + .292 |
204 |
IF (NN.GT.64) YN(NN) = YN(NN) + .292 |
205 |
ENDDO |
206 |
C |
207 |
C CONVERSIONE XVIEW-CM, YVIEW-CM |
208 |
C |
209 |
DO NN = 1, NPLA |
210 |
RN = FLOAT(NN) |
211 |
X(NN) = (RN - 1.)*.862 + .169 |
212 |
X2(NN) = (RN - 1.)*.862 + .843 |
213 |
ENDDO |
214 |
C |
215 |
PRINT *,'DATA FILENAME : ' |
216 |
READ (*,5050)FILENAME |
217 |
C |
218 |
OPEN (UNIT=83,FILE=FILENAME,STATUS='OLD',FORM='UNFORMATTED') |
219 |
C |
220 |
DO JJ=1, NUMEV |
221 |
C |
222 |
READ (83,END=900),LENGHT |
223 |
CALL VZERO(ZXY,2*LENSEV) |
224 |
CALL VZERO(VALORI,NWORD) |
225 |
CALL VZERO(TRASX,22) |
226 |
CALL VZERO(TRASY,22) |
227 |
NEVENT = NEVENT + 1 |
228 |
CALL LEGGE(VALORI,LENGHT) |
229 |
C |
230 |
C VALORI IS A VECTOR OF VARIABLE LENGHT (LENGHT) WHICH CONTAINS ALL THE NEEDED |
231 |
C INFORMATION FROM THE SIMULATION |
232 |
C |
233 |
IF (LENGHT.LT.17) GO TO 449 |
234 |
C |
235 |
C PUT IN THE VECTOR ZXY ALL THE ENERGIES RELEASED BY THE |
236 |
C PARTICLE IN EVERY STRIP |
237 |
C |
238 |
DO I = 6,LENGHT-17,2 |
239 |
KK = VALORI(I) |
240 |
IF (KK.LT.5000) THEN |
241 |
ZXY(1,KK) = VALORI(I+1) |
242 |
C |
243 |
C SATURAZIONE DELLE STRIP A 1100MIP |
244 |
C |
245 |
IF (ZXY(1,KK).GT.1100.) ZXY(1,KK) = 1100. |
246 |
ELSE |
247 |
KS = KK - 5000 |
248 |
ZXY(2,KS) = VALORI(I+1) |
249 |
IF (ZXY(2,KS).GT.1100.) ZXY(2,KS) = 1100. |
250 |
ENDIF |
251 |
ENDDO |
252 |
EVENTO = EVENTO + 1 |
253 |
C |
254 |
C RICAVA GLI ANGOLI VERI DI INCIDENZA... |
255 |
C |
256 |
ZETA = ABS(VALORI(LENGHT)) |
257 |
THETA = DBLE(VALORI(LENGHT-10)) |
258 |
IF (THETA.EQ.0.) THEN |
259 |
CPHI = 0. |
260 |
SPHI = 1. |
261 |
ELSE |
262 |
CPHI = DBLE(VALORI(3))/DSIND(THETA) |
263 |
SPHI = DBLE(VALORI(4))/DSIND(THETA) |
264 |
ENDIF |
265 |
IF (CPHI.GT.1.) CPHI = 1. |
266 |
IF (CPHI.LT.-1.) CPHI = -1. |
267 |
IF (SPHI.GE.0.) PHI = DACOSD(CPHI) |
268 |
IF (SPHI.LT.0.) PHI =-DACOSD(CPHI)+DBLE(360.) |
269 |
C |
270 |
C SE LA PARTICELLA ARRIVA DALLA FACCIA X LATERALE GLI ANGOLI CHE CI DA` |
271 |
C LA SIMULAZIONE SONO RELATIVI AD UN SIS DI RIF RUOTATO E LI SI DEVONO |
272 |
C CONVERTIRE NEGLI ANGOLI GIUSTI |
273 |
C |
274 |
IF (ZETA.NE.18.964) THEN |
275 |
IF (PHI.LT.90.) THEN |
276 |
TTHY = PHI |
277 |
ELSE |
278 |
TTHY = PHI - DBLE(360.) |
279 |
ENDIF |
280 |
IF (PHI.LT.90.OR.PHI.GT.270.) THEN |
281 |
TTHX = DATAND(DABS(DTAND(THETA))/(DSQRT( |
282 |
& DBLE(1.)+DTAND(PHI)**2.)))-DBLE(90.) |
283 |
ENDIF |
284 |
IF (PHI.GT.90..AND.PHI.LT.270.) THEN |
285 |
TTHX = DATAND(-DABS(DTAND(THETA))/(DSQRT( |
286 |
& DBLE(1.)+DTAND(PHI)**2.)))-DBLE(90.) |
287 |
ENDIF |
288 |
MMX = DTAND(TTHX) |
289 |
MMY = DTAND(TTHY) |
290 |
THETA = DATAND(DSQRT((MMX*MMX) + (MMY*MMY))) |
291 |
C |
292 |
IF (MMX.EQ.0..AND.MMY.GT.0.) PHI = DBLE(90.) |
293 |
IF (MMX.EQ.0..AND.MMY.LT.0.) PHI = DBLE(270.) |
294 |
IF (MMY.EQ.0..AND.MMX.GE.0.) PHI = DBLE(0.) |
295 |
IF (MMY.EQ.0..AND.MMX.LT.0.) PHI = DBLE(180.) |
296 |
C |
297 |
IF (MMY.NE.0..AND.MMX.NE.0.) THEN |
298 |
PHI = DATAND(MMY/MMX) |
299 |
IF (MMY.LT.0..AND.MMX.GT.0.) PHI = PHI + DBLE(360.) |
300 |
IF (MMX.LT.0.) PHI = PHI + DBLE(180.) |
301 |
ENDIF |
302 |
CPHI = DCOSD(PHI) |
303 |
SPHI = DSIND(PHI) |
304 |
ENDIF |
305 |
C |
306 |
C ...E LA PROIEZIONE DELLA TRAIETTORIA SULL'ULTIMO PIANO |
307 |
C |
308 |
POSIXD=DBLE(VALORI(LENGHT-13))+DBLE(ZETA)* |
309 |
& DTAND(THETA)*CPHI |
310 |
POSIYD=DBLE(VALORI(LENGHT-12))+DBLE(ZETA)* |
311 |
& DTAND(THETA)*SPHI |
312 |
C |
313 |
C CONTAINED OR NOT? |
314 |
C |
315 |
IF (ALFA.EQ.'C'.OR.ALFA.EQ.'c') THEN |
316 |
IF (POSIXD.GT.LIMLA.OR.POSIXD.LT.-LIMLA.OR.POSIYD.GT. |
317 |
& LIMLA.OR.POSIYD.LT.-LIMLA) GO TO 449 |
318 |
ZETA2 = 18.964 - LIMLA2*.862 |
319 |
ZETA3 = 18.964 - (LIMLA2+LIMLA5)*.862 |
320 |
IF (ZETA.GT.ZETA2.OR.ZETA.LT.ZETA3) GO TO 449 |
321 |
ENDIF |
322 |
C |
323 |
C EVENTI REALMENTE CONTENUTI |
324 |
C |
325 |
IF (POSIXD.LE.LIMLA.AND.POSIXD.GE.-LIMLA.AND.POSIYD.LE. |
326 |
& LIMLA.AND.POSIYD.GE.-LIMLA) EVCON = EVCON + 1 |
327 |
C### |
328 |
c BANDI = 0. |
329 |
c QQZ = ZETA |
330 |
c PX = SNGL(POSIXD) |
331 |
c PY = SNGL(POSIYD) |
332 |
c THET = SNGL(THETA) |
333 |
c PH = SNGL(PHI) |
334 |
c QQX = SNGL(VALORI(LENGHT-13)) |
335 |
c IF ((QQZ.GE.12.964.AND.QQZ.LE.18.964.AND. |
336 |
c & (ABS(PX).GE.0.AND.ABS(PX).LT.(12.1+1.35*4)).AND. |
337 |
c & (ABS(PY).GE.0.AND.ABS(PY).LT.(12.1+1.35*4)).AND. |
338 |
c & (THET.GE.0.AND.THET.LT.63.5))) THEN |
339 |
cC |
340 |
c FFLACO3=(ABS(QQZ)+ABS((12.1-QQX)-18.964)/(TAND(THET)*COSD(PH))) |
341 |
cC |
342 |
c X0ATT = .76*(FFLACO3/(COSD(THET))) |
343 |
cC |
344 |
c IF (X0ATT.GT.10.) THEN |
345 |
c TOTTO2 = TOTTO2 + 1. |
346 |
c BANDI = 1. |
347 |
c ENDIF |
348 |
cC |
349 |
c ENDIF |
350 |
ccc IF (BANDI.NE.1.) GO TO 449 |
351 |
C### |
352 |
C TROVA IL PIANO SU CUI E` STATA RILASCIATA LA MASSIMA ENERGIA |
353 |
C (FINPAR) E REGISTRA LE ENERGIE RILASCIATE SU CIASCUN PIANO (PARZ) |
354 |
C |
355 |
ENET = 0. |
356 |
PARZE = 0. |
357 |
FINPAR = 1 |
358 |
FINPAR1 = 1 |
359 |
CALL VZERO(PARZ,2*NPLA) |
360 |
CALL VZERO(ESTMA,3) |
361 |
DO NN = 1, NPLA |
362 |
DO KK = 1, NCHA |
363 |
NLK = (KK - 1)*NPLA + NN |
364 |
IF (ZXY(1,NLK).GE.ESTMA(1).AND.ZXY(1,NLK).LT.1100.) THEN |
365 |
ESTMA(1) = ZXY(1,NLK) |
366 |
ESTMA(2) = NN |
367 |
ESTMA(3) = KK |
368 |
ENDIF |
369 |
IF (ZXY(2,NLK).GE.ESTMA(1).AND.ZXY(2,NLK).LT.1100.) THEN |
370 |
ESTMA(1) = ZXY(2,NLK) |
371 |
ESTMA(2) = NN |
372 |
ESTMA(3) = KK |
373 |
ENDIF |
374 |
PARZ(1,NN) = PARZ(1,NN) + ZXY(1,NLK) |
375 |
PARZ(2,NN) = PARZ(2,NN) + ZXY(2,NLK) |
376 |
ENET = ENET + ZXY(1,NLK) + ZXY(2,NLK) |
377 |
ENDDO |
378 |
IF (PARZ(1,NN).GE.PARZE) THEN |
379 |
PARZE = PARZ(1,NN) |
380 |
FINPAR = NN |
381 |
FINPAR1 = 1 |
382 |
ENDIF |
383 |
IF (PARZ(2,NN).GE.PARZE) THEN |
384 |
PARZE = PARZ(2,NN) |
385 |
FINPAR = NN |
386 |
FINPAR1 = 2 |
387 |
ENDIF |
388 |
ENDDO |
389 |
cc FINPAR = FINPAR + 3 |
390 |
cc IF (FINPAR.GT.NPLA) FINPAR = NPLA |
391 |
C$$$$ |
392 |
C |
393 |
C CERCA IL CLUSTER DI TRE STRIP CHE HA VISTO LA MAGGIORE ENERGIA |
394 |
C PER CIASCUN PIANO |
395 |
C |
396 |
PO3S = 0. |
397 |
DO I = 1, NPLA |
398 |
STRIMAX = 0. |
399 |
STRIMAY = 0. |
400 |
DO M = 1, 94 |
401 |
NLK = (M - 1) * NPLA + I |
402 |
NLK1 = M * NPLA + I |
403 |
NLK2 = (M + 1) * NPLA + I |
404 |
IF (ZXY(1,NLK).GT.TR.AND.ZXY(1,NLK1).GT.TR |
405 |
& .AND.ZXY(1,NLK2).GT.TR) THEN |
406 |
PPX = ZXY(1,NLK) + ZXY(1,NLK1) + ZXY(1,NLK2) |
407 |
ELSE |
408 |
PPX = 0. |
409 |
ENDIF |
410 |
IF (ZXY(2,NLK).GT.TR.AND.ZXY(2,NLK1).GT.TR |
411 |
& .AND.ZXY(2,NLK2).GT.TR) THEN |
412 |
PPY = ZXY(2,NLK) + ZXY(2,NLK1) + ZXY(2,NLK2) |
413 |
ELSE |
414 |
PPY = 0. |
415 |
ENDIF |
416 |
IF (PPX.GT.STRIMAX) THEN |
417 |
STRIMAX = PPX |
418 |
POSTRIX = M + 1 |
419 |
ENDIF |
420 |
IF (PPY.GT.STRIMAY) THEN |
421 |
STRIMAY = PPY |
422 |
POSTRIY = M + 1 |
423 |
ENDIF |
424 |
ENDDO |
425 |
C |
426 |
C REGISTRA LA POSIZIONE DELLA STRIP CENTRALE |
427 |
C |
428 |
IF (STRIMAX.GT.0.) THEN |
429 |
TRASX(I) = REAL(POSTRIX) + .122 |
430 |
ELSE |
431 |
TRASX(I) = -10. |
432 |
ENDIF |
433 |
IF (STRIMAY.GT.0.) THEN |
434 |
TRASY(I) = REAL(POSTRIY) + .122 |
435 |
ELSE |
436 |
TRASY(I) = -10. |
437 |
ENDIF |
438 |
C |
439 |
IF (STRIMAX.GT.PO3S) THEN |
440 |
PO3S = STRIMAX |
441 |
FFLA2 = FLOAT(I) |
442 |
ENDIF |
443 |
C |
444 |
ENDDO |
445 |
C |
446 |
C LA PRIMA VOLTA IL FIT VIENE ESEGUITO SULLE STRIP CENTRALI |
447 |
C |
448 |
DO F = 1, 2 |
449 |
cc F = 1 |
450 |
C |
451 |
C PREPARA I VETTORI PER IL SECONDO FIT; ADESSO LA STRIP "CENTRALE" |
452 |
C DIVENTA QUELLA DOVE PASSA LA TRAIETTORIA DELLA PRIMA RETTA |
453 |
C |
454 |
IF (F.EQ.2) THEN |
455 |
DO NN = 1, NPLA |
456 |
PY = X(NN)*SNGL(MX) + SNGL(QX) |
457 |
CALL TRAIE(PY,TRASX(NN)) |
458 |
C |
459 |
PY = X2(NN)*SNGL(MY) + SNGL(QY) |
460 |
CALL TRAIE(PY,TRASY(NN)) |
461 |
C |
462 |
C CORR. |
463 |
C |
464 |
IF (NN.GT.(FINPAR+1)) THEN |
465 |
IF (MX.LT.-0.58) |
466 |
& TRASX(NN) = TRASX(NN) - .4 |
467 |
IF (MX.GT.0.58) |
468 |
& TRASX(NN) = TRASX(NN) + .4 |
469 |
IF (MY.LT.-0.58) |
470 |
& TRASY(NN) = TRASY(NN) - .4 |
471 |
IF (MY.GT.0.58) |
472 |
& TRASY(NN) = TRASY(NN) + .4 |
473 |
ENDIF |
474 |
IF (NN.LT.(FINPAR-1)) THEN |
475 |
IF (MX.LT.-0.58) |
476 |
& TRASX(NN) = TRASX(NN) + .2 |
477 |
IF (MX.GT.0.58) |
478 |
& TRASX(NN) = TRASX(NN) - .2 |
479 |
IF (MY.LT.-0.58) |
480 |
& TRASY(NN) = TRASY(NN) + .2 |
481 |
IF (MY.GT.0.58) |
482 |
& TRASY(NN) = TRASY(NN) - .2 |
483 |
ENDIF |
484 |
C |
485 |
ENDDO |
486 |
ENDIF |
487 |
C |
488 |
C NUMERO DI STRIP SOPRA LA SOGLIA SOG ATTORNO (48 A SX E 48 A DX) |
489 |
C ALLA POSIZIONE DELLA TRAIETTORIA |
490 |
C |
491 |
IF (F.EQ.2) THEN |
492 |
CALL VZERO(NSTR,4*NPLA) |
493 |
CALL VZERO(ENCEN,12*NPLA) |
494 |
CALL VZERO(NNPLO,2*NPLA) |
495 |
VERT = -1. |
496 |
DO NN = 1, NPLA |
497 |
NPLOTX = 0 |
498 |
NPLOTY = 0 |
499 |
NPLOTX2 = 0 |
500 |
NPLOTY2 = 0 |
501 |
DO NJ = 1, NCHA |
502 |
NLK = (NJ - 1.)*NPLA + NN |
503 |
IF (ZXY(1,NLK).GT.0.) NSTR(3,NN) = NSTR(3,NN)+1 |
504 |
IF (ZXY(2,NLK).GT.0.) NSTR(4,NN) = NSTR(4,NN)+1 |
505 |
IF (ZXY(1,NLK).GT.0.) THEN |
506 |
NPLOTX = NPLOTX + 1 |
507 |
ELSE |
508 |
NPLOTX = 0 |
509 |
ENDIF |
510 |
IF (ZXY(2,NLK).GT.0.) THEN |
511 |
NPLOTY = NPLOTY + 1 |
512 |
ELSE |
513 |
NPLOTY = 0 |
514 |
ENDIF |
515 |
IF (NPLOTX.GT.NNPLO(1,NN)) NNPLO(1,NN) = NPLOTX |
516 |
IF (NPLOTY.GT.NNPLO(2,NN)) NNPLO(2,NN) = NPLOTY |
517 |
IF (NPLOTX.GE.2) NPLOTX2 = 1 |
518 |
IF (NPLOTY.GE.3) NPLOTY2 = 1 |
519 |
ENDDO |
520 |
IF (VERT.LT.0..AND.NPLOTX2.EQ.1.AND.NPLOTY2.EQ.1) |
521 |
& VERT = FLOAT(NN) |
522 |
IF (NINT(TRASX(NN)).GT.0.) THEN |
523 |
NINF = NINT(TRASX(NN)) - 48 |
524 |
NSUP = NINT(TRASX(NN)) + 48 |
525 |
IF (NINF.LT.1) NINF = 1 |
526 |
IF (NSUP.GT.96) NSUP = 96 |
527 |
DO NJ = NINF, NSUP |
528 |
NLK = (NJ - 1.)*NPLA + NN |
529 |
IF (ZXY(1,NLK).GT.SOG.AND.ZXY(1,NLK).GT.0.) |
530 |
& NSTR(1,NN) = NSTR(1,NN) + 1 |
531 |
C |
532 |
C 1 RM |
533 |
C |
534 |
IF (NJ.GE.NINT(TRASX(NN)-3.).AND. |
535 |
& NJ.LE.NINT(TRASX(NN)+3.)) THEN |
536 |
ENCEN(1,NN) = ENCEN(1,NN)+ZXY(1,NLK) |
537 |
IF (ZXY(1,NLK).GT.0.) THEN |
538 |
ENCEN(3,NN) = ENCEN(3,NN)+1. |
539 |
ENDIF |
540 |
ENDIF |
541 |
C |
542 |
C 3 RM |
543 |
C |
544 |
IF (NJ.GE.NINT(TRASX(NN)-12.).AND. |
545 |
& NJ.LE.NINT(TRASX(NN)+12.)) THEN |
546 |
ENCEN(5,NN) = ENCEN(5,NN)+ZXY(1,NLK) |
547 |
IF (ZXY(1,NLK).GT.0.) THEN |
548 |
ENCEN(7,NN) = ENCEN(7,NN)+1. |
549 |
ENDIF |
550 |
ENDIF |
551 |
C |
552 |
C 6 RM |
553 |
C |
554 |
IF (NJ.GE.NINT(TRASX(NN)-25.).AND. |
555 |
& NJ.LE.NINT(TRASX(NN)+25.)) THEN |
556 |
ENCEN(9,NN) = ENCEN(9,NN)+ZXY(1,NLK) |
557 |
IF (ZXY(1,NLK).GT.0.) THEN |
558 |
ENCEN(11,NN) = ENCEN(11,NN)+1. |
559 |
ENDIF |
560 |
ENDIF |
561 |
ENDDO |
562 |
ELSE |
563 |
NSTR(1,NN) = -1 |
564 |
C ENCEN(1,NN) = -1. |
565 |
DO NJ = 1, 12, 2 |
566 |
ENCEN(NJ,NN) = 0. |
567 |
ENDDO |
568 |
ENDIF |
569 |
IF (NINT(TRASY(NN)).GT.0.) THEN |
570 |
NINF = NINT(TRASY(NN)) - 48 |
571 |
NSUP = NINT(TRASY(NN)) + 48 |
572 |
IF (NINF.LT.1) NINF = 1 |
573 |
IF (NSUP.GT.96) NSUP = 96 |
574 |
DO NJ = NINF, NSUP |
575 |
NLK = (NJ - 1.)*NPLA + NN |
576 |
IF (ZXY(2,NLK).GT.SOG.AND.ZXY(2,NLK).GT.0.) |
577 |
& NSTR(2,NN) = NSTR(2,NN) + 1 |
578 |
C |
579 |
C 1 RM |
580 |
C |
581 |
IF (NJ.GE.NINT(TRASY(NN)-3.).AND. |
582 |
& NJ.LE.NINT(TRASY(NN)+3.)) THEN |
583 |
ENCEN(2,NN) = ENCEN(2,NN)+ZXY(2,NLK) |
584 |
IF (ZXY(2,NLK).GT.0.) THEN |
585 |
ENCEN(4,NN) = ENCEN(4,NN)+1. |
586 |
ENDIF |
587 |
ENDIF |
588 |
C |
589 |
C 3 RM |
590 |
C |
591 |
IF (NJ.GE.NINT(TRASY(NN)-12.).AND. |
592 |
& NJ.LE.NINT(TRASY(NN)+12.)) THEN |
593 |
ENCEN(6,NN) = ENCEN(6,NN)+ZXY(2,NLK) |
594 |
IF (ZXY(2,NLK).GT.0.) THEN |
595 |
ENCEN(8,NN) = ENCEN(8,NN)+1. |
596 |
ENDIF |
597 |
ENDIF |
598 |
C |
599 |
C 6 RM |
600 |
C |
601 |
IF (NJ.GE.NINT(TRASY(NN)-25.).AND. |
602 |
& NJ.LE.NINT(TRASY(NN)+25.)) THEN |
603 |
ENCEN(10,NN) = ENCEN(10,NN)+ZXY(2,NLK) |
604 |
IF (ZXY(2,NLK).GT.0.) THEN |
605 |
ENCEN(12,NN) = ENCEN(12,NN)+1. |
606 |
ENDIF |
607 |
ENDIF |
608 |
ENDDO |
609 |
ELSE |
610 |
NSTR(2,NN) = -1 |
611 |
C ENCEN(2,NN) = -1. |
612 |
DO NJ = 1, 12, 2 |
613 |
ENCEN(NJ+1,NN) = 0. |
614 |
ENDDO |
615 |
ENDIF |
616 |
ENDDO |
617 |
ENDIF |
618 |
C |
619 |
C CALCOLA IL BARICENTRO, PER CIASCUN PIANO, RISPETTO 7 STRIP |
620 |
C NEL PRIMO FIT, RISPETTO 3 STRIP NEL SECONDO |
621 |
C |
622 |
DO NN = 1, NPLA |
623 |
CALL BARIC(TRASX(NN),Y(NN),EX(NN),NN,F,1) |
624 |
C |
625 |
CALL BARIC(TRASY(NN),Y2(NN),EY(NN),NN,F,2) |
626 |
ENDDO |
627 |
C |
628 |
C --- PER LE VISTE X --- |
629 |
C METTE NEI VETTORI YV XV EXV I DATI DA FITTARE (ESCLUDE GLI ZERI) |
630 |
C |
631 |
CALL VZERO(XV,NPLA) |
632 |
CALL VZERO(YV,NPLA) |
633 |
CALL VZERO(EXV,NPLA) |
634 |
CALL VZERO(DELTX,NPLA) |
635 |
FLA2 = 0 |
636 |
DO KK = 1, NPLA !!finpar |
637 |
IF (Y(KK).GT.0.) THEN |
638 |
FLA2 = FLA2 + 1 |
639 |
XV(FLA2) = X(KK) |
640 |
YV(FLA2) = Y(KK) |
641 |
EXV(FLA2) = EX(KK) |
642 |
PY = X(KK)*TAN(THX) + QQX + 12.1 |
643 |
DELTX(KK) = PY - YV(FLA2) |
644 |
ENDIF |
645 |
ENDDO |
646 |
IF (FLA2.LT.2) THEN |
647 |
c IF (BANDI.EQ.1.) SCARTA = SCARTA + 1. |
648 |
GO TO 449 |
649 |
ENDIF |
650 |
C |
651 |
C RICHIAMA SOLO NEL PRIMO PASSAGGIO UN ANALISI TOPOLOGICA |
652 |
C DELL'EVENTO PER ESCLUDERE I PUNTI CHE CORRONO SUL BORDO |
653 |
C |
654 |
IF (F.EQ.1) THEN |
655 |
CALL BORDO(FLA2) |
656 |
IF (FLA2.LT.2) THEN |
657 |
c IF (BANDI.EQ.1.) SCARTA = SCARTA + 1. |
658 |
GO TO 449 |
659 |
ENDIF |
660 |
ENDIF |
661 |
C |
662 |
C ESEGUE IL FIT DELLA RETTA, LA PRIMA VOLTA A PARTIRE DA UNA RETTA |
663 |
C A PENDENZA 0, LA SECONDA DALLA RETTA DEL PRIMO FIT |
664 |
C |
665 |
CHOPT = 'QE' |
666 |
IF (F.EQ.1) THEN |
667 |
PA(1) = DBLE(12.1) |
668 |
PA(2) = DBLE(0.) |
669 |
ELSE |
670 |
PA(1) = QX |
671 |
PA(2) = MX |
672 |
ENDIF |
673 |
C |
674 |
CALL HFITV(FLA2,NPLA,NVAR,XV,YV,EXV,FUN,CHOPT,NPAR,PA,STEP, |
675 |
& PMIN,PMAX,SIGPAR,CHI2) |
676 |
C |
677 |
C MX COEFF. ANGOLARE, QX QUOTA |
678 |
C |
679 |
MX = PA(2) |
680 |
QX = PA(1) |
681 |
C |
682 |
C REGISTRA I PULL |
683 |
C |
684 |
IF (F.EQ.2) THEN |
685 |
DO I = 1, NPLA |
686 |
ZZ(1,I) = -100. |
687 |
ZZ(2,I) = -100. |
688 |
ENDDO |
689 |
C |
690 |
DO I = 1, FLA2 |
691 |
NDEP = NINT(1.+((XV(I)-.169)/.862)) |
692 |
YVM = XV(I)*SNGL(MX) + SNGL(QX) |
693 |
EYM = SQRT((XV(I)*SNGL(SIGPAR(2)))**2+(SNGL(SIGPAR(1)))**2) |
694 |
ZZ(1,NDEP) = (YV(I)-YVM)/SQRT(ABS((EXV(I))**2-EYM**2)) |
695 |
ENDDO |
696 |
ENDIF |
697 |
C |
698 |
C --- PER LE VISTE Y --- |
699 |
C METTE NEI VETTORI XV YV EYV I DATI DA FITTARE (ESCLUDE GLI ZERI) |
700 |
C |
701 |
CALL VZERO(XV,NPLA) |
702 |
CALL VZERO(YV,NPLA) |
703 |
CALL VZERO(EYV,NPLA) |
704 |
CALL VZERO(DELTY,NPLA) |
705 |
FLA = 0 |
706 |
DO KK = 1, NPLA !!finpar |
707 |
IF (Y2(KK).GT.0.) THEN |
708 |
FLA = FLA + 1 |
709 |
XV(FLA) = X2(KK) |
710 |
YV(FLA) = Y2(KK) |
711 |
EYV(FLA) = EY(KK) |
712 |
PY = X2(KK)*TAN(THY) + QQY + 12.1 |
713 |
DELTY(KK) = PY - YV(FLA) |
714 |
ENDIF |
715 |
ENDDO |
716 |
IF (FLA.LT.2) THEN |
717 |
c IF (BANDI.EQ.1.) SCARTA = SCARTA + 1. |
718 |
GO TO 449 |
719 |
ENDIF |
720 |
C |
721 |
C COME PER LE VISTE X |
722 |
C |
723 |
IF (F.EQ.1) THEN |
724 |
CALL BORDO(FLA) |
725 |
IF (FLA.LT.2) THEN |
726 |
c IF (BANDI.EQ.1.) SCARTA = SCARTA + 1. |
727 |
GO TO 449 |
728 |
ENDIF |
729 |
ENDIF |
730 |
C |
731 |
C ESEGUE IL FIT DELLA RETTA |
732 |
C |
733 |
CHOPT = 'QE' |
734 |
IF (F.EQ.1) THEN |
735 |
PA(1) = DBLE(12.1) |
736 |
PA(2) = DBLE(0.) |
737 |
ELSE |
738 |
PA(1) = QY |
739 |
PA(2) = MY |
740 |
ENDIF |
741 |
C |
742 |
CALL HFITV(FLA,NPLA,NVAR,XV,YV,EYV,FUN,CHOPT,NPAR,PA,STEP, |
743 |
& PMIN,PMAX,SIGPAR,CHI2) |
744 |
C |
745 |
C MY COEFF. ANGOLARE, QY QUOTA |
746 |
C |
747 |
MY = PA(2) |
748 |
QY = PA(1) |
749 |
C |
750 |
C REGISTRA I PULL |
751 |
C |
752 |
IF (F.EQ.2) THEN |
753 |
DO I = 1, FLA2 |
754 |
NDEP = NINT(1.+((XV(I)-.843)/.862)) |
755 |
YVM = XV(I)*SNGL(MY) + SNGL(QY) |
756 |
EYM = SQRT((XV(I)*SNGL(SIGPAR(2)))**2+(SNGL(SIGPAR(1)))**2) |
757 |
ZZ(2,NDEP) = (YV(I)-YVM)/SQRT(ABS((EYV(I))**2-EYM**2)) |
758 |
ENDDO |
759 |
ENDIF |
760 |
C |
761 |
C RIPETE UNA SECONDA VOLTA IL FIT |
762 |
C |
763 |
ENDDO |
764 |
C |
765 |
C CORREZIONE DOVUTA ALLA NON SIMMETRIA DELLO SCIAME |
766 |
C |
767 |
cc QX = QX - DBLE(.1) * DATAN(MX) !.038818 |
768 |
cc QY = QY - DBLE(.1) * DATAN(MY) !.243585 |
769 |
cc MX = DTAN(DATAN(MX)/DBLE(1.-.0064651)) !.00506855)) |
770 |
cc MY = DTAN(DATAN(MY)/DBLE(1.-.0064651)) !.00506855)) |
771 |
C |
772 |
C CALCOLA THETA (analoga e` DACOSD(1./SQRT(1.+(MX**2)+(MY**2))) ) |
773 |
C |
774 |
TMISD = DATAND(DSQRT((MX*MX) + (MY*MY))) |
775 |
C |
776 |
C CALCOLA PHI |
777 |
C |
778 |
IF (MX.EQ.0..AND.MY.GT.0.) PMISD = DBLE(90.) |
779 |
IF (MX.EQ.0..AND.MY.LT.0.) PMISD = DBLE(270.) |
780 |
IF (MY.EQ.0..AND.MX.GE.0.) PMISD = DBLE(0.) |
781 |
IF (MY.EQ.0..AND.MX.LT.0.) PMISD = DBLE(180.) |
782 |
C |
783 |
IF (MY.NE.0..AND.MX.NE.0.) THEN |
784 |
PMID = DATAND(MY/MX) |
785 |
IF (MY.LT.0..AND.MX.GT.0.) PMISD = PMID + DBLE(360.) |
786 |
IF (MX.LT.0.) PMISD = PMID + DBLE(180.) |
787 |
IF (MY.GT.0..AND.MX.GT.0.) PMISD = PMID |
788 |
ENDIF |
789 |
C |
790 |
C TROVA IL VERTICE DI INTERAZIONE SECONDO IL FIT |
791 |
C |
792 |
QXP = QX - DBLE(12.1) |
793 |
QYP = QY - DBLE(12.1) |
794 |
C |
795 |
ZETAM = 18.964 |
796 |
C |
797 |
IF (QXP.GT.12.1) THEN |
798 |
ZETAM = SNGL(DBLE(18.964)+(QXP-12.1)/MX) |
799 |
QXP = DBLE(12.1) |
800 |
ENDIF |
801 |
IF (QXP.LT.-12.1) THEN |
802 |
ZETAM = SNGL(DBLE(18.964)+(QXP+12.1)/MX) |
803 |
QXP = DBLE(-12.1) |
804 |
ENDIF |
805 |
IF (QYP.GT.12.1) THEN |
806 |
ZETAM = SNGL(DBLE(18.964)+(QYP-12.1)/MY) |
807 |
QYP = DBLE(12.1) |
808 |
ENDIF |
809 |
IF (QYP.LT.-12.1) THEN |
810 |
ZETAM = SNGL(DBLE(18.964)+(QYP+12.1)/MY) |
811 |
QYP = DBLE(-12.1) |
812 |
ENDIF |
813 |
C |
814 |
C |
815 |
C EVENTI CONTENUTI SECONDO IL FIT |
816 |
C |
817 |
POSIXMD = QXP + DBLE(ZETAM)*DTAND(TMISD)*DCOSD(PMISD) |
818 |
POSIYMD = QYP + DBLE(ZETAM)*DTAND(TMISD)*DSIND(PMISD) |
819 |
C |
820 |
IF (POSIXMD.LE.LIMLA.AND.POSIXMD.GE.-LIMLA.AND.POSIYMD.LE. |
821 |
& LIMLA.AND.POSIYMD.GE.-LIMLA) EVCONMIS = EVCONMIS + 1 |
822 |
C |
823 |
C EVENTI CONTENUTI SECONDO SIA SECONDO IL FIT CHE REALMENTE CONTENUTI |
824 |
C |
825 |
IF ((POSIXMD.LE.LIMLA.AND.POSIXMD.GE.-LIMLA.AND.POSIYMD.LE. |
826 |
& LIMLA.AND.POSIYMD.GE.-LIMLA).AND. |
827 |
& (POSIXD.LE.LIMLA.AND.POSIXD.GE.-LIMLA.AND.POSIYD.LE. |
828 |
& LIMLA.AND.POSIYD.GE.-LIMLA)) EVCOINC = EVCOINC + 1 |
829 |
C |
830 |
C CONVERTE IN RADIANTI E IN SINGOLA PRECISIONE PER L'ENTUPLA |
831 |
C |
832 |
C ANGOLI THETA E PHI |
833 |
C |
834 |
THET = SNGL(THETA*PI/DBLE(180.)) |
835 |
PH = SNGL(PHI*PI/DBLE(180.)) |
836 |
TMIS = SNGL(TMISD*PI/DBLE(180.)) |
837 |
PMIS = SNGL(PMISD*PI/DBLE(180.)) |
838 |
C |
839 |
C PROIEZIONE DEGLI ANGOLI SUI PIANI XZ E YZ |
840 |
C |
841 |
IF (PHI.LT.90..OR.PHI.GT.270.) THEN |
842 |
THX = SNGL(DATAN2(DABS(DTAND(THETA)),DSQRT( |
843 |
& DBLE(1.)+DTAND(PHI)**2.))) |
844 |
THY = SNGL(DATAN2(DTAND(PHI)*DABS(DTAND(THETA)), |
845 |
& DSQRT(DBLE(1.)+DTAND(PHI)**2.))) |
846 |
ENDIF |
847 |
IF (PHI.GT.90..AND.PHI.LT.270.) THEN |
848 |
THX = SNGL(DATAN2(-DABS(DTAND(THETA)),DSQRT( |
849 |
& DBLE(1.)+DTAND(PHI)**2.))) |
850 |
THY = SNGL(DATAN2(-DTAND(PHI)*DABS(DTAND(THETA)), |
851 |
& DSQRT(DBLE(1.)+DTAND(PHI)**2.))) |
852 |
ENDIF |
853 |
C |
854 |
C PROIEZIONE DEGLI ANGOLI MISURATI |
855 |
C |
856 |
THXM = SNGL(DATAN(MX)) |
857 |
THYM = SNGL(DATAN(MY)) |
858 |
C |
859 |
C LA PROIEZIONE SULL'ULTIMO PIANO |
860 |
C |
861 |
POSIXM = SNGL(POSIXMD) |
862 |
POSIYM = SNGL(POSIYMD) |
863 |
POSIX = SNGL(POSIXD) |
864 |
POSIY = SNGL(POSIYD) |
865 |
C |
866 |
C LE QUOTE REALI E MISURATE |
867 |
C |
868 |
QQX = SNGL(VALORI(LENGHT-13)) |
869 |
QQY = SNGL(VALORI(LENGHT-12)) |
870 |
QQZ = ZETA |
871 |
QQXM = SNGL(QXP) |
872 |
QQYM = SNGL(QYP) |
873 |
QQZM = ZETAM |
874 |
C |
875 |
C REGISTRA IN ENMEZX(Y) LE ENERGIE RILASCIATE SU TUTTE LE STRIP |
876 |
C A DX E A SX DELLA VERA TRAIETTORIA DELLA PARTICELLA PER I PIANI X(Y) |
877 |
C |
878 |
CALL VZERO(ENMEZX,2*22) |
879 |
CALL VZERO(ENMEZY,2*22) |
880 |
CALL VZERO(QPLX,6*22) |
881 |
CALL VZERO(QPLY,6*22) |
882 |
C |
883 |
DO LL = 1, NPLA |
884 |
PMETX = 0. |
885 |
PMETY = 0. |
886 |
PY = X(LL)*TAN(THX) + QQX + 12.1 |
887 |
CALL TRAIE(PY,PMETX) |
888 |
PY = X2(LL)*TAN(THY) + QQY + 12.1 |
889 |
CALL TRAIE(PY,PMETY) |
890 |
DO JG = 1, NCHA |
891 |
NLK = (JG - 1)*NPLA + LL |
892 |
IF (PMETX.LT.0.) GO TO 446 |
893 |
IF (JG.LT.NINT(PMETX)) |
894 |
& ENMEZX(1,LL) = ENMEZX(1,LL) + ZXY(1,NLK) |
895 |
IF (JG.EQ.NINT(PMETX)) THEN |
896 |
ENMEZX(1,LL) = ENMEZX(1,LL) + ZXY(1,NLK)/2. |
897 |
ENMEZX(2,LL) = ENMEZX(2,LL) + ZXY(1,NLK)/2. |
898 |
ENDIF |
899 |
IF (JG.GT.NINT(PMETX)) |
900 |
& ENMEZX(2,LL) = ENMEZX(2,LL) + ZXY(1,NLK) |
901 |
446 CONTINUE |
902 |
IF (PMETY.LT.0.) GO TO 447 |
903 |
IF (JG.LT.NINT(PMETY)) |
904 |
& ENMEZY(1,LL) = ENMEZY(1,LL) + ZXY(2,NLK) |
905 |
IF (JG.EQ.NINT(PMETY)) THEN |
906 |
ENMEZY(1,LL) = ENMEZY(1,LL) + ZXY(2,NLK)/2. |
907 |
ENMEZY(2,LL) = ENMEZY(2,LL) + ZXY(2,NLK)/2. |
908 |
ENDIF |
909 |
IF (JG.GT.NINT(PMETY)) |
910 |
& ENMEZY(2,LL) = ENMEZY(2,LL) + ZXY(2,NLK) |
911 |
447 CONTINUE |
912 |
ENDDO |
913 |
C |
914 |
DO RT = 1, 6 |
915 |
DO RE = 1, 16 |
916 |
NLL = (RE+(RT-1)*16 - 1)*NPLA + LL |
917 |
QPLX(RT,LL) = ZXY(1,NLL) + QPLX(RT,LL) |
918 |
QPLY(RT,LL) = ZXY(2,NLL) + QPLY(RT,LL) |
919 |
ENDDO |
920 |
ENDDO |
921 |
C |
922 |
ENDDO |
923 |
C |
924 |
C MISURA DELL'ENERGIA: SI INTEGRA L'ENERGIA RIVELATA FINO AL PIANO DI |
925 |
C MASSIMO |
926 |
C |
927 |
PARZEN = 0. |
928 |
DO LL = 1, FINPAR |
929 |
PARZEN = PARZEN + PARZ(1,LL) + PARZ(2,LL) |
930 |
ENDDO |
931 |
FFLA = FLOAT(FINPAR) |
932 |
IF (FINPAR1.EQ.1) THEN |
933 |
PARZEN = PARZEN - PARZ(2,FINPAR) |
934 |
FFLA = FFLA - 1. |
935 |
ENDIF |
936 |
C |
937 |
C CORREZIONE DELL'ENERGIA PER EVENTI USCITI DAL CALORIMETRO |
938 |
C |
939 |
cc IF ((FFLA/COS(TMIS)).LT.28..AND.(POSIXMD.LT.-12.1.OR.POSIXMD |
940 |
cc & .GT.12.1.OR.POSIYMD.LT.-12.1.OR.POSIYMD.GT.12.1)) THEN |
941 |
cc PARZEN = PARZEN + PARZEN*ENCORR(FFLA/COS(TMIS)) |
942 |
cc ENDIF |
943 |
C |
944 |
C |
945 |
ERZEN = ENERGER(PARZEN) |
946 |
ER3 = ER3 + PARZEN*(ERZEN**(-2)) |
947 |
ER4 = ER4 + ERZEN**(-2) |
948 |
C |
949 |
ERENET = ENERGER(ENET) |
950 |
VENET = VALORI(LENGHT-14) |
951 |
C |
952 |
ER1 = ER1 + ENET*(ERENET**(-2)) |
953 |
ER2 = ER2 + ERENET**(-2) |
954 |
C |
955 |
ENETM = ER1/ER2 |
956 |
C |
957 |
ERENETM = ER2**(-.5) |
958 |
C |
959 |
ENM = ENM + ENET |
960 |
ERENM = ERENM + ERENET |
961 |
CCC |
962 |
CALL VZERO(XEN,NPL) |
963 |
CALL VZERO(YEN,NPL) |
964 |
CALL VZERO(EYEN,NPL) |
965 |
CALL VZERO(XEN2,NPL) |
966 |
CALL VZERO(YEN2,NPL) |
967 |
CALL VZERO(EYEN2,NPL) |
968 |
IFIN2 = FINPAR !NINT(FFLA2) |
969 |
DO NN = 1, NPLA |
970 |
c NPIA = NN - 3 + IFIN2 |
971 |
NPIA = NN |
972 |
IF (NPIA.GT.22) NPIA = 22 |
973 |
IF (NPIA.LT.0) NPIA = 0 |
974 |
RN = FLOAT(NPIA) |
975 |
c |
976 |
c XEN(NN) = (RN - 1.)*0.75/COS(TMIS) ! (RN - 1.)*.862 + .169 |
977 |
c XEN(NN+22) = RN*0.75/COS(TMIS) !(RN - 1.)*.862 + .843 |
978 |
c |
979 |
YEN(NN) = PARZ(1,NPIA) |
980 |
YEN(NN+22) = PARZ(2,NPIA) |
981 |
EYEN(NN) = 5.*SQRT(YEN(NN)) |
982 |
EYEN(NN+22) = 5.*SQRT(YEN(NN+22)) |
983 |
ENDDO |
984 |
NPLL = 0 |
985 |
DO NN = 1, NPLA |
986 |
IF (YEN(NN).NE.0..AND.NN.EQ.1) THEN |
987 |
NPLL = NPLL + 1 |
988 |
XEN2(NPLL) = NPLL*.76/COS(TMIS) !XEN(NN) |
989 |
YEN2(NPLL) = YEN(NN) |
990 |
EYEN2(NPLL) = EYEN(NN) |
991 |
ENDIF |
992 |
IF (YEN(NN+22).NE.0..AND.NN.EQ.22) THEN |
993 |
NPLL = NPLL + 1 |
994 |
XEN2(NPLL) = NPLL*.76/COS(TMIS) !XEN(NN+22) |
995 |
YEN2(NPLL) = YEN(NN+22) |
996 |
EYEN2(NPLL) = EYEN(NN+22) |
997 |
ENDIF |
998 |
IF (NN.GT.1.AND.NN.LT.22.AND.YEN(NN+21). |
999 |
& NE.0..AND.YEN(NN).NE.0.) THEN |
1000 |
NPLL = NPLL + 1 |
1001 |
XEN2(NPLL) = NPLL*.76/COS(TMIS) !XEN(NN) |
1002 |
YEN2(NPLL) = (YEN(NN) + YEN(NN+21))/2. |
1003 |
EYEN2(NPLL) = (EYEN(NN) + EYEN(NN+21))/2. |
1004 |
ENDIF |
1005 |
IF (NN.GT.1.AND.NN.LT.22.AND.YEN(NN+21). |
1006 |
& NE.0..AND.YEN(NN).EQ.0.) THEN |
1007 |
NPLL = NPLL + 1 |
1008 |
XEN2(NPLL) = NPLL*.76/COS(TMIS) !XEN(NN) |
1009 |
YEN2(NPLL) = YEN(NN+21) |
1010 |
EYEN2(NPLL) = EYEN(NN+21) |
1011 |
ENDIF |
1012 |
IF (NN.GT.1.AND.NN.LT.22.AND.YEN(NN+21). |
1013 |
& EQ.0..AND.YEN(NN).NE.0.) THEN |
1014 |
NPLL = NPLL + 1 |
1015 |
XEN2(NPLL) = NPLL*.76/COS(TMIS) !XEN(NN) |
1016 |
YEN2(NPLL) = YEN(NN) |
1017 |
EYEN2(NPLL) = EYEN(NN) |
1018 |
ENDIF |
1019 |
ENDDO |
1020 |
IF (NPLL.LT.13..OR.BETA.NE.'Y') THEN |
1021 |
PAR(1) = -1. |
1022 |
PAR(2) = -1. |
1023 |
POSMAX = -1. |
1024 |
PAR1 = -1. |
1025 |
GO TO 448 |
1026 |
ENDIF |
1027 |
CHOPT = 'BEQ' |
1028 |
PAR(1) = 1. |
1029 |
PAR(2) = 1. |
1030 |
PMI2(1) = 2. |
1031 |
PMI2(2) = -0.5 !.4 |
1032 |
PMA(1) = 20. !11. |
1033 |
PMA(2) = 6. |
1034 |
STE(1) = 0.01 |
1035 |
STE(2) = 0.01 |
1036 |
C |
1037 |
CALL HFITV(NPLL,NPL,NVAR,XEN2,YEN2,EYEN2,ENERGIF,CHOPT, |
1038 |
& NPAR2,PAR,STE,PMI2,PMA,SIGPA,CHI) |
1039 |
CCC |
1040 |
IF (PAR(2).NE.0.) THEN |
1041 |
POSMAX = PAR(1)/PAR(2) |
1042 |
PAR1 = PAR(1) |
1043 |
ELSE |
1044 |
PAR1 = PAR(1) |
1045 |
POSMAX = 0. |
1046 |
ENDIF |
1047 |
448 CONTINUE |
1048 |
C |
1049 |
C RIEMPIE LE ENTUPLE |
1050 |
C |
1051 |
CALL HFNT(1) |
1052 |
C |
1053 |
TRIG = TRIG + 1 |
1054 |
C |
1055 |
c IF ((THET-TMIS).LT.-.3) THEN |
1056 |
c PRINT *,'FLA2 ',FLA2 |
1057 |
c PRINT *,'FLA ',FLA |
1058 |
c PRINT *,'QXP ',DBLE(VALORI(LENGHT-13)) |
1059 |
c PRINT *,'QXPM ',QXP |
1060 |
c PRINT *,'QYP ',DBLE(VALORI(LENGHT-12)) |
1061 |
c PRINT *,'QYPM ',QYP |
1062 |
c PRINT *,'THETA: ',THETA |
1063 |
c PRINT *,'TMIS ',DATAND(DSQRT((MX**2) + (MY**2))) |
1064 |
C PRINT *,'thet: ',thet |
1065 |
c PRINT *,'TMIS: ',TMISD |
1066 |
c PRINT *,'PHI: ',PHI |
1067 |
c PRINT *,'PMIS ',PMISD |
1068 |
c PRINT *,'TRIG ',TRIG |
1069 |
c PRINT *,'***********' |
1070 |
c ENDIF |
1071 |
C |
1072 |
449 CONTINUE |
1073 |
C |
1074 |
C ESEGUE PER TUTTI GLI EVENTI |
1075 |
C |
1076 |
450 ENDDO |
1077 |
C CALL HFNT(2) |
1078 |
500 CONTINUE |
1079 |
900 CONTINUE |
1080 |
C |
1081 |
CLOSE(83) |
1082 |
C |
1083 |
C PRINT RESULTS |
1084 |
C |
1085 |
PRINT *,'FILE PROCESSATO: ',FILENAME |
1086 |
PRINT *,'FILE DATI SALVATO: ',NAME |
1087 |
PRINT *,'EVENTI TOTALI: ',EVENTO |
1088 |
PRINT *,'EVENTI CONTENUTI REALEMENTE: ',EVCON |
1089 |
PRINT *,'EVENTI ANALIZZATI: ',TRIG |
1090 |
PRINT *,'EVENTI CONTENUTI MISURATI: ',EVCONMIS |
1091 |
PRINT *,'EVENTI CONTENUTI COINCIDENTI: ',EVCOINC |
1092 |
PRINT *,'EVENTI BUONI: ',TOTTO2 |
1093 |
PRINT *,'EVENTI SCARTATI: ',SCARTA |
1094 |
c PRINT *,'ENERGIA MISURATA PESATA: ',ENETM |
1095 |
c PRINT *,'ERRORE EN. PESATA: ',ERENETM |
1096 |
C PRINT *,'ENERGIA REALE ',VENET/CALIB |
1097 |
C PRINT *,'ENERGIA MEDIA MISURATA: ',ENM/TRIG |
1098 |
C PRINT *,'ERRORE MEDIO MISURATO: ',ERENM/TRIG |
1099 |
c PRINT *,'EN. MIS. PARZIALE PESATA: ',ER3/ER4 |
1100 |
c PRINT *,'ERRORE EN. PARZ. PESATA: ',ER4**(-.5) |
1101 |
C |
1102 |
950 CONTINUE |
1103 |
C |
1104 |
CALL HROUT(0,ICYCLE,' ') |
1105 |
CALL HREND('SIMULATION') |
1106 |
1000 CONTINUE |
1107 |
C |
1108 |
5050 FORMAT(A40) |
1109 |
5060 FORMAT('Y',I2) |
1110 |
5070 FORMAT('NUEV') |
1111 |
5080 FORMAT(1X,I5) |
1112 |
6001 FORMAT(A1) |
1113 |
6005 FORMAT(F3.3) |
1114 |
6006 FORMAT(F4.3) |
1115 |
C |
1116 |
END |
1117 |
C |
1118 |
C--------------------------------------------------------------------- |
1119 |
SUBROUTINE LEGGE(VALORI,LENGHT) |
1120 |
C--------------------------------------------------------------------- |
1121 |
C |
1122 |
REAL VALORI(LENGHT) |
1123 |
C |
1124 |
READ(83),VALORI |
1125 |
C |
1126 |
RETURN |
1127 |
END |
1128 |
C |
1129 |
C--------------------------------------------------------------------- |
1130 |
SUBROUTINE BORDO(NB) |
1131 |
C--------------------------------------------------------------------- |
1132 |
C |
1133 |
PARAMETER (STMA=23.616,STMI=.584) |
1134 |
PARAMETER (NPLA=22) |
1135 |
INTEGER CONF, CONF2, NB, STORE |
1136 |
REAL YV |
1137 |
C |
1138 |
COMMON/VETFIT/YV(NPLA) |
1139 |
C |
1140 |
C ROUTINE PER SCARTARE DAL FIT I PUNTI CHE SCORRONO SUL BORDO DEL |
1141 |
C CALORIMETRO E CHE NON SONO IL PROSEGUIMENTO DELLA VERA TRAIETTORIA |
1142 |
C |
1143 |
CONF = 0 |
1144 |
CONF2 = 0 |
1145 |
DO KK = 1, NB |
1146 |
IF (YV(KK).LT.STMA.OR.YV(KK).GT.STMI) CONF = 1 |
1147 |
IF ((YV(KK).GE.STMA.OR.YV(KK).LE.STMI).AND.CONF.EQ.1 |
1148 |
& .AND.CONF2.EQ.0) STORE = KK |
1149 |
IF (YV(KK).GE.STMA.OR.YV(KK).LE.STMI.AND.CONF.EQ.1) |
1150 |
& CONF2 = CONF2 + 1 |
1151 |
IF (CONF2.GE.2) THEN |
1152 |
NB = STORE |
1153 |
GO TO 7000 |
1154 |
ENDIF |
1155 |
ENDDO |
1156 |
C |
1157 |
7000 CONTINUE |
1158 |
RETURN |
1159 |
END |
1160 |
C |
1161 |
C--------------------------------------------------------------------- |
1162 |
SUBROUTINE TRAIE(YH,TRH) |
1163 |
C--------------------------------------------------------------------- |
1164 |
C |
1165 |
REAL YH, TRH |
1166 |
C |
1167 |
IF (YH.GE.0..AND.YH.LE.24.2) THEN |
1168 |
IF (YH.LE..218) TRH = 1. |
1169 |
IF (YH.GT..218.AND.YH.LT.7.782) |
1170 |
& TRH = 1. + ( YH - .218) / .244 |
1171 |
IF (YH.LT.8.05.AND.YH.GE.7.782) |
1172 |
& TRH = 32. |
1173 |
IF (YH.LE.8.318.AND.YH.GT.8.05) |
1174 |
& TRH = 33. |
1175 |
IF (YH.GT.8.318.AND.YH.LT.15.882) |
1176 |
& TRH = 1. + (YH - .51) / .244 |
1177 |
IF (YH.LT.16.15.AND.YH.GE.15.882) |
1178 |
& TRH = 64. |
1179 |
IF (YH.LE.16.418.AND.YH.GT.16.15) |
1180 |
& TRH = 65. |
1181 |
IF (YH.GT.16.418.AND.YH.LT.23.982) |
1182 |
& TRH = 1. + (YH - .802) / .244 |
1183 |
IF (YH.GE.23.982) TRH = 96. |
1184 |
ELSE |
1185 |
TRH = -10. |
1186 |
ENDIF |
1187 |
RETURN |
1188 |
END |
1189 |
C |
1190 |
C--------------------------------------------------------------------- |
1191 |
SUBROUTINE BARIC(TRA,YY,EYY,NNF,FF,NX) |
1192 |
C--------------------------------------------------------------------- |
1193 |
C |
1194 |
PARAMETER (STRIP1=7.,STRIP2=5.) !!7 5 |
1195 |
PARAMETER (NPLA=22, NCHA=96) |
1196 |
PARAMETER (LENSEV=NPLA*NCHA) |
1197 |
C |
1198 |
REAL TRA, YY, EYY |
1199 |
INTEGER ST0, ST1, ST2, ST3 |
1200 |
REAL UFF, ENTO |
1201 |
REAL ZXY, YN |
1202 |
C |
1203 |
INTEGER NLK, NNF, FF, NX, NPOS |
1204 |
C |
1205 |
COMMON/ENSTRIP/ZXY(2,LENSEV),YN(NCHA) |
1206 |
C |
1207 |
ST0 = NINT(STRIP1) |
1208 |
ST1 = NINT(STRIP1-STRIP2) |
1209 |
ST2 = NINT((STRIP1-1.)/2. + 1.) |
1210 |
ST3 = NINT(FLOAT(ST1)/2.) |
1211 |
C |
1212 |
IF (TRA.GT.0.) THEN |
1213 |
UFF = 0. |
1214 |
ENTO = 0. |
1215 |
DO P = 1, (ST0-(FF-1)*ST1) |
1216 |
IF (NINT(TRA).LT.(ST2-ST3*(FF-1))) THEN |
1217 |
IF (P.GT.(ST0-ABS(NINT(TRA)-ST2-ST3*(FF-1)))) GO TO 7100 |
1218 |
NPOS = P |
1219 |
ELSE |
1220 |
IF ((NINT(TRA)+P-ST2+ST3*(FF-1)).GT.96) GO TO 7100 |
1221 |
NPOS = NINT(TRA)+P-ST2+ST3*(FF-1) |
1222 |
ENDIF |
1223 |
NLK = (NPOS-1)*NPLA + NNF |
1224 |
UFF = UFF + YN(NPOS)*ZXY(NX,NLK) |
1225 |
ENTO = ENTO + ZXY(NX,NLK) |
1226 |
ENDDO |
1227 |
7100 CONTINUE |
1228 |
IF (ENTO.GT.0.) THEN |
1229 |
YY = UFF/ENTO |
1230 |
c EYY = 1. |
1231 |
EYY = YY/(ENTO**.79) |
1232 |
c IF (NNF.GT.8.AND.NNF.LT.20) THEN |
1233 |
cc EYY = YY/(ENTO**(.85 + 2.*EXP(-NNF*1.5))) |
1234 |
c EYY = YY/(ENTO**(.71)) |
1235 |
c ELSE |
1236 |
c EYY = YY/(ENTO**(.89)) |
1237 |
c ENDIF |
1238 |
ELSE |
1239 |
YY = -1. |
1240 |
EYY = 1E5 |
1241 |
ENDIF |
1242 |
IF (NNF.LT.2.AND.YY.GT.0.) EYY = .244/SQRT(12.) |
1243 |
ccc EYY = .244/SQRT(12.) |
1244 |
ccc |
1245 |
ELSE |
1246 |
YY = -1. |
1247 |
EYY = 1E5 |
1248 |
ENDIF |
1249 |
RETURN |
1250 |
END |
1251 |
C |
1252 |
C--------------------------------------------------------------------------- |
1253 |
REAL FUNCTION FUN(X) |
1254 |
C--------------------------------------------------------------------------- |
1255 |
C |
1256 |
PARAMETER (NPAR=2, NPAR2=2) |
1257 |
PARAMETER (NWPAWC=3000000) |
1258 |
COMMON/HCFITD/PA(NPAR),SIGPAR(NPAR),CHI2,PMIN(NPAR),PMAX(NPAR), |
1259 |
& STEP(NPAR), |
1260 |
& PAR(NPAR2),SIGPA(NPAR2),CHI,PMI2(NPAR2),PMA(NPAR2), |
1261 |
& STE(NPAR2) |
1262 |
COMMON/PAWC/HMEMOR(NWPAWC) |
1263 |
C |
1264 |
FUN = PA(1) + X*PA(2) |
1265 |
C |
1266 |
END |
1267 |
C |
1268 |
C--------------------------------------------------------------------------- |
1269 |
REAL FUNCTION ENERGER(X) |
1270 |
C--------------------------------------------------------------------------- |
1271 |
C |
1272 |
PARAMETER (CALIB=0.0001059994) |
1273 |
C |
1274 |
ENERGER = .0205 * SQRT(7.4 * 6.76/(X * CALIB)) * X |
1275 |
C |
1276 |
END |
1277 |
C |
1278 |
C--------------------------------------------------------------------------- |
1279 |
REAL FUNCTION ENCORR(X) |
1280 |
C--------------------------------------------------------------------------- |
1281 |
C |
1282 |
REAL FFUN |
1283 |
C |
1284 |
FFUN = .5+.353*ATAN((X-21.972)*.374) |
1285 |
c FFUN = .5466+.3831*ATAN((X-21.386)*.343) |
1286 |
IF (X.LT.28.) THEN |
1287 |
IF (FFUN.LT.0.) FFUN = 1.E-5 |
1288 |
ENCORR = .45/FFUN |
1289 |
ELSE |
1290 |
ENCORR = 0. |
1291 |
ENDIF |
1292 |
C |
1293 |
END |
1294 |
C |
1295 |
C--------------------------------------------------------------------------- |
1296 |
REAL FUNCTION ENERGIF(X) |
1297 |
C--------------------------------------------------------------------------- |
1298 |
PARAMETER (NPAR=2, NPAR2=2) |
1299 |
PARAMETER (NWPAWC=3000000) |
1300 |
COMMON/HCFITD/PA(NPAR),SIGPAR(NPAR),CHI2,PMIN(NPAR),PMAX(NPAR), |
1301 |
& STEP(NPAR), |
1302 |
& PAR(NPAR2),SIGPA(NPAR2),CHI,PMI2(NPAR2),PMA(NPAR2), |
1303 |
& STE(NPAR2) |
1304 |
COMMON/PAWC/HMEMOR(NWPAWC) |
1305 |
C |
1306 |
C Rossi equation B |
1307 |
C |
1308 |
ENERGIF = (0.8*(X**(PAR(1)))+0.2)*EXP(-X*PAR(2)) |
1309 |
C |
1310 |
END |
1311 |
|