1 | * Program to calculate differenced data for forward scattering |
---|
2 | CHARACTER*80 runid |
---|
3 | character runo*4,run1*5,fout*80 |
---|
4 | REAL T(10000),CT(10000),err(10000),SUM1(9),SUM2(9),SUM3(9) |
---|
5 | REAL COUT(10000),CTHIN(10000),CTHICK(10000),NFR(9),sum4(9) |
---|
6 | REAL ERROUT(10000),ERRTHIN(10000),ERRTHICK(10000) |
---|
7 | REAL OUT(2000),THIN(2000),THICK(2000) |
---|
8 | REAL MOUT(10000),MTHIN(10000),MTHICK(10000) |
---|
9 | INTEGER IARRAY(100),ERRCODE,LOUT,NCHAN,NDET,IFSN,NOS |
---|
10 | INTEGER rpb1(32),ICT(10000),ICTA(9,10000),MICTA(9,10000),NI(9) |
---|
11 | real xin(10000),yin(10000),ein(10000) |
---|
12 | real xout(10000),yout(10000),eout(10000) |
---|
13 | |
---|
14 | write(6,*) ' First and last run numbers?' |
---|
15 | read(5,*) IRUN1,IRUN2 |
---|
16 | write(6,*) ' First and last spectrum numbers?' |
---|
17 | READ(5,*) ISPECT1,ISPECT2 |
---|
18 | * WRITE(6,*) ' File name for output?' |
---|
19 | * read(5,'(a)') fout |
---|
20 | fout='temp.dat' |
---|
21 | open(unit=4,file=fout,status='new') |
---|
22 | |
---|
23 | DO 1000 ISPECT=ISPECT1,ISPECT2 |
---|
24 | WRITE(6,*) ' spectrum: ',ISPECT |
---|
25 | IOPT=1 |
---|
26 | IF(IOPT.EQ.1) THEN |
---|
27 | WRITE(6,*) ' THIN DIFFERENCE' |
---|
28 | ELSE IF(IOPT.EQ.2) THEN |
---|
29 | WRITE(6,*) ' DOUBLE DIFFERENCE' |
---|
30 | BETA=ST_VAR_IN(5) |
---|
31 | WRITE(6,*) ' BETA=',BETA |
---|
32 | ELSE IF(IOPT.EQ.3) THEN |
---|
33 | WRITE(6,*) ' THICK DIFFERENCE' |
---|
34 | END IF |
---|
35 | |
---|
36 | DO I=1,9 |
---|
37 | DO J=1,10000 |
---|
38 | ICTA(I,J)=0 |
---|
39 | MICTA(I,J)=0 |
---|
40 | END DO |
---|
41 | END DO |
---|
42 | NGOOD_FRAMES=0 |
---|
43 | |
---|
44 | DO 101 IRUN=IRUN1,IRUN2 |
---|
45 | * Define name of data file. |
---|
46 | if(irun.lt.10000) then |
---|
47 | write(runo(1:4),'(I4.4)')irun |
---|
48 | runid='evs_data:evs0'//runo//'.raw' |
---|
49 | else |
---|
50 | write(run1(1:5),'(I5.5)')irun |
---|
51 | runid='evs_data:evs'//run1//'.raw' |
---|
52 | end if |
---|
53 | * write(6,5) runid |
---|
54 | 5 Format(' ',' data file is;',50a) |
---|
55 | |
---|
56 | * Get relevant run information. |
---|
57 | CALL GETPARI(runid,'RPB',rpb1,32,length_out,errcode) |
---|
58 | ngood_frames =ngood_frames+rpb1(10) |
---|
59 | CALL GETPARI(runid,'NPER',NPD,100,LOUT,ERRCODE) |
---|
60 | CALL GETPARI(RUNID,'NSP1',IARRAY,100,LOUT,ERRCODE) |
---|
61 | NDET=IARRAY(1) |
---|
62 | |
---|
63 | CALL GETPARR(RUNID,'TCB1',T,10000,LOUT,ERRCODE) |
---|
64 | CALL GETPARI(RUNID,'NTC1',IARRAY,100,LOUT,ERRCODE) |
---|
65 | NCHAN=IARRAY(1) |
---|
66 | NMON=MON(RUNID) |
---|
67 | DT=T(3)-T(2) |
---|
68 | DO I=1,NCHAN |
---|
69 | T(I)=T(I)-DT |
---|
70 | END DO |
---|
71 | |
---|
72 | DO IP=1,NPD |
---|
73 | IS=ISPECT+(IP-1)*(NDET+1) |
---|
74 | CALL GETDAT(RUNID,IS,1,ICT,10000,ERRCODE) |
---|
75 | DO IT=2,NCHAN |
---|
76 | ICTA(IP,IT)=ICTA(IP,IT)+ICT(IT) |
---|
77 | END DO |
---|
78 | IS=NMON+(IP-1)*(NDET+1) |
---|
79 | CALL GETDAT(RUNID,IS,1,ICT,10000,ERRCODE) |
---|
80 | DO IT=2,NCHAN |
---|
81 | MICTA(IP,IT)=MICTA(IP,IT)+ICT(IT) |
---|
82 | END DO |
---|
83 | END DO |
---|
84 | |
---|
85 | 101 CONTINUE |
---|
86 | |
---|
87 | write(6,*) ' Total number of good frames=',ngood_frames |
---|
88 | WRITE(6,*) ' NUMBER OF PERIODS=',NPD |
---|
89 | WRITE(6,*) ' NUMBER OF SPECTRA=',NDET |
---|
90 | WRITE(6,*) ' Spectrum number=',ISPECT |
---|
91 | WRITE(6,*) ' Number of monitor spectrum=',nmon |
---|
92 | * DETERMINE WHICH SPECTRA IS OUT,THIN,THICK |
---|
93 | IOUT=1 |
---|
94 | ITHIN=2 |
---|
95 | ITHICK=3 |
---|
96 | |
---|
97 | DO IP=1,NPD |
---|
98 | DO IT=2,NCHAN |
---|
99 | IF(T(IT).GT.100.0.AND.T(IT).LT.400.0) THEN |
---|
100 | SUM1(IP)=SUM1(IP)+ICTA(IP,IT) |
---|
101 | ELSE IF(T(IT).GT.550.0.AND.T(IT).LT.600.0) THEN |
---|
102 | SUM2(IP)=SUM2(IP)+ICTA(IP,IT) |
---|
103 | END IF |
---|
104 | END DO |
---|
105 | write(6,*) ' sum1=',sum1(ip),' sum2=',sum2(ip) |
---|
106 | END DO |
---|
107 | |
---|
108 | DO IP=1,NPD |
---|
109 | NI(IP)=IP |
---|
110 | IF(SUM2(IP).NE.0.0) THEN |
---|
111 | SUM1(IP)=SUM1(IP)/SUM2(IP) |
---|
112 | SUM4(IP)=SUM1(IP) |
---|
113 | END IF |
---|
114 | END DO |
---|
115 | |
---|
116 | CALL ORDER(SUM1,NI,NPD) |
---|
117 | CALL WRITESP(NI,ISPECT,NDET,NPD) ! Write, out,thin,thick spectra. |
---|
118 | |
---|
119 | * WRITE(6,*) ' 1' |
---|
120 | * Calculate total foil out, thin, thick and convert to Cts/micro-sec |
---|
121 | |
---|
122 | IF(ISPECT.GE.135.AND.ISPECT.LE.142) THEN |
---|
123 | NI(1)=ISPECT |
---|
124 | NI(2)=ISPECT+2.0*(NDET+1) |
---|
125 | NI(3)=ISPECT+4.0*(NDET+1) |
---|
126 | NI(4)=ISPECT+1.0*(NDET+1) |
---|
127 | NI(5)=ISPECT+3.0*(NDET+1) |
---|
128 | NI(6)=ISPECT+5.0*(NDET+1) |
---|
129 | ELSE IF(ISPECT.GE.143.AND.ISPECT.LE.150) THEN |
---|
130 | NI(1)=ISPECT+1.0*(NDET+1) |
---|
131 | NI(2)=ISPECT+3.0*(NDET+1) |
---|
132 | NI(3)=ISPECT+5.0*(NDET+1) |
---|
133 | NI(4)=ISPECT |
---|
134 | NI(5)=ISPECT+2.0*(NDET+1) |
---|
135 | NI(6)=ISPECT+4.0*(NDET+1) |
---|
136 | ELSE IF(ISPECT.GE.151.AND.ISPECT.LE.158) THEN |
---|
137 | NI(1)=ISPECT |
---|
138 | NI(2)=ISPECT+2.0*(NDET+1) |
---|
139 | NI(3)=ISPECT+4.0*(NDET+1) |
---|
140 | NI(4)=ISPECT+1.0*(NDET+1) |
---|
141 | NI(5)=ISPECT+3.0*(NDET+1) |
---|
142 | NI(6)=ISPECT+5.0*(NDET+1) |
---|
143 | ELSE IF(ISPECT.GE.159.AND.ISPECT.LE.166) THEN |
---|
144 | NI(1)=ISPECT+1.0*(NDET+1) |
---|
145 | NI(2)=ISPECT+3.0*(NDET+1) |
---|
146 | NI(3)=ISPECT+5.0*(NDET+1) |
---|
147 | NI(4)=ISPECT |
---|
148 | NI(5)=ISPECT+2.0*(NDET+1) |
---|
149 | NI(6)=ISPECT+4.0*(NDET+1) |
---|
150 | ELSE IF(ISPECT.GE.167.AND.ISPECT.LE.174) THEN |
---|
151 | NI(1)=ISPECT |
---|
152 | NI(2)=ISPECT+2.0*(NDET+1) |
---|
153 | NI(3)=ISPECT+4.0*(NDET+1) |
---|
154 | NI(4)=ISPECT+1.0*(NDET+1) |
---|
155 | NI(5)=ISPECT+3.0*(NDET+1) |
---|
156 | NI(6)=ISPECT+5.0*(NDET+1) |
---|
157 | ELSE IF(ISPECT.GE.175.AND.ISPECT.LE.182) THEN |
---|
158 | NI(1)=ISPECT+1.0*(NDET+1) |
---|
159 | NI(2)=ISPECT+3.0*(NDET+1) |
---|
160 | NI(3)=ISPECT+5.0*(NDET+1) |
---|
161 | NI(4)=ISPECT |
---|
162 | NI(5)=ISPECT+2.0*(NDET+1) |
---|
163 | NI(6)=ISPECT+4.0*(NDET+1) |
---|
164 | ELSE IF(ISPECT.GE.183.AND.ISPECT.LE.190) THEN |
---|
165 | NI(1)=ISPECT |
---|
166 | NI(2)=ISPECT+2.0*(NDET+1) |
---|
167 | NI(3)=ISPECT+4.0*(NDET+1) |
---|
168 | NI(4)=ISPECT+1.0*(NDET+1) |
---|
169 | NI(5)=ISPECT+3.0*(NDET+1) |
---|
170 | NI(6)=ISPECT+5.0*(NDET+1) |
---|
171 | ELSE IF(ISPECT.GE.191.AND.ISPECT.LE.198) THEN |
---|
172 | NI(1)=ISPECT+1.0*(NDET+1) |
---|
173 | NI(2)=ISPECT+3.0*(NDET+1) |
---|
174 | NI(3)=ISPECT+5.0*(NDET+1) |
---|
175 | NI(4)=ISPECT |
---|
176 | NI(5)=ISPECT+2.0*(NDET+1) |
---|
177 | NI(6)=ISPECT+4.0*(NDET+1) |
---|
178 | END IF |
---|
179 | WRITE(6,*) ' FOIL OUT SPECTRA ARE' |
---|
180 | WRITE(6,*) NI(4),NI(5),NI(6) |
---|
181 | WRITE(6,*) ' FOIL IN SPECTRA ARE' |
---|
182 | WRITE(6,*) NI(1),NI(2),NI(3) |
---|
183 | |
---|
184 | |
---|
185 | IF(ISPECT.GE.135.AND.ISPECT.LE.142) THEN |
---|
186 | NI(1)=1 |
---|
187 | NI(2)=3 |
---|
188 | NI(3)=5 |
---|
189 | NI(4)=2 |
---|
190 | NI(5)=4 |
---|
191 | NI(6)=6 |
---|
192 | ELSE IF(ISPECT.GE.143.AND.ISPECT.LE.150) THEN |
---|
193 | NI(1)=2 |
---|
194 | NI(2)=4 |
---|
195 | NI(3)=6 |
---|
196 | NI(4)=1 |
---|
197 | NI(5)=3 |
---|
198 | NI(6)=5 |
---|
199 | ELSE IF(ISPECT.GE.151.AND.ISPECT.LE.158) THEN |
---|
200 | NI(1)=1 |
---|
201 | NI(2)=3 |
---|
202 | NI(3)=5 |
---|
203 | NI(4)=2 |
---|
204 | NI(5)=4 |
---|
205 | NI(6)=6 |
---|
206 | ELSE IF(ISPECT.GE.159.AND.ISPECT.LE.166) THEN |
---|
207 | NI(1)=2 |
---|
208 | NI(2)=4 |
---|
209 | NI(3)=6 |
---|
210 | NI(4)=1 |
---|
211 | NI(5)=3 |
---|
212 | NI(6)=5 |
---|
213 | ELSE IF(ISPECT.GE.167.AND.ISPECT.LE.174) THEN |
---|
214 | NI(1)=1 |
---|
215 | NI(2)=3 |
---|
216 | NI(3)=5 |
---|
217 | NI(4)=2 |
---|
218 | NI(5)=4 |
---|
219 | NI(6)=6 |
---|
220 | ELSE IF(ISPECT.GE.175.AND.ISPECT.LE.182) THEN |
---|
221 | NI(1)=2 |
---|
222 | NI(2)=4 |
---|
223 | NI(3)=6 |
---|
224 | NI(4)=1 |
---|
225 | NI(5)=3 |
---|
226 | NI(6)=5 |
---|
227 | ELSE IF(ISPECT.GE.183.AND.ISPECT.LE.190) THEN |
---|
228 | NI(1)=1 |
---|
229 | NI(2)=3 |
---|
230 | NI(3)=5 |
---|
231 | NI(4)=2 |
---|
232 | NI(5)=4 |
---|
233 | NI(6)=6 |
---|
234 | ELSE IF(ISPECT.GE.191.AND.ISPECT.LE.198) THEN |
---|
235 | NI(1)=2 |
---|
236 | NI(2)=4 |
---|
237 | NI(3)=6 |
---|
238 | NI(4)=1 |
---|
239 | NI(5)=3 |
---|
240 | NI(6)=5 |
---|
241 | END IF |
---|
242 | ! WRITE(6,*) ' FOIL OUT SPECTRA ARE' |
---|
243 | ! WRITE(6,*) NI(4),NI(5),NI(6) |
---|
244 | ! WRITE(6,*) ' FOIL IN SPECTRA ARE' |
---|
245 | ! WRITE(6,*) NI(1),NI(2),NI(3) |
---|
246 | |
---|
247 | DO I=2,NCHAN |
---|
248 | COUT(I)=FLOAT(ICTA(NI(4),I)+ICTA(NI(5),I)+ICTA(NI(6),I)) |
---|
249 | $/(T(I)-T(I-1)) |
---|
250 | ERROUT(I)=SQRT(FLOAT(ICTA(NI(4),I)+ICTA(NI(5),I)+ICTA(NI(6),I))) |
---|
251 | $/(T(I)-T(I-1)) |
---|
252 | CTHIN(I)=FLOAT(ICTA(NI(1),I)+ICTA(NI(2),I)+ICTA(NI(3),I) |
---|
253 | $)/(T(I)-T(I-1)) |
---|
254 | ERRTHIN(I)=SQRT(FLOAT(ICTA(NI(1),I)+ICTA(NI(2),I) |
---|
255 | $ +ICTA(NI(3),I)))/(T(I)-T(I-1)) |
---|
256 | CTHICK(I)=FLOAT(ICTA(NI(1),I)+ICTA(NI(2),I) |
---|
257 | $)/(T(I)-T(I-1)) |
---|
258 | ERRTHICK(I)=SQRT(FLOAT(ICTA(NI(1),I)+ICTA(NI(2),I) |
---|
259 | $))/(T(I)-T(I-1)) |
---|
260 | END DO |
---|
261 | SUM3(IOUT)=SUM2(NI(4))+SUM2(NI(5))+SUM2(NI(6)) |
---|
262 | SUM3(ITHIN)=SUM2(NI(1))+SUM2(NI(2))+SUM2(NI(3)) |
---|
263 | SUM3(ITHICK)=SUM2(NI(1))+SUM2(NI(2)) |
---|
264 | |
---|
265 | * Calculate total foil out, thin, thick monitor cts |
---|
266 | IF(NPD.EQ.2) THEN |
---|
267 | DO I=2,NCHAN |
---|
268 | MOUT(I)=FLOAT(MICTA(NI(2),I))/(T(I)-T(I-1)) |
---|
269 | MTHIN(I)=FLOAT(MICTA(NI(1),I))/(T(I)-T(I-1)) |
---|
270 | END DO |
---|
271 | ELSE IF(NPD.EQ.3) THEN |
---|
272 | DO I=2,NCHAN |
---|
273 | MOUT(I)=FLOAT(MICTA(NI(3),I))/(T(I)-T(I-1)) |
---|
274 | MTHIN(I)=FLOAT(MICTA(NI(2),I))/(T(I)-T(I-1)) |
---|
275 | MTHICK(I)=FLOAT(MICTA(NI(1),I))/(T(I)-T(I-1)) |
---|
276 | END DO |
---|
277 | ELSE IF(NPD.EQ.6) THEN |
---|
278 | DO I=2,NCHAN |
---|
279 | MOUT(I)=FLOAT(MICTA(NI(5),I)+MICTA(NI(6),I) |
---|
280 | $)/(T(I)-T(I-1)) |
---|
281 | MTHIN(I)=FLOAT(MICTA(NI(3),I)+MICTA(NI(4),I) |
---|
282 | $)/(T(I)-T(I-1)) |
---|
283 | MTHICK(I)=FLOAT(MICTA(NI(1),I)+MICTA(NI(2),I) |
---|
284 | $)/(T(I)-T(I-1)) |
---|
285 | END DO |
---|
286 | ELSE IF(NPD.EQ.9) THEN |
---|
287 | DO I=2,NCHAN |
---|
288 | MOUT(I)=FLOAT(MICTA(NI(7),I)+MICTA(NI(8),I) |
---|
289 | $ +MICTA(NI(9),I))/(T(I)-T(I-1)) |
---|
290 | MTHIN(I)=FLOAT(MICTA(NI(4),I)+MICTA(NI(5),I) |
---|
291 | $ +MICTA(NI(6),I))/(T(I)-T(I-1)) |
---|
292 | MTHICK(I)=FLOAT(MICTA(NI(1),I)+MICTA(NI(2),I) |
---|
293 | $ +MICTA(NI(3),I))/(T(I)-T(I-1)) |
---|
294 | END DO |
---|
295 | END IF |
---|
296 | |
---|
297 | |
---|
298 | * Calculate number of frames in each period. |
---|
299 | SUM3T=0.0 |
---|
300 | DO I=1,3 |
---|
301 | SUM3T=SUM3T+SUM3(I) |
---|
302 | END DO |
---|
303 | DO I=1,3 |
---|
304 | if(sum3t.ne.0.0)NFR(I)=NGOOD_FRAMES*SUM3(I)/SUM3T |
---|
305 | IF(NFR(I).EQ.0.0) THEN |
---|
306 | write(6,*) ' nframes=0 in period',i |
---|
307 | nfr(i)=1 |
---|
308 | end if |
---|
309 | END DO |
---|
310 | |
---|
311 | |
---|
312 | * Start of routines for dead-time correction |
---|
313 | * CONVERT TO COUNTS/FRAME |
---|
314 | * DO I=2,NCHAN |
---|
315 | * COUT(I)=COUT(I)/NFR(IOUT) |
---|
316 | * ERROUT(I)=ERROUT(I)/NFR(IOUT) |
---|
317 | * CTHIN(I)=CTHIN(I)/NFR(ITHIN) |
---|
318 | * ERRTHIN(I)=ERRTHIN(I)/NFR(ITHIN) |
---|
319 | * IF(NPD.NE.2) THEN |
---|
320 | * CTHICK(I)=CTHICK(I)/NFR(ITHICK) |
---|
321 | * ERRTHICK(I)=ERRTHICK(I)/NFR(ITHICK) |
---|
322 | * END IF |
---|
323 | * END DO |
---|
324 | * |
---|
325 | * Correct for dead times in detector |
---|
326 | * TAU=1.0 |
---|
327 | * DO I=2,NCHAN |
---|
328 | * COUT(I)=COUT(I)/(1.0-TAU*COUT(I)) |
---|
329 | * ERROUT(I)=ERROUT(I)/(1.0-TAU*ERROUT(I)) |
---|
330 | * CTHIN(I)=CTHIN(I)/(1.0-TAU*CTHIN(I)) |
---|
331 | * ERRTHIN(I)=ERRTHIN(I)/(1.0-TAU*ERRTHIN(I)) |
---|
332 | * IF(NPD.NE.2) THEN |
---|
333 | * ERRTHICK(I)=ERRTHICK(I)/(1.0-TAU*ERRTHICK(I)) |
---|
334 | * END IF |
---|
335 | * END DO |
---|
336 | * |
---|
337 | * CONVERT BACK TO TOTAL COUNTS. |
---|
338 | * DO I=2,NCHAN |
---|
339 | * COUT(I)=COUT(I)*NFR(IOUT) |
---|
340 | * ERROUT(I)=ERROUT(I)*NFR(IOUT) |
---|
341 | * CTHIN(I)=CTHIN(I)*NFR(ITHIN) |
---|
342 | * ERRTHIN(I)=ERRTHIN(I)*NFR(ITHIN) |
---|
343 | * IF(NPD.NE.2) THEN |
---|
344 | * CTHICK(I)=CTHICK(I)*NFR(ITHICK) |
---|
345 | * ERRTHICK(I)=ERRTHICK(I)*NFR(ITHICK) |
---|
346 | * END IF |
---|
347 | * END DO |
---|
348 | |
---|
349 | |
---|
350 | |
---|
351 | * Calculate monitor cts between 600 and 700 usec. |
---|
352 | RMOUT=0.0 |
---|
353 | RMTHIN=0.0 |
---|
354 | RMTHICK=0.0 |
---|
355 | DO I=2,NCHAN |
---|
356 | IF(T(I).GE.600.0.AND.T(I).LT.700.0) THEN |
---|
357 | RMOUT=RMOUT+MOUT(I) |
---|
358 | RMTHIN=RMTHIN+MTHIN(I) |
---|
359 | RMTHICK=RMTHICK+MTHICK(I) |
---|
360 | * WRITE(6,*) ' T=',T(I),' MOUT=',MOUT(I) |
---|
361 | END IF |
---|
362 | END DO |
---|
363 | if(rmthick.eq.0.0)rmthick=1e-6 |
---|
364 | |
---|
365 | * WRITE(6,*) ' RMOUT=',RMOUT,' RMTHIN=',RMTHIN,' RMTHICK=',RMTHICK |
---|
366 | |
---|
367 | * Normalise spectra to monitor counter. |
---|
368 | const=1000.0 ! Constant to give suitable scale in tof spectra. |
---|
369 | DO I=2,NCHAN |
---|
370 | COUT(I)=COUT(I)*const/RMOUT |
---|
371 | CTHIN(I)=CTHIN(I)*CONST/RMTHIN |
---|
372 | CTHICK(I)=CTHICK(I)*CONST/RMTHICK |
---|
373 | ERROUT(I)=ERROUT(I)*CONST/RMOUT |
---|
374 | ERRTHIN(I)=ERRTHIN(I)*CONST/RMTHIN |
---|
375 | ERRTHICK(I)=ERRTHICK(I)*CONST/RMTHICK |
---|
376 | END DO |
---|
377 | |
---|
378 | |
---|
379 | * Calculate cts between 400 and 450 usec. |
---|
380 | SUMOUT=0.0 |
---|
381 | SUMTHIN=0.0 |
---|
382 | SUMTHICK=0.0 |
---|
383 | DO I=2,NCHAN |
---|
384 | IF(T(I).GE.400.0.AND.T(I).LT.450.0) THEN |
---|
385 | SUMOUT=SUMOUT+COUT(I) |
---|
386 | SUMTHIN=SUMTHIN+CTHIN(I) |
---|
387 | SUMTHICK=SUMTHICK+CTHICK(I) |
---|
388 | END IF |
---|
389 | END DO |
---|
390 | |
---|
391 | |
---|
392 | if(sumthin.eq.0.0) then |
---|
393 | write(6,*) ' no cts in thin foil spectra' |
---|
394 | sumthin=1.0 |
---|
395 | end if |
---|
396 | if(sumout.eq.0.0) then |
---|
397 | write(6,*) ' no cts in foil out spectra' |
---|
398 | sumout=1.0 |
---|
399 | end if |
---|
400 | if(sumthick.eq.0.0) then |
---|
401 | write(6,*) ' no cts in thick foil spectra' |
---|
402 | sumthick=1.0 |
---|
403 | end if |
---|
404 | |
---|
405 | * Normalise spectra to same area as foil out between 400 and 450 usec. |
---|
406 | DO I=2,NCHAN |
---|
407 | CTHIN(I)=CTHIN(I)*SUMOUT/SUMTHIN |
---|
408 | ERRTHIN(I)=ERRTHIN(I)*SUMOUT/SUMTHIN |
---|
409 | CTHICK(I)=CTHICK(I)*SUMOUT/SUMTHICK |
---|
410 | ERRTHICK(I)=ERRTHICK(I)*SUMOUT/SUMTHICK |
---|
411 | END DO |
---|
412 | |
---|
413 | IT=0.0 |
---|
414 | IF(IOPT.EQ.1)THEN ! THIN DIFFERENCE |
---|
415 | DT=T(3)-T(2) |
---|
416 | DO I=2,NCHAN |
---|
417 | IF(T(I).LT.600.0.and.t(i).gt.1.0) THEN |
---|
418 | IT=IT+1 |
---|
419 | XOUT(IT)=T(I) |
---|
420 | YOUT(IT)=-COUT(I)+CTHIN(I) |
---|
421 | EOUT(IT)=SQRT(ERROUT(I)**2+ERRTHIN(I)**2) |
---|
422 | END IF |
---|
423 | END DO |
---|
424 | ELSE IF(IOPT.EQ.2) THEN ! DOUBLE DIFFERENCE |
---|
425 | DO I=2,NCHAN |
---|
426 | IF(T(I).LT.600.0) THEN |
---|
427 | IT=IT+1 |
---|
428 | XOUT(IT)=T(I) |
---|
429 | YOUT(IT)=COUT(I)*(1-BETA)-CTHIN(I)+BETA*CTHICK(I) |
---|
430 | EOUT(I)=SQRT((1-BETA)**2*ERROUT(I)**2+ERRTHIN(I)**2 |
---|
431 | $ +BETA**2*ERRTHICK(I)**2) |
---|
432 | END IF |
---|
433 | END DO |
---|
434 | ELSE IF(IOPT.EQ.3) THEN ! THICK DIFFERENCE |
---|
435 | DO I=2,NCHAN |
---|
436 | IF(T(I).LT.600.0) THEN |
---|
437 | IT=IT+1 |
---|
438 | XOUT(IT)=T(I) |
---|
439 | YOUT(IT)=COUT(I)-CTHICK(I) |
---|
440 | EOUT(IT)=SQRT(ERROUT(I)**2+ERRTHICK(I)**2) |
---|
441 | END IF |
---|
442 | END DO |
---|
443 | END IF |
---|
444 | |
---|
445 | |
---|
446 | |
---|
447 | LPTOUT=IT-1 |
---|
448 | write(4,*) ispect,lptout |
---|
449 | do i=1,lptout |
---|
450 | if(t(i).lt.600.0) then |
---|
451 | write(4,*) xout(i),yout(i),eout(i) |
---|
452 | end if |
---|
453 | end do |
---|
454 | |
---|
455 | 1000 Continue |
---|
456 | write(4,*) irun1,irun2 |
---|
457 | close(4) |
---|
458 | write(6,*) ' Data is in TEMP.DAT' |
---|
459 | END |
---|
460 | |
---|
461 | c Sort in order of increasing x values. |
---|
462 | SUBROUTINE ORDER(X,Y,N) |
---|
463 | REAL X(N) |
---|
464 | INTEGER Y(N) |
---|
465 | |
---|
466 | DO IS=1,N |
---|
467 | IMIN=IS |
---|
468 | XMIN=X(IS) |
---|
469 | DO I=IS,N |
---|
470 | IF(X(I).LE.XMIN) THEN |
---|
471 | IMIN=I |
---|
472 | XMIN=X(I) |
---|
473 | END IF |
---|
474 | END DO |
---|
475 | * Put smallest x value and corresponding y value in X(IS) and Y(IS). |
---|
476 | XTEMP=X(IMIN) |
---|
477 | YTEMP=Y(IMIN) |
---|
478 | X(IMIN)=X(IS) |
---|
479 | Y(IMIN)=Y(IS) |
---|
480 | X(IS)=XTEMP |
---|
481 | Y(IS)=YTEMP |
---|
482 | END DO |
---|
483 | |
---|
484 | END |
---|
485 | |
---|
486 | SUBROUTINE WRITESP(NI,ISPECT,NDET,NPD) |
---|
487 | INTEGER NI(9) |
---|
488 | If(npd.eq.9) then |
---|
489 | write(6,*) ' Foil out Spectra are:' |
---|
490 | do i=7,9 |
---|
491 | IV=ISPECT+(ni(i)-1)*(NDET+1) |
---|
492 | write(6,*) iv |
---|
493 | end do |
---|
494 | |
---|
495 | write(6,*) ' Thin Foil Spectra are:' |
---|
496 | do i=4,6 |
---|
497 | IV=ISPECT+(ni(i)-1)*(NDET+1) |
---|
498 | write(6,*) iv |
---|
499 | end do |
---|
500 | |
---|
501 | write(6,*) ' Thick foil Spectra are:' |
---|
502 | do i=1,3 |
---|
503 | IV=ISPECT+(ni(i)-1)*(NDET+1) |
---|
504 | write(6,*) iv |
---|
505 | end do |
---|
506 | end if |
---|
507 | |
---|
508 | If(npd.eq.3) then |
---|
509 | write(6,*) ' Foil out Spectrum is:' |
---|
510 | IV=ISPECT+(ni(3)-1)*(NDET+1) |
---|
511 | write(6,*) iv |
---|
512 | write(6,*) ' Thin foil Spectrum is:' |
---|
513 | IV=ISPECT+(ni(2)-1)*(NDET+1) |
---|
514 | write(6,*) iv |
---|
515 | write(6,*) ' Thick foil Spectrum is:' |
---|
516 | IV=ISPECT+(ni(1)-1)*(NDET+1) |
---|
517 | write(6,*) iv |
---|
518 | end if |
---|
519 | |
---|
520 | If(npd.eq.2) then |
---|
521 | write(6,*) ' Foil out Spectrum is:' |
---|
522 | IV=ISPECT+(ni(2)-1)*(NDET+1) |
---|
523 | write(6,*) iv |
---|
524 | write(6,*) ' Thin foil Spectrum is:' |
---|
525 | IV=ISPECT+(ni(1)-1)*(NDET+1) |
---|
526 | write(6,*) iv |
---|
527 | end if |
---|
528 | RETURN |
---|
529 | END |
---|
530 | |
---|
531 | * FUNCTION TO READ NUMBER OF MONITOR COUNTER |
---|
532 | function mon(runid) |
---|
533 | integer nmon*4,nrun*4 |
---|
534 | CHARACTER runid*80 |
---|
535 | |
---|
536 | CALL GETPARI(runid,'MDET',NMON,100,LOUT,ERRCODE) |
---|
537 | MON=NMON |
---|
538 | END |
---|