source: trunk/MagicSoft/Mars/msignal/MExtractPINDiode.cc@ 4829

Last change on this file since 4829 was 3962, checked in by gaug, 21 years ago
*** empty log message ***
File size: 7.6 KB
Line 
1/* ======================================================================== *\
2!
3! *
4! * This file is part of MARS, the MAGIC Analysis and Reconstruction
5! * Software. It is distributed to you in the hope that it can be a useful
6! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
7! * It is distributed WITHOUT ANY WARRANTY.
8! *
9! * Permission to use, copy, modify and distribute this software and its
10! * documentation for any purpose is hereby granted without fee,
11! * provided that the above copyright notice appear in all copies and
12! * that both that copyright notice and this permission notice appear
13! * in supporting documentation. It is provided "as is" without express
14! * or implied warranty.
15! *
16!
17!
18! Author(s): Markus Gaug, 02/2004 <mailto:markus@ifae.es>
19!
20! Copyright: MAGIC Software Development, 2000-2004
21!
22!
23\* ======================================================================== */
24
25//////////////////////////////////////////////////////////////////////////////
26//
27// MExtractPINDiode
28//
29// Extracts the signal from a fixed window in a given range.
30//
31// Call: SetRange(fHiGainFirst, fHiGainLast, fLoGainFirst, fLoGainLast)
32// to modify the ranges. Ranges have to be an even number. In case of odd
33// ranges, the last slice will be reduced by one.
34// Defaults are:
35//
36// fHiGainFirst = fgHiGainFirst = 3
37// fHiGainLast = fgHiGainLast = 14
38// fLoGainFirst = fgLoGainFirst = 3
39// fLoGainLast = fgLoGainLast = 14
40//
41//////////////////////////////////////////////////////////////////////////////
42#include "MExtractPINDiode.h"
43
44#include <fstream>
45
46#include "MLog.h"
47#include "MLogManip.h"
48
49#include "MParList.h"
50
51#include "MRawEvtData.h"
52#include "MRawEvtPixelIter.h"
53
54#include "MPedestalCam.h"
55#include "MPedestalPix.h"
56
57#include "MExtractedSignalPINDiode.h"
58
59ClassImp(MExtractPINDiode);
60
61using namespace std;
62
63const UInt_t MExtractPINDiode::fgPINDiodeIdx = 100;
64const Byte_t MExtractPINDiode::fgHiGainFirst = 0;
65const Byte_t MExtractPINDiode::fgHiGainLast = 14;
66const Byte_t MExtractPINDiode::fgLoGainFirst = 0;
67const Byte_t MExtractPINDiode::fgLoGainLast = 14;
68// --------------------------------------------------------------------------
69//
70// Default constructor.
71//
72// Calls:
73// - SetRange(fgHiGainFirst, fgHiGainLast, fgLoGainFirst, fgLoGainLast)
74// - SetPINDiodeIdx()
75//
76MExtractPINDiode::MExtractPINDiode(const char *name, const char *title)
77{
78
79 fName = name ? name : "MExtractPINDiode";
80 fTitle = title ? title : "Task to extract the signal from the FADC slices";
81
82 SetRange(fgHiGainFirst, fgHiGainLast, fgLoGainFirst, fgLoGainLast);
83 SetPINDiodeIdx();
84}
85
86// --------------------------------------------------------------------------
87//
88// SetRange:
89//
90// Checks:
91// - if the window defined by (fHiGainLast-fHiGainFirst-1) are odd, subtract one
92// - if the window defined by (fLoGainLast-fLoGainFirst-1) are odd, subtract one
93// - if the Hi Gain window is smaller than 2, set fHiGainLast to fHiGainFirst+1
94// - if the Lo Gain window is smaller than 2, set fLoGainLast to fLoGainFirst+1
95//
96// Calls:
97// - MExtractor::SetRange(hifirst,hilast,lofirst,lolast);
98//
99// Sets:
100// - fNumHiGainSamples to: (Float_t)(fHiGainLast-fHiGainFirst+1)
101// - fNumLoGainSamples to: (Float_t)(fLoGainLast-fLoGainFirst+1)
102// - fSqrtHiGainSamples to: TMath::Sqrt(fNumHiGainSamples)
103// - fSqrtLoGainSamples to: TMath::Sqrt(fNumLoGainSamples)
104//
105void MExtractPINDiode::SetRange(Byte_t hifirst, Byte_t hilast, Byte_t lofirst, Byte_t lolast)
106{
107
108 const Byte_t window = hilast-hifirst+1+lolast-lofirst+1;
109 const Byte_t weven = window & ~1;
110
111 if (weven != window)
112 {
113 *fLog << warn << GetDescriptor()
114 << Form("%s%2i%s%2i",": Total window size has to be even, set last slice from "
115 ,(int)lolast," to ",(int)(lolast-1)) << endl;
116 lolast -= 1;
117 }
118
119 if (window<2)
120 {
121 *fLog << warn << GetDescriptor()
122 << Form("%s%2i%s%2i",": Total window is smaller than 2 FADC sampes, set last slice from"
123 ,(int)lolast," to ",(int)(lofirst+1)) << endl;
124 hilast = hifirst+1;
125 }
126
127
128 MExtractor::SetRange(hifirst,hilast,lofirst,lolast);
129
130
131 fNumHiGainSamples = (Float_t)(fHiGainLast-fHiGainFirst+1);
132 fNumLoGainSamples = (Float_t)(fLoGainLast-fLoGainFirst+1);
133
134 fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples);
135 fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples);
136
137 fNumSamples = fHiGainLast-fHiGainFirst+1+fLoGainLast-fLoGainFirst+1;
138 fSqrtSamples = TMath::Sqrt((Float_t)fNumSamples);
139
140}
141
142// --------------------------------------------------------------------------
143//
144// Calls:
145// - MExtractor::PreProcess
146//
147// The following output containers are also searched and created if
148// they were not found:
149//
150// - MExtractedPINDiode
151//
152Int_t MExtractPINDiode::PreProcess(MParList *pList)
153{
154
155 if (!MExtractor::PreProcess(pList))
156 return kFALSE;
157
158 fPINDiode = (MExtractedSignalPINDiode*)pList->FindCreateObj(AddSerialNumber("MExtractedSignalPINDiode"));
159 if (!fPINDiode)
160 return kFALSE;
161
162 const MPedestalPix &ped = (*fPedestals)[fPINDiodeIdx];
163
164 if (&ped)
165 {
166 fPedestal = ped.GetPedestal();
167 fPedRms = ped.GetPedestalRms();
168 }
169 else
170 {
171 *fLog << err << " Cannot find MPedestalPix of the PIN Diode (idx="
172 << fPINDiodeIdx << ")" << endl;
173 return kFALSE;
174 }
175
176 return kTRUE;
177}
178
179// --------------------------------------------------------------------------
180//
181// The ReInit calls:
182// - MExtractor::ReInit()
183// - fPINDiode->SetUsedFADCSlices(fHiGainFirst, fLoGainLast);
184//
185Bool_t MExtractPINDiode::ReInit(MParList *pList)
186{
187
188 MExtractor::ReInit(pList);
189
190 fPINDiode->SetUsedFADCSlices(fHiGainFirst, fLoGainLast);
191
192 return kTRUE;
193
194}
195
196
197
198void MExtractPINDiode::FindSignalandVarianceHiGain(Byte_t *ptr, Int_t &sum, Int_t &sum2, Byte_t &sat) const
199{
200
201 Byte_t *end = ptr + fHiGainLast - fHiGainFirst + 1;
202
203 while (ptr<end)
204 {
205 sum += *ptr;
206 sum2 += *ptr * *ptr;
207
208 if (*ptr++ >= fSaturationLimit)
209 sat++;
210 }
211}
212
213void MExtractPINDiode::FindSignalandVarianceLoGain(Byte_t *ptr, Int_t &sum, Int_t &sum2, Byte_t &sat) const
214{
215
216 Byte_t *end = ptr + fLoGainLast - fLoGainFirst + 1;
217
218 while (ptr<end)
219 {
220 sum += *ptr;
221 sum2 += *ptr * *ptr;
222
223 if (*ptr++ >= fSaturationLimit)
224 sat++;
225 }
226}
227
228
229
230// --------------------------------------------------------------------------
231//
232// Calculate the integral of the FADC time slices and store them as a new
233// pixel in the MExtractedPINDiode container.
234//
235Int_t MExtractPINDiode::Process()
236{
237
238
239 MRawEvtPixelIter pixel(fRawEvt);
240
241 fPINDiode->Clear();
242
243 pixel.Jump(fPINDiodeIdx);
244
245 Int_t sum = 0;
246 Int_t sum2 = 0;
247 Byte_t sat = 0;
248 Int_t max = 0;
249
250 //
251 // Calculate first the time:
252 //
253 const Int_t maxhi = pixel.GetIdxMaxHiGainSample();
254 const Int_t maxlo = pixel.GetIdxMaxLoGainSample();
255
256 if (maxhi > maxlo)
257 max = maxhi;
258 else
259 max = maxlo + pixel.GetNumHiGainSamples();
260
261 FindSignalandVarianceHiGain(pixel.GetHiGainSamples()+fHiGainFirst,sum,sum2,sat);
262 FindSignalandVarianceLoGain(pixel.GetLoGainSamples()+fLoGainFirst,sum,sum2,sat);
263
264 const Float_t var = ((Float_t)sum2 - (Float_t)sum*sum/fNumSamples)/(fNumSamples-1);
265 const Float_t rms = TMath::Sqrt(var);
266
267 //
268 // FIXME: The following formulae have to be revised!!
269 //
270 fPINDiode->SetExtractedSignal(sum - fPedestal*fNumSamples, fPedRms*fSqrtSamples);
271 fPINDiode->SetExtractedRms (rms, rms/2./fSqrtSamples);
272 fPINDiode->SetExtractedTime (max, rms/fSqrtSamples);
273 fPINDiode->SetSaturation(sat);
274 fPINDiode->SetReadyToSave();
275
276 return kTRUE;
277}
278
Note: See TracBrowser for help on using the repository browser.