/[PAMELA software]/DarthVader/CalorimeterLevel2/src/calol2tr.for
ViewVC logotype

Annotation of /DarthVader/CalorimeterLevel2/src/calol2tr.for

Parent Directory Parent Directory | Revision Log Revision Log

Revision 1.16 - (hide annotations) (download)
Fri Jul 20 08:24:53 2007 UTC (17 years, 4 months ago) by mocchiut
Branch: MAIN
CVS Tags: v4r00
Changes since 1.15: +54 -29 lines 
Formal changes to use fortran routines in the presampler analysis
1 mocchiut 1.1 *****************************************************************************
2     INTEGER FUNCTION CALOL2TR()
3     c
4     IMPLICIT NONE
5     C
6     INCLUDE 'INTEST.TXT'
7     C
8     DOUBLE PRECISION al_p(5),
9 mocchiut 1.16 & xout(nplav),yout(nplav),zin(nplav)
10 mocchiut 1.1 C
11 mocchiut 1.16 REAL PIANO(NPLAV), VARFIT(2)
12 mocchiut 1.1 REAL TX, TY, SHIFT
13 mocchiut 1.16 REAL BAR(2,NPLAV), DISTY
14     REAL DISTX, Y(NPLAV), YY(NPLAV)
15 mocchiut 1.1 REAL RIG, PPLANEMAX, RMASS
16     REAL RNSS, QTOTT, RQT, MX, MY
17     REAL CHECK, ENER, CX, CY
18     REAL EINF, ESUP, RPIANO(2)
19     REAL hmemor(9000000), X01PL
20     C
21     REAL ax,bx,eax,ebx,chi2x
22     REAL ay,by,eay,eby,chi2y
23     REAL parzen3, TMISD
24     INTEGER Nfitx,Nfity
25     C
26     INTEGER INDEX, NTOT(2), NPIANI, GTR
27     INTEGER j, m, i, IWPL(2), timpx, timpy, T, nn
28     INTEGER IPLANE, NNX, NNY, INFX, INFY, ISUPX, ISUPY
29 mocchiut 1.16 INTEGER IBAR(2,NPLAV), NPFIT(2), CHTRACK,IWPLU
30 mocchiut 1.1 INTEGER Iquest(100), ICONTROL5, nin, IFAIL
31     C
32     PARAMETER (X01PL=0.74)
33     C
34    
35     C
36     COMMON / slftrig / tmisd,ax,bx,eax,ebx,chi2x,Nfitx,ay,by,eay,eby,
37     & chi2y,Nfity,parzen3
38     SAVE / slftrig /
39     C
40     COMMON / TAGLIOEN / EINF, ESUP, ENER(2)
41     SAVE / TAGLIOEN /
42     C
43     COMMON / SHIFT / SHIFT
44     SAVE / SHIFT /
45     C
46     COMMON / ANGOLO / BAR, IBAR
47     SAVE / ANGOLO /
48     C
49     COMMON / WHERE / CX, CY, PIANO
50     SAVE / WHERE /
51     C
52     COMMON / GENERAL / RIG, RMASS
53     SAVE / GENERAL /
54     C
55     COMMON / CH / CHECK
56     SAVE / CH /
57     C
58 mocchiut 1.8 COMMON / CALOFIT / VARFIT, NPFIT, IWPL,CHTRACK
59 mocchiut 1.1 SAVE / CALOFIT /
60     C
61     COMMON / pawcd / hmemor
62     save / pawcd /
63     C
64     Common / QUESTd / Iquest
65     save / questd /
66     C
67     C Begin !
68     C
69 mocchiut 1.16 c print *,' sono qui'
70 mocchiut 1.1 CALOL2TR = 0;
71     NCORE = 0.
72     QCORE = 0.
73     NOINT = 0.
74     QCYL = 0.
75     NCYL = 0.
76     QLOW = 0.
77     NLOW = 0.
78     QTR = 0.
79     NTR = 0.
80     QLAST = 0.
81     QTRACK = 0.
82     QPRESH = 0.
83     NPRESH = 0.
84     QTRACKX = 0.
85     QTRACKY = 0.
86     DXTRACK = 0.
87     DYTRACK = 0.
88     QPRE = 0.
89     NPRE = 0.
90     NLAST = 0.
91     PLANETOT = 0.
92     QMEAN = 0.
93 mocchiut 1.8 C SELFTRIGGER = 0
94 mocchiut 1.1 C
95     C BEGIN WITH THE FISRT TRACK IF WE HAVE A TRACK FROM TRACKER
96     C
97     T = 1
98     C
99     10 CONTINUE
100     C
101     IF (GOOD2.EQ.1) THEN
102     C
103     CHTRACK = 0
104     C
105     CALL VZERO(IWPL,2)
106 mocchiut 1.16 CALL VZERO(BAR,2*NPLAV)
107     CALL VZERO(IBAR,2*NPLAV)
108     CALL VZERO(TBAR,2*NPLAV)
109     CALL VZERO(TIBAR,2*NPLAV)
110     CALL VZERO(Y,NPLAV)
111     CALL VZERO(YY,NPLAV)
112     CALL VZERO(XOUT,NPLAV)
113     CALL VZERO(YOUT,NPLAV)
114 mocchiut 1.1 do m = 1, 5
115     al_p(m) = al_pp(t,m)
116 mocchiut 1.16 c print *,' al_p(',m,') = ',al_p(m)
117 mocchiut 1.1 enddo
118     if (al_p(5).eq.0.) THEN
119 mocchiut 1.4 PRINT *,' CALORIMETER - WARNING F77: track with R = 0, discarded'
120 mocchiut 1.1 GOOD2 = 0
121     GOTO 969
122     ENDIF
123     DO M = 1,2
124     DO I = 1,NPLA
125     XOUT(I) = 0.
126     YOUT(I) = 0.
127     IF (MOD(M,2).EQ.0) THEN
128     DISTX = PIANO(I) + ZALIG
129     ELSE
130     DISTX = PIANO(I) - 5.81 + ZALIG
131     ENDIF
132     ZIN(I) = distx / 10.
133 mocchiut 1.16 c print *,' ZIN(',I,') = ',ZIN(I)
134 mocchiut 1.1 TBAR(M,I) = 0.
135     TIBAR(M,I) = 0
136     enddo
137     IFAIL = 0
138     call DOTRACK(NPLA,ZIN,XOUT,YOUT,AL_P,IFAIL)
139     if(IFAIL.ne.0)then
140     GOOD2 = 0
141 mocchiut 1.3 c print *,' CALORIMETER - WARNING F77: tracking failed '
142 mocchiut 1.1 goto 969
143     endif
144     TX = TAN(ASIN(AL_P(3))) * COS(AL_P(4))
145     TY = TAN(ASIN(AL_P(3))) * SIN(AL_P(4))
146     DO I = 1, NPLA
147     NN = 0
148 mocchiut 1.14 IF (M.EQ.2) NN = 1
149 mocchiut 1.1 IF (MOD(I,2).EQ.NN) THEN
150 mocchiut 1.16 IF (REVERSE.EQ.0) THEN
151     SHIFT = +0.5
152     ELSE
153     SHIFT = -0.5
154     ENDIF
155 mocchiut 1.1 ELSE
156 mocchiut 1.16 IF (REVERSE.EQ.0) THEN
157     SHIFT = -0.5
158     ELSE
159     SHIFT = +0.5
160     ENDIF
161 mocchiut 1.1 ENDIF
162     C
163     C CHECK IF XOUT OR YOUT ARE NaN
164     C
165     IF (XOUT(I).NE.XOUT(I).OR.YOUT(I).NE.YOUT(I)) THEN
166 mocchiut 1.3 c print *,
167     c & ' CALORIMETER - WARNING F77: tracking error (NaN values)'
168 mocchiut 1.1 GOOD2 = 0
169     GOTO 969
170     ENDIF
171     C
172     CX = XOUT(I)*10. + XALIG
173 mocchiut 1.5 CY = YOUT(I)*10. + YALIG
174 mocchiut 1.1 C
175     IF (I.EQ.1) THEN
176     TIMPX = CX
177     TIMPY = CY
178     ENDIF
179     IF (M.EQ.1) THEN
180     Y(I) = CX
181     BAR(M,I) = Y(I)
182     TBAR(M,I) = (Y(I) - XALIG)/10.
183 mocchiut 1.16 IF (I.EQ.NPLA) MX=ABS(Y(1)-Y(NPLA))/
184     & ABS(ZIN(1)-ZIN(NPLA))
185 mocchiut 1.1 ELSE
186     YY(I) = CY
187     BAR(M,I) = YY(I)
188 mocchiut 1.5 TBAR(M,I) = (-YALIG + YY(I))/10.
189 mocchiut 1.16 IF (I.EQ.NPLA) MY=ABS(Y(1)-Y(NPLA))/
190     & ABS(ZIN(1)-ZIN(NPLA))
191 mocchiut 1.1 ENDIF
192     CALL LASTRISCIA(BAR(M,I),IBAR(M,I))
193     tibar(M,I) = ibar(m,i)
194     IF (ibar(m,i).EQ.-1) THEN
195     CHTRACK = CHTRACK + 1
196     ELSE
197     IWPL(M) = IWPL(M) + 1
198     ENDIF
199     ENDDO
200     ENDDO
201     969 continue
202     cC
203     cC IF WE HAVE A GOOD CALORIMETER FIT DOES IT MATCH WITH TRACKER FIT?
204     cC
205     c IF (GOOD2.EQ.1.AND.NPFIT(2).GT.15.AND.VARFIT(2).LT.1000
206     c & .AND.TRKCHI2.EQ.1) THEN
207     c IF (ABS(TBAR(2,1)-CBAR(2,1))<40.) THEN
208     cC
209     cC GOOD, THE TWO TRACKS COINCIDE
210     cC
211     c IF (T.EQ.2) TRKCHI2 = 2
212     c GOTO 6996
213     c ELSE
214     cC
215     cC IT IS NOT A GOOD FIT BUT WE HAVE AN IMAGE AND IT IS THE FIRST TRACK
216     cC
217     c IF (T.EQ.1) THEN
218     c T = 2
219     c GOTO 10
220     c ENDIF
221     c IF (T.EQ.2) THEN
222     c TRKCHI2 = -1
223     c T = 1
224     c GOTO 10
225     c ENDIF
226     c ENDIF
227     c ENDIF
228     C
229     IF (GOOD2.EQ.0) THEN
230     c IF (T.EQ.1.AND.TRKCHI2.EQ.1) THEN
231     c GOOD2 = 1
232     c T = 2
233     c GOTO 10
234     c ENDIF
235     GOTO 50
236     ENDIF
237     C
238     GOTO 6996
239     C
240     ENDIF
241     C
242     C WE MUST PROCESS A SELFTRIGGER EVENT
243     C
244 mocchiut 1.7 IF (TRIGTY.GE.2.AND.HZN.EQ.0) THEN
245 mocchiut 1.1 C
246     C CALL SELFTRIGGER SUBROUTINE
247     C
248 mocchiut 1.8 CALL VZERO(IWPL,2)
249 mocchiut 1.7 CALL VZERO(VARCFIT,2)
250     CALL VZERO(NPCFIT,2)
251 mocchiut 1.16 CALL VZERO(TBAR,2*NPLAV)
252     CALL VZERO(TIBAR,2*NPLAV)
253     CALL VZERO(BAR,2*NPLAV)
254     CALL VZERO(IBAR,2*NPLAV)
255     CALL VZERO(Y,NPLAV)
256     CALL VZERO(YY,NPLAV)
257     CALL VZERO(XOUT,NPLAV)
258     CALL VZERO(YOUT,NPLAV)
259 mocchiut 1.7 C
260 mocchiut 1.1 CALL SELFTRIG()
261     ELEN = PARZEN3
262     SELEN = ABS(ELEN) * (11.98*1E-2 + 7.6 * EXP(-5736/ABS(ELEN)))
263     C
264     NPCFIT(1) = NFITX
265     NPCFIT(2) = NFITY
266     C
267     DO M = 1,2
268     C
269 mocchiut 1.13 c print *,' ax ',ax,' ay ',ay
270     c print *,' bx ',bx,' by ',by
271 mocchiut 1.1 IF (NPCFIT(M).GE.2) THEN
272     IF (M.EQ.1) THEN
273     VARCFIT(1) = CHI2X
274 mocchiut 1.13 IMPX = AX ! PAMELA REF
275 mocchiut 1.1 TANX = BX
276     ELSE
277     VARCFIT(2) = CHI2Y
278 mocchiut 1.13 IMPY = AY ! PAMELA REF
279 mocchiut 1.1 TANY = BY
280     ENDIF
281     C
282     DO I = 1,NPLA
283     NN = 0
284 mocchiut 1.14 IF (M.EQ.2) NN = 1
285 mocchiut 1.1 IF (MOD(I,2).EQ.NN) THEN
286 mocchiut 1.16 IF (REVERSE.EQ.0) THEN
287     SHIFT = +0.5
288     ELSE
289     SHIFT = -0.5
290     ENDIF
291 mocchiut 1.1 ELSE
292 mocchiut 1.16 IF (REVERSE.EQ.0) THEN
293     SHIFT = -0.5
294     ELSE
295     SHIFT = +0.5
296     ENDIF
297 mocchiut 1.1 ENDIF
298     C
299     IF (M.EQ.1) THEN
300     DISTX = PIANO(I) - 5.81
301 mocchiut 1.13 Y(I) = (DISTX * TANX) + AX - XALIG
302 mocchiut 1.10 c CBAR(M,I) = Y(I)
303 mocchiut 1.1 BAR(M,I) = Y(I)
304 mocchiut 1.13 CBAR(M,I) = (Y(I) + XALIG)/10.
305 mocchiut 1.16 IF (I.EQ.NPLA) MX=ABS(Y(1)-Y(NPLA))/
306     & ABS(ZIN(1)-ZIN(NPLA))
307 mocchiut 1.1 C
308     ELSE
309     DISTY = PIANO(I)
310 mocchiut 1.13 YY(I) = (DISTY * TANY) + AY - YALIG
311 mocchiut 1.10 c CBAR(M,I) = YY(I)
312 mocchiut 1.1 BAR(M,I) = YY(I)
313 mocchiut 1.13 CBAR(M,I) = (YY(I) + YALIG)/10.
314 mocchiut 1.16 IF (I.EQ.NPLA) MY=ABS(Y(1)-Y(NPLA))/
315     & ABS(ZIN(1)-ZIN(NPLA))
316 mocchiut 1.1 C
317     ENDIF
318     CALL LASTRISCIA(BAR(M,I),IBAR(M,I))
319     cibar(M,I) = ibar(m,i)
320 mocchiut 1.8 IF (ibar(m,i).EQ.-1) THEN
321     CHTRACK = CHTRACK + 1
322     ELSE
323     IWPL(M) = IWPL(M) + 1
324     ENDIF
325 mocchiut 1.1 ENDDO
326     ENDIF
327     C
328     ENDDO
329     C
330     ENDIF
331 mocchiut 1.8 C
332     IF (TRIGTY.GE.2.AND.HZN.NE.0) THEN
333     IF (GOOD2.EQ.1) THEN
334     PRINT *,' CALORIMETER - WARNING F77: unknown request'
335     GOOD2 = 1
336     GOTO 50
337     ENDIF
338     IF ( NPCFIT(1).EQ.0.OR.NPCFIT(2).EQ.0 ) THEN
339     GOOD2 = 1
340     GOTO 50
341     ENDIF
342     ENDIF
343 mocchiut 1.1 C
344     6996 CONTINUE
345     C
346     DX0L = 0.
347     C
348     C IF THE TRACK IS OUTSIDE THE CALORIMETER GO OUT, IF NOT CALCULATE DX0L
349     C
350 mocchiut 1.7 IF (CHTRACK.EQ.44) THEN ! CHTRACK is the number of planes not touched by the track
351 mocchiut 1.1 GOOD2 = 0
352 mocchiut 1.3 c PRINT *,' CALORIMETER - WARNING F77: track outside calorimeter'
353 mocchiut 1.1 GOTO 50
354     ELSE
355     IF ( IWPL(1).LE.IWPL(2) ) THEN
356     IWPLU = IWPL(1)
357     ELSE
358     IWPLU = IWPL(2)
359     ENDIF
360     C
361     DX0L = IWPLU * SQRT((BAR(2,1)-(2.66*MY+BAR(2,1)))**2
362     & + (BAR(1,1)-(2.66*MX+BAR(1,1)))**2 + 2.66**2) /
363     & 3.6
364     C
365     ENDIF
366     C
367     C
368     C RIG IS RIGIDITY AS DETERMINED BY THE TRACKER
369     C OR by CALORIMETER IF IN SELFTRIGGER MODE
370     C
371     IF (GOOD2.EQ.1) THEN
372     GTR = 1
373     IF (TRKCHI2.LT.0) GTR = 2
374     IF ( AL_PP(GTR,5).NE.0. ) THEN
375     RIG = 1./(AL_PP(GTR,5))
376     ELSE
377     GOOD2 = 0
378     PRINT *,' CALORIMETER - WARNING F77: track with R = 0'
379     GOTO 50
380     ENDIF
381     ENDIF
382 mocchiut 1.7 IF (TRIGTY.GE.2.AND.HZN.EQ.0.AND.GOOD2.EQ.0) THEN
383 mocchiut 1.1 RIG = ELEN ! SELFTRIGGER RIGIDITY
384     IF ( RIG.EQ.0. ) THEN
385 mocchiut 1.8 GOOD2 = 1
386 mocchiut 1.1 PRINT *,' CALORIMETER - WARNING F77: ST track with R = 0'
387     GOTO 50
388     ENDIF
389     ENDIF
390     C
391 mocchiut 1.7 IF (GOOD2.EQ.0.AND.(TRIGTY.LT.2.OR.HZN.EQ.1)) THEN
392     RIG = RIGINPUT
393     ENDIF
394     C
395 mocchiut 1.1 RNSS = 0.
396     QTOTT = 0.
397     C
398 mocchiut 1.15 PPLANEMAX = 1.01*(LOG(ABS(RIG)/0.0081)-1.) / 0.74
399 mocchiut 1.1 C
400     IPLANE = INT(ANINT(PPLANEMAX)) + 5
401     C
402     IF (IPLANE.GT.NPLA) IPLANE=NPLA
403     IF (IPLANE.LT.1) IPLANE = 1
404 mocchiut 1.16 c print *,' calcolo...'
405 mocchiut 1.1 C
406     C CALCULATE QLOW AND NLOW
407     C
408     DO J = IPLANE,NPLA
409     DO I = 1,NCHA
410     IF (DEXY(1,J,I).GE.EMIN) THEN
411     NLOW = NLOW + 1
412     QLOW = QLOW + DEXY(1,J,I)
413     ENDIF
414     IF (DEXY(2,J,I).GE.EMIN) THEN
415     NLOW = NLOW + 1
416     QLOW = QLOW + DEXY(2,J,I)
417     ENDIF
418     ENDDO
419     ENDDO
420     C
421     C CALCULATE QCORE AND NCORE
422     C
423     C
424     C 8 STRIPS ARE 2.88 cm , A MOLIERE RADIUS IS ABOUT 0.7 cm .
425     C
426     DO J = 1,IPLANE
427     NNX = IBAR(1,J)
428 mocchiut 1.15 RNSS = 0. ! BACO!!
429     QTOTT = 0. ! BACO!!
430 mocchiut 1.1 IF (NNX.NE.-1) THEN
431     IF (NNX.LT.9) NNX = 9
432     IF (NNX.GT.88) NNX = 88
433     INFX = NNX - 8
434     ISUPX = NNX + 8
435     DO I = INFX,ISUPX
436     IF (DEXY(1,J,I).GE.EMIN) THEN
437     RNSS = RNSS + 1
438     QTOTT = QTOTT + DEXY(1,J,I)
439     ENDIF
440     ENDDO
441     ENDIF
442     C
443     NNY = IBAR(2,J)
444     IF (NNY.NE.-1) THEN
445     IF (NNY.LT.9) NNY = 9
446     IF (NNY.GT.88) NNY = 88
447     INFY = NNY - 8
448     ISUPY = NNY + 8
449     DO I = INFY,ISUPY
450     IF (DEXY(2,J,I).GE.EMIN) THEN
451     RNSS = RNSS + 1
452     QTOTT = QTOTT + DEXY(2,J,I)
453     ENDIF
454     ENDDO
455     ENDIF
456     NCORE = RNSS * FLOAT(J) + NCORE
457     QCORE = QTOTT * FLOAT(J) + QCORE
458     ENDDO
459     C
460     C CALCULATE NOINT
461     C
462     CALL NOINTER(NIN)
463     NOINT = FLOAT(NIN)
464     C
465     C
466     C QCYL = DETECTED ENERGY AND NCYL = NUMBER OF HIT STRIPS IN A CYLINDER oF
467     C RADIUS 8.5 STRIPS WITH AXIS DEFINED BY THE DIRECTION OF THE INCOMING
468     C PARTICLE .
469     C
470     C 8 STRIPS ARE 2.88 cm , A MOLIERE RADIUS IS ABOUT 0.7 cm .
471     C
472     DO J = 1,NPLA
473     NNX = IBAR(1,J)
474     IF (NNX.NE.-1) THEN
475     IF (NNX.LT.9) NNX = 9
476     IF (NNX.GT.88) NNX = 88
477     INFX = NNX - 8
478     ISUPX = NNX + 8
479     DO I = INFX,ISUPX
480     IF (DEXY(1,J,I).LT.EMIN) GO TO 710
481     NCYL = NCYL + 1
482     QCYL = QCYL + DEXY(1,J,I)
483     710 ENDDO
484     ENDIF
485     NNY = IBAR(2,J)
486     IF (NNY.NE.-1) THEN
487     IF (NNY.LT.9) NNY = 9
488     IF (NNY.GT.88) NNY = 88
489     INFY = NNY - 8
490     ISUPY = NNY + 8
491     DO I=INFY,ISUPY
492     IF (DEXY(2,J,I).LT.EMIN) GO TO 810
493     NCYL = NCYL + 1
494     QCYL = QCYL + DEXY(2,J,I)
495     810 ENDDO
496     ENDIF
497     C
498     C QTR = DETECTED ENERGY AND NTR = NUMBER OF HIT STRIPS IN A CYLINDER oF
499     C RADIUS 4.5 STRIPS WITH AXIS DEFINED BY THE DIRECTION OF THE INCOMING
500     C PARTICLE .
501     C
502     NNX = IBAR(1,J)
503     IF (NNX.NE.-1) THEN
504     IF (NNX.LT.5) NNX = 5
505     IF (NNX.GT.92) NNX = 92
506     INFX = NNX - 4
507     ISUPX = NNX + 4
508     DO I = INFX,ISUPX
509     IF (DEXY(1,J,I).GT.EMIN) THEN
510     NTR = NTR + 1
511     QTR = QTR + DEXY(1,J,I)
512     ENDIF
513     ENDDO
514     ENDIF
515     C
516     NNY = IBAR(2,J)
517     IF (NNY.NE.-1) THEN
518     IF (NNY.LT.5) NNY = 5
519     IF (NNY.GT.92) NNY = 92
520     INFY = NNY - 4
521     ISUPY = NNY + 4
522     DO I = INFY, ISUPY
523     IF (DEXY(2,J,I).GT.EMIN) THEN
524     NTR = NTR + 1
525     QTR = QTR + DEXY(2,J,I)
526     ENDIF
527     ENDDO
528     ENDIF
529     ENDDO
530     C
531     C CALCULATE QTRACK
532     C
533     CALL LATERALE(QTRACK,RQT)
534    
535     C
536     C CALCULATE NPRESH AND QPRESH
537     C
538     DO I = 1,4
539     NNX = IBAR(1,I)
540     IF (NNX.NE.-1) THEN
541     IF (NNX.LT.3) NNX = 3
542     IF (NNX.GT.94) NNX = 94
543     INFX = NNX - 2
544     ISUPX = NNX + 2
545     DO J = INFX,ISUPX
546     IF (DEXY(1,I,J).GE.EMIN) THEN
547     NPRESH = NPRESH + 1
548     QPRESH = QPRESH + DEXY(1,I,J)
549     ENDIF
550     ENDDO
551     ENDIF
552     C
553     NNY = IBAR(2,I)
554     IF (NNY.NE.-1) THEN
555     IF (NNY.LT.3) NNY = 3
556     IF (NNY.GT.94) NNY = 94
557     INFY = NNY - 2
558     ISUPY = NNY + 2
559     DO J = INFY,ISUPY
560     IF (DEXY(2,I,J).GE.EMIN) THEN
561     NPRESH = NPRESH + 1
562     QPRESH = QPRESH + DEXY(2,I,J)
563     ENDIF
564     ENDDO
565     ENDIF
566     ENDDO
567     C
568     C CALCULATE DXTRACK, DYTRACK, QTRACKX AND QTRACKY
569     C
570     ICONTROL5 = 0
571     CALL NSHOWER(ICONTROL5,DXTRACK,DYTRACK,QTRACKX,QTRACKY)
572     C
573     C CALCULATE QPRE AND NPRE
574     C
575     DO J = 1,3
576     NNX = IBAR(1,J)
577     IF (NNX.NE.-1) THEN
578     IF (NNX.LT.9) NNX = 9
579     IF (NNX.GT.88) NNX = 88
580     INFX = NNX - 8
581     ISUPX = NNX + 8
582     DO I = INFX,ISUPX
583     IF (DEXY(1,J,I).GE.EMIN) THEN
584     NPRE = NPRE + 1
585     QPRE = QPRE + DEXY(1,J,I)
586     ENDIF
587     ENDDO
588     ENDIF
589     C
590     NNY = IBAR(2,J)
591     IF (NNY.NE.-1) THEN
592     IF (NNY.LT.9) NNY = 9
593     IF (NNY.GT.88) NNY = 88
594     INFY = NNY - 8
595     ISUPY = NNY + 8
596     DO I=INFY,ISUPY
597     IF (DEXY(2,J,I).GE.EMIN) THEN
598     NPRE = NPRE + 1
599     QPRE = QPRE + DEXY(2,J,I)
600     ENDIF
601     ENDDO
602     ENDIF
603     ENDDO
604     C
605     C CALCULATE NLAST AND QLAST
606     C
607     DO J = NPLA-4,NPLA
608     NNX = IBAR(1,J)
609     IF (NNX.NE.-1) THEN
610 mocchiut 1.2 IF (NNX.LT.5) NNX = 5
611     IF (NNX.GT.92) NNX = 92
612     c IF (NNX.LT.9) NNX = 9
613     c IF (NNX.GT.88) NNX = 88
614     INFX = NNX - 4
615     ISUPX = NNX + 4
616     c INFX = NNX - 8
617     c ISUPX = NNX + 8
618 mocchiut 1.1 DO I = INFX,ISUPX
619     IF (DEXY(1,J,I).GE.EMIN) THEN
620     NLAST = NLAST + 1
621     QLAST = QLAST + DEXY(1,J,I)
622     ENDIF
623     ENDDO
624     ENDIF
625     C
626     NNY = IBAR(2,J)
627     IF (NNY.NE.-1) THEN
628 mocchiut 1.2 IF (NNY.LT.5) NNY = 5
629     IF (NNY.GT.92) NNY = 92
630     c IF (NNY.LT.9) NNY = 9
631     c IF (NNY.GT.88) NNY = 88
632     INFY = NNY - 4
633     ISUPY = NNY + 4
634     c INFY = NNY - 8
635     c ISUPY = NNY + 8
636 mocchiut 1.1 DO I=INFY,ISUPY
637     IF (DEXY(2,J,I).GE.EMIN) THEN
638     NLAST = NLAST + 1
639     QLAST = QLAST + DEXY(2,J,I)
640     ENDIF
641     ENDDO
642     ENDIF
643     ENDDO
644     C
645     C
646     C CALCULATE PLANETOT AND QMEAN
647     C
648     DO M = 1,2
649     RPIANO(M) = 0.
650     NTOT(M) = 0
651     ENDDO
652     NPIANI = 5
653     QMEAN = 0.
654     INDEX = 0
655 mocchiut 1.12 C
656     IF (TRIGTY.GE.2.AND.HZN.NE.0) THEN
657     EINF = 50.
658     ESUP = 15000.
659     CALL NUCLEI(RPIANO,NPIANI,QMEAN,NTOT,INDEX)
660     PLANETOT = RPIANO(1) + RPIANO(2)
661     ELSE
662     EINF = EMIN
663     ESUP = 15000.
664     CALL ELIO(RPIANO,NPIANI,QMEAN,NTOT,INDEX)
665     PLANETOT = RPIANO(1) + RPIANO(2)
666     ENDIF
667 mocchiut 1.1 C
668     50 CONTINUE
669     C
670 mocchiut 1.16 c print *,' esco'
671 mocchiut 1.1 RETURN
672     END
673    
674    
  ViewVC Help
Powered by ViewVC 1.1.23