source: trunk/MagicSoft/Mars/mtemp/mifae/library/MFDisp.cc@ 5879

Last change on this file since 5879 was 5671, checked in by domingo, 20 years ago
*** empty log message ***
File size: 8.9 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): Eva Domingo, 12/2004 <mailto:domingo@ifae.es>
19! Wolfgang Wittek, 12/2004 <mailto:wittek@mppmu.mpg.de>
20!
21!
22! Copyright: MAGIC Software Development, 2000-2005
23!
24!
25\* ======================================================================== */
26
27//////////////////////////////////////////////////////////////////////////////
28// //
29// MFDisp //
30// //
31// This is a class to set cuts to select an event sample //
32// for the Disp optimization //
33// //
34// //
35//////////////////////////////////////////////////////////////////////////////
36
37#include "MFDisp.h"
38
39#include "MGeomCam.h"
40#include "MHillas.h"
41#include "MHillasSrc.h"
42#include "MNewImagePar.h"
43
44#include "MParList.h"
45
46#include "MLog.h"
47#include "MLogManip.h"
48
49
50ClassImp(MFDisp);
51
52using namespace std;
53
54// --------------------------------------------------------------------------
55//
56// Default constructor.
57//
58MFDisp::MFDisp(const char *name, const char *title)
59{
60 fName = name ? name : "MFDisp";
61 fTitle = title ? title : "Cuts for Disp optimization";
62
63 // default values of cuts
64 SetCuts(0, 1000, 0, 1000, 0, 1000,
65 0.0, 1.e10, 0.0, 1.0, 0.0, 1.0, 0.0, 0.0);
66}
67
68
69// --------------------------------------------------------------------------
70//
71// Set the values for the cuts
72//
73void MFDisp::SetCuts(Int_t islandsmin, Int_t islandsmax,
74 Int_t usedpixelsmin, Int_t usedpixelsmax,
75 Int_t corepixelsmin, Int_t corepixelsmax,
76 Float_t sizemin, Float_t sizemax,
77 Float_t leakage1min, Float_t leakage1max,
78 Float_t leakage2min, Float_t leakage2max,
79 Float_t lengthmin, Float_t widthmin)
80{
81 fIslandsMin = islandsmin;
82 fIslandsMax = islandsmax;
83
84 fUsedPixelsMin = usedpixelsmin;
85 fUsedPixelsMax = usedpixelsmax;
86
87 fCorePixelsMin = corepixelsmin;
88 fCorePixelsMax = corepixelsmax;
89
90 fSizeMin = sizemin;
91 fSizeMax = sizemax;
92
93 fLeakage1Min = leakage1min;
94 fLeakage1Max = leakage1max;
95
96 fLeakage2Min = leakage2min;
97 fLeakage2Max = leakage2max;
98
99 fLengthMin = lengthmin;
100 fWidthMin = widthmin;
101}
102
103
104// --------------------------------------------------------------------------
105//
106Int_t MFDisp::PreProcess(MParList *pList)
107{
108 MGeomCam *cam = (MGeomCam*)pList->FindObject("MGeomCam");
109 if (!cam)
110 {
111 *fLog << err << "MGeomCam (Camera Geometry) not found... aborting." << endl;
112 return kFALSE;
113 }
114
115 fMm2Deg = cam->GetConvMm2Deg();
116
117
118 fHil = (MHillas*)pList->FindObject("MHillas");
119 if (!fHil)
120 {
121 *fLog << err << "MHillas not found... aborting." << endl;
122 return kFALSE;
123 }
124
125 fHilSrc = (MHillasSrc*)pList->FindObject("MHillasSrc");
126 if (!fHilSrc)
127 {
128 *fLog << err << "MHillasSrc not found... aborting." << endl;
129 return kFALSE;
130 }
131
132 fNewImgPar = (MNewImagePar*)pList->FindObject("MNewImagePar");
133 if (!fNewImgPar)
134 {
135 *fLog << err << "MNewImagePar not found... aborting." << endl;
136 return kFALSE;
137 }
138
139
140 memset(fCut, 0, sizeof(fCut));
141
142
143 //--------------------
144 *fLog << inf << "MFDisp cut values :" << endl;
145
146 *fLog << inf << " fIslandsMin, fIslandsMax = "
147 << fIslandsMin << ", " << fIslandsMax << endl;
148
149 *fLog << inf << " fUsedPixelsMin, fUsedPixelsMax = "
150 << fUsedPixelsMin << ", " << fUsedPixelsMax << endl;
151
152 *fLog << inf << " fCorePixelsMin, fCorePixelsMax = "
153 << fCorePixelsMin << ", " << fCorePixelsMax << endl;
154
155 *fLog << inf << " fSizeMin, fSizeMax = "
156 << fSizeMin << ", " << fSizeMax << endl;
157
158 *fLog << inf << " fLeakage1Min, fLeakage1Max = "
159 << fLeakage1Min << ", " << fLeakage1Max << endl;
160
161 *fLog << inf << " fLeakage2Min, fLeakage2Max = "
162 << fLeakage2Min << ", " << fLeakage2Max << endl;
163
164 *fLog << inf << " fLengthMin, fWidthMin = "
165 << fLengthMin << ", " << fWidthMin << endl;
166 //--------------------
167
168 return kTRUE;
169}
170
171
172// --------------------------------------------------------------------------
173//
174// Sum up the number of events passing each cut
175//
176Bool_t MFDisp::Set(Int_t rc)
177{
178 fResult = kTRUE;
179 fCut[rc]++;
180 return kTRUE;
181}
182
183
184// --------------------------------------------------------------------------
185//
186// Evaluate cuts for Disp optimization
187// if selections are fulfilled set fResult = kTRUE;
188//
189Int_t MFDisp::Process()
190{
191 const Double_t length = fHil->GetLength() * fMm2Deg;
192 const Double_t width = fHil->GetWidth() * fMm2Deg;
193 //const Double_t dist = fHilSrc->GetDist()* fMm2Deg;
194 //const Double_t delta = fHil->GetDelta() * kRad2Deg;
195 const Double_t size = fHil->GetSize();
196
197 const Double_t leakage1 = fNewImgPar->GetLeakage1();
198 const Double_t leakage2 = fNewImgPar->GetLeakage2();
199
200 const Int_t numusedpixels = fNewImgPar->GetNumUsedPixels();
201 const Int_t numcorepixels = fNewImgPar->GetNumCorePixels();
202 //const Int_t numislands = fNewImgPar->GetNumIslands();
203
204 fResult = kFALSE;
205
206 //if (numislands<fIslandsMin || numislands>fIslandsMax )
207 // return Set(1);
208
209 if (numusedpixels<fUsedPixelsMin || numusedpixels>fUsedPixelsMax )
210 return Set(2);
211
212 if (numcorepixels<fCorePixelsMin || numcorepixels>fCorePixelsMax )
213 return Set(3);
214
215 if (size<fSizeMin || size>fSizeMax)
216 return Set(4);
217
218 if (leakage1<fLeakage1Min || leakage1>fLeakage1Max)
219 return Set(5);
220
221 if (leakage2<fLeakage2Min || leakage2>fLeakage2Max)
222 return Set(6);
223
224 if (length<fLengthMin || width<fWidthMin)
225 return Set(7);
226
227 fCut[0]++;
228 return kTRUE;
229}
230
231
232// --------------------------------------------------------------------------
233//
234// Prints some statistics about the cuts selections
235//
236Int_t MFDisp::PostProcess()
237{
238 if (GetNumExecutions()==0)
239 return kTRUE;
240
241 *fLog << inf << endl;
242 *fLog << GetDescriptor() << " execution statistics:" << endl;
243
244 *fLog << dec << setfill(' ');
245 *fLog << " " << setw(7) << fCut[1] << " (" << setw(3);
246 *fLog << (int)(fCut[1]*100/GetNumExecutions()+0.5) ;
247 *fLog << "%) Evts skipped due to: No.of islands < " << fIslandsMin ;
248 *fLog << " or > " << fIslandsMax << endl;
249
250 *fLog << " " << setw(7) << fCut[2] << " (" << setw(3);
251 *fLog << (int)(fCut[2]*100/GetNumExecutions()+0.5) ;
252 *fLog << "%) Evts skipped due to: No.of used pixels < " << fUsedPixelsMin ;
253 *fLog << " or > " << fUsedPixelsMax << endl;
254
255 *fLog << " " << setw(7) << fCut[3] << " (" << setw(3);
256 *fLog << (int)(fCut[3]*100/GetNumExecutions()+0.5) ;
257 *fLog << "%) Evts skipped due to: No.of core pixels < " << fCorePixelsMin ;
258 *fLog << " or > " << fCorePixelsMax << endl;
259
260 *fLog << " " << setw(7) << fCut[4] << " (" << setw(3) ;
261 *fLog << (int)(fCut[4]*100/GetNumExecutions()+0.5) ;
262 *fLog << "%) Evts skipped due to: SIZE < " << fSizeMin;
263 *fLog << " or > " << fSizeMax << endl;
264
265 *fLog << " " << setw(7) << fCut[5] << " (" << setw(3) ;
266 *fLog << (int)(fCut[5]*100/GetNumExecutions()+0.5) ;
267 *fLog << "%) Evts skipped due to: Leakage1 < " << fLeakage1Min;
268 *fLog << " or > " << fLeakage1Max << endl;
269
270 *fLog << " " << setw(7) << fCut[6] << " (" << setw(3) ;
271 *fLog << (int)(fCut[6]*100/GetNumExecutions()+0.5) ;
272 *fLog << "%) Evts skipped due to: Leakage2 < " << fLeakage2Min;
273 *fLog << " or > " << fLeakage2Max << endl;
274
275 *fLog << " " << setw(7) << fCut[7] << " (" << setw(3) ;
276 *fLog << (int)(fCut[7]*100/GetNumExecutions()+0.5) ;
277 *fLog << "%) Evts skipped due to: LENGTH <= " << fLengthMin;
278 *fLog << " or WIDTH <= " << fWidthMin << endl;
279
280 *fLog << " " << fCut[0] << " (" ;
281 *fLog << (int)(fCut[0]*100/GetNumExecutions()+0.5) ;
282 *fLog << "%) Evts survived the Disp cuts!" << endl;
283 *fLog << endl;
284
285 return kTRUE;
286}
287
288
289
290
291
Note: See TracBrowser for help on using the repository browser.