/[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.12 - (hide annotations) (download)
Tue Jan 23 11:52:26 2007 UTC (17 years, 10 months ago) by mocchiut
Branch: MAIN
Changes since 1.11: +14 -11 lines
OrbitalInfo updated, small changes in calo code

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

  ViewVC Help
Powered by ViewVC 1.1.23