source: trunk/Mars/manalysis/MMcTriggerLvl2.cc@ 10083

Last change on this file since 10083 was 8106, checked in by tbretz, 18 years ago
*** empty log message ***
File size: 32.3 KB
Line 
1/* ======================================================================== *\
2! $Name: not supported by cvs2svn $:$Id: MMcTriggerLvl2.cc,v 1.14 2006-10-17 17:15:58 tbretz Exp $
3! --------------------------------------------------------------------------
4!
5! *
6! * This file is part of MARS, the MAGIC Analysis and Reconstruction
7! * Software. It is distributed to you in the hope that it can be a useful
8! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
9! * It is distributed WITHOUT ANY WARRANTY.
10! *
11! * Permission to use, copy, modify and distribute this software and its
12! * documentation for any purpose is hereby granted without fee,
13! * provided that the above copyright notice appear in all copies and
14! * that both that copyright notice and this permission notice appear
15! * in supporting documentation. It is provided "as is" without express
16! * or implied warranty.
17! *
18!
19!
20! Author(s): Antonio Stamerra 1/2003 <mailto:antono.stamerra@pi.infn.it>
21! Author(s): Marcos Lopez 1/2003 <mailto:marcos@gae.ucm.es>
22! Author(s): Nicola Galante 7/2003 <mailto:nicola.galante@pi.infn.it>
23!
24! Copyright: MAGIC Software Development, 2000-2003
25!
26!
27\* ======================================================================== */
28
29
30/////////////////////////////////////////////////////////////////////////////
31// //
32// MMcTriggerLvl2 //
33// Storage container for the 2nd level trigger selection parameters //
34// as part of the 2nd level trigger simulation //
35// //
36// input parameter: //
37// fCompactNN number of next neighboors that define a compact pixel //
38// //
39// //
40// Basic L2T Selection Parameters: //
41// //
42// fLutPseudoSize number of compact pixels in one LUT //
43// fPseudoSize Multiplicity of the bigger cluster //
44// fSizeBiggerCell higher number of fired pixel in cell //
45// //
46/////////////////////////////////////////////////////////////////////////////
47
48#include "MMcTriggerLvl2.h"
49
50#include "MGeomCam.h"
51#include "MGeomPix.h"
52#include "MGeomCamMagic.h"
53
54#include "MMcTrig.hxx"
55
56#include "MMcEvt.hxx"
57
58#include "MLog.h"
59
60#include <TCanvas.h>
61
62ClassImp(MMcTriggerLvl2);
63
64using namespace std;
65
66// ---------------------------
67// Initialization of static const data members pixel_in_cell and pixel_in_lut
68
69//
70// Correspondence TABLE between pixel numbering in the trigger cells and
71// the standard spiral counting
72// (*Note*: Pixels start to count ** from 1 ** instead of 0)
73//
74// This correspondence is valid only for MAGIC-like geometries!
75//
76// FIXME! These definitions should be included in a GeomTrig class
77//
78const Int_t MMcTriggerLvl2::gsPixelsInCell[gsNPixInCell][gsNCells] = {
79 {26, 91, 66, 71, 76, 81, 86, 269, 224, 233, 242, 251, 260, 391, 336, 347, 358, 369, 380},
80 {25, 61, 41, 45, 49, 53, 57, 215, 175, 183, 191, 199, 207, 325, 275, 285, 295, 305, 315},
81 {24, 37, 22, 25, 28, 31, 34, 167, 132, 139, 146, 153, 160, 265, 220, 229, 238, 247, 256},
82 {23, 19, 9, 11, 13, 15, 17, 125, 95, 101, 107, 113, 119, 211, 171, 179, 187, 195, 203},
83 {27, 62, 67, 72, 77, 82, 87, 270, 225, 234, 243, 252, 261, 392, 337, 348, 359, 370, 381},
84 {12, 38, 42, 46, 50, 54, 58, 216, 176, 184, 192, 200, 208, 326, 276, 286, 296, 306, 316},
85 {11, 20, 23, 26, 29, 32, 35, 168, 133, 140, 147, 154, 161, 266, 221, 230, 239, 248, 257},
86 {10, 8, 10, 12, 14, 16, 18, 126, 96, 102, 108, 114, 120, 212, 172, 180, 188, 196, 204},
87 {22, 2, 3, 4, 5, 6, 7, 90, 65, 70, 75, 80, 85, 164, 129, 136, 143, 150, 157},
88 {28, 93, 99, 105, 111, 117, 123, 271, 226, 235, 244, 253, 262, 393, 338, 349, 360, 371, 382},
89 {13, 63, 68, 73, 78, 83, 88, 217, 177, 185, 193, 201, 209, 327, 277, 287, 297, 307, 317},
90 { 4, 39, 43, 47, 51, 55, 59, 169, 134, 141, 148, 155, 162, 267, 222, 231, 240, 249, 258},
91 { 3, 21, 24, 27, 30, 33, 36, 127, 97, 103, 109, 115, 121, 213, 173, 181, 189, 197, 205},
92 { 9, 9, 11, 13, 15, 17, 19, 91, 66, 71, 76, 81, 86, 165, 130, 137, 144, 151, 158},
93 {21, 3, 4, 5, 6, 7, 2, 61, 41, 45, 49, 53, 57, 123, 93, 99, 105, 111, 117},
94 {29, 130, 137, 144, 151, 158, 165, 218, 227, 236, 245, 254, 263, 394, 339, 350, 361, 372, 383},
95 {14, 94, 100, 106, 112, 118, 124, 170, 178, 186, 194, 202, 210, 328, 278, 288, 298, 308, 318},
96 { 5, 64, 69, 74, 79, 84, 89, 128, 135, 142, 149, 156, 163, 268, 223, 232, 241, 250, 259},
97 { 1, 40, 44, 48, 52, 56, 60, 92, 98, 104, 110, 116, 122, 214, 174, 182, 190, 198, 206},
98 { 2, 22, 25, 28, 31, 34, 37, 62, 67, 72, 77, 82, 87, 166, 131, 138, 145, 152, 159},
99 { 8, 10, 12, 14, 16, 18, 8, 38, 42, 46, 50, 54, 58, 124, 94, 100, 106, 112, 118},
100 {20, 11, 13, 15, 17, 19, 9, 20, 23, 26, 29, 32, 35, 88, 63, 68, 73, 78, 83},
101 {30, 131, 138, 145, 152, 159, 166, 219, 228, 237, 246, 255, 264, 329, 279, 289, 299, 309, 319},
102 {15, 95, 101, 107, 113, 119, 125, 171, 179, 187, 195, 203, 211, 269, 224, 233, 242, 251, 260},
103 { 6, 65, 70, 75, 80, 85, 90, 129, 136, 143, 150, 157, 164, 215, 175, 183, 191, 199, 207},
104 { 7, 41, 45, 49, 53, 57, 61, 93, 99, 105, 111, 117, 123, 167, 132, 139, 146, 153, 160},
105 {19, 23, 26, 29, 32, 35, 20, 63, 68, 73, 78, 83, 88, 125, 95, 101, 107, 113, 119},
106 {37, 24, 27, 30, 33, 36, 21, 39, 43, 47, 51, 55, 59, 89, 64, 69, 74, 79, 84},
107 {31, 132, 139, 146, 153, 160, 167, 220, 229, 238, 247, 256, 265, 270, 225, 234, 243, 252, 261},
108 {16, 96, 102, 108, 114, 120, 126, 172, 180, 188, 196, 204, 212, 216, 176, 184, 192, 200, 208},
109 {17, 66, 71, 76, 81, 86, 91, 130, 137, 144, 151, 158, 165, 168, 133, 140, 147, 154, 161},
110 {18, 42, 46, 50, 54, 58, 38, 94, 100, 106, 112, 118, 124, 126, 96, 102, 108, 114, 120},
111 {36, 43, 47, 51, 55, 59, 39, 64, 69, 74, 79, 84, 89, 90, 65, 70, 75, 80, 85},
112 {32, 133, 140, 147, 154, 161, 168, 221, 230, 239, 248, 257, 266, 217, 177, 185, 193, 201, 209},
113 {33, 97, 103, 109, 115, 121, 127, 173, 181, 189, 197, 205, 213, 169, 134, 141, 148, 155, 162},
114 {35, 68, 73, 78, 83, 88, 63, 95, 101, 107, 113, 119, 125, 91, 66, 71, 76, 81, 86}
115 };
116
117//
118// corrispondence between pixels in cell and lut (pix numbering starts from 1)
119//
120const Int_t MMcTriggerLvl2::gsPixelsInLut[gsNLutInCell][gsNPixInLut] = {
121 { 1, 2, 3, 4, 6, 7, 8, 9, 12, 13, 14, 15},
122 {34, 29, 23, 16, 30, 24, 17, 10, 25, 18, 11, 5},
123 {35, 31, 26, 20, 19, 32, 27, 21, 36, 33, 28, 22},
124 };
125
126
127// --------------------------------------------------------------------------
128//
129// Default constructor
130//
131MMcTriggerLvl2::MMcTriggerLvl2(const char *name, const char *title):fCompactNN(2),fTriggerPattern(0)
132{
133 fName = name ? name : ClassName();
134 fTitle = title;
135
136 *fLog << "created MMcTriggerLvl2" << endl;
137
138 //
139 // Initialization of the fPixels array to zero
140 //
141 memset (fPixels,0,gsNPixInCell*gsNCells*sizeof(Int_t));
142 memset (fFiredPixel,0,gsNTrigPixels*sizeof(Int_t));
143
144 // create a new camera
145 SetNewCamera(new MGeomCamMagic);
146
147}
148
149// --------------------------------------------------------------------------
150//
151// Destructor
152//
153MMcTriggerLvl2::~MMcTriggerLvl2()
154{
155 delete fGeomCam;
156}
157
158// --------------------------------------------------------------------------
159//
160// Print functions for debugging and testing
161//
162// Options available:
163//
164// "cell":
165// Print out the pixel number of a given trigger cell
166//
167// "status":
168// Print a table with the status (pixel is fired or not after Lvl1) for
169// all the pixels in all the trigger cells
170//
171// "check"
172// Print a warning message when no starting pixel is found in the
173// CalcPseudoSize method.
174//
175// no option:
176// Print the value of the selection parameters
177//
178//
179void MMcTriggerLvl2::Print(Option_t *opt) const
180{
181 TString str(opt);
182
183 if (str.Contains("status", TString::kIgnoreCase))
184 {
185 *fLog << " Status of cells 1-9" <<endl;
186 for(int i=0; i<gsNPixInCell; i++)
187 {
188 for(int j=0; j<9; j++)
189 {
190 // *fLog.width(3);
191 *fLog <<gsPixelsInCell[i][j]-1 << ":" << fPixels[i][j] << "\t ";
192 }
193 *fLog << endl;
194 }
195 }
196 else if (str.Contains("cell", TString::kIgnoreCase))
197 {
198 *fLog << " Pixel numbering in cells 1-9" <<endl;
199 for (int i=0;i<gsNPixInCell;i++)
200 {
201 for(int j=0; j<9; j++)
202 {
203 *fLog << gsPixelsInCell[i][j]-1 << "\t";
204 }
205 *fLog << endl;
206 }
207 }
208 else if (str.Contains("check", TString::kIgnoreCase))
209 {
210 // check when no starting pixels was found (<<--to be fixed)
211 if (fPseudoSize < 0)
212 *fLog << " Warning: starting pixel not found. Pseudosize set to -1." << endl;
213 }
214 else
215 {
216 *fLog << " L2T selection parameters: " << endl;
217 *fLog << " - LutPseudoSize = " << fLutPseudoSize << endl;
218 *fLog << " - CellPseudoSize = " << fCellPseudoSize << endl;
219 *fLog << " - PseudoSize = " << fPseudoSize << endl;
220 *fLog << " - BiggerCellSize = " << fSizeBiggerCell << endl;
221 *fLog << " - TriggerPattern = " << fTriggerPattern << endl;
222 }
223
224}
225
226// --------------------------------------------------------------------------
227//
228// Take the information supplied by the Lvl1 (it reads from the object MMcTrig
229// the status of all pixels after Lvl1) and pass it to the Lvl2 (fill the
230// array fPixels)
231//
232//
233void MMcTriggerLvl2::SetLv1(MMcTrig *trig)
234{
235 if (!trig)
236 return;
237
238 fMcTrig = trig;
239
240 for(int i=0; i<gsNPixInCell; i++)
241 {
242 for(int j=0; j<gsNCells; j++)
243 {
244 int pixel = gsPixelsInCell[i][j]-1;
245 fPixels[i][j] = (fMcTrig->IsPixelFired(pixel,0)) ? 1 : 0;
246 fFiredPixel[pixel]=(fMcTrig->IsPixelFired(pixel,0)) ? 1 : 0;
247 //if (fFiredPixel[pixel]==1)
248 //*fLog << pixel<<",";
249 }
250 }
251 //*fLog<<endl<<"Fine evento"<<endl;
252}
253
254
255// --------------------------------------------------------------------------
256//
257// Set the trigger status ( pixel fired(=1) or not(=0) ) manually for a given
258// pixel of a given cell
259//
260//
261void MMcTriggerLvl2::SetPixelFired(Int_t pixel, Int_t fired)
262{
263 for(int i=0; i<gsNPixInCell; i++)
264 {
265 for(int j=0; j<gsNCells; j++)
266 {
267 if(gsPixelsInCell[i][j]-1==pixel) fPixels[i][j]=fired;
268 }
269 }
270
271}
272
273// --------------------------------------------------------------------------
274//
275// Calculate the Level 2 Trigger (L2T) parameters. They can be used
276// to select the events that have passed the L2T selection rule.
277//
278//
279void MMcTriggerLvl2::Calc()
280{
281
282 // Find the Lut and cell with the higher LutPseudoSize.
283 int lutcell = CalcBiggerLutPseudoSize();
284 fMaxCell = lutcell/10;
285 int maxlut = lutcell-fMaxCell*10;
286 fLutPseudoSize = GetLutCompactPixel(fMaxCell,maxlut);
287
288 fMaxCell = CalcBiggerCellPseudoSize(); // fCellPseudoSize
289 // fMaxCell is used by the PseudoSize to find the starting pixel
290 CalcPseudoSize();
291
292 fSizeBiggerCell = GetCellNumberFired(CalcBiggerFiredCell());
293
294 //*fLog << "fLPS="<<fLutPseudoSize<<endl;
295}
296
297
298
299// --------------------------------------------------------------------------
300//
301// For a given cell, just count how many pixels have been fired after Lvl1
302//
303Int_t MMcTriggerLvl2::GetCellNumberFired(int cell)
304{
305 int size=0;
306
307 for(int i=0; i<gsNPixInCell; i++)
308 {
309 size += fPixels[i][cell];
310 }
311
312 return size;
313
314}
315
316
317// --------------------------------------------------------------------------
318//
319// Find the cell which the bigger number of fired pixels
320//
321//
322Int_t MMcTriggerLvl2::CalcBiggerFiredCell()
323{
324 int size=-1;
325 int cell=-1;
326
327 for(int j=0; j<gsNCells; j++)
328 {
329 if (GetCellNumberFired(j) > size)
330 {
331 size = GetCellNumberFired(j);
332 cell = j;
333 }
334 }
335
336 // *fLog << "size " <<size <<" in cell " << cell << endl;
337
338 return cell;
339
340}
341
342// --------------------------------------------------------------------------
343//
344// Calculate the higher LutPseudoSize of one event defined as the higher number
345// of compact pixel in the LUTs of the trigger cells.
346// neighpix is the number of Next-Neighbors which defines the compact pixel
347// topology (it can be 2-NN or 3-NN)
348//
349// It returns the cell and the LUT with the bigger LutPseudoSize, coded
350// accordingly to: cell*10 + LUT
351//
352Int_t MMcTriggerLvl2::CalcBiggerLutPseudoSize()
353{
354 int size=0;
355 int cell=-1;
356 int lut=-1;
357
358 for(int j=0; j<gsNCells; j++)
359 {
360 for(int i=0; i<gsNLutInCell; i++)
361 {
362 if (GetLutCompactPixel(j,i) >= size)
363 {
364 size = GetLutCompactPixel(j,i);
365 cell = j;
366 lut = i;
367 }
368 }
369 }
370
371 //*fLog <<"Max cell: " << cell+1 << " Max Lut: " << lut+1 << " PseudoSize: " << size <<endl;
372
373 return cell*10+lut;
374}
375
376
377// --------------------------------------------------------------------------
378//
379// Identify compact pixels in one LUT and calculate the LutPseudoSize
380// (=number of compact pixels in one LUT)
381// neighpix is the number of Next-Neighbors which defines the compact pixel
382// topology.
383// Up to now only 2NN or 3NN are considered
384// <!if changed: change check made by Preprocess in MMcTriggerLvl2Calc>
385//
386// Returns:
387// -1 wrong neighpix
388// -2 wrong cell (<1 or >19)
389// -3 wrong lut (<1 or >3)
390// else:
391// the number of compact pixels in one LUT.
392//
393Int_t MMcTriggerLvl2::GetLutCompactPixel(int cell, int lut)
394{
395 int size=0;
396 int lutpix, a[gsNPixInLut];
397 int neighpix= (*this).fCompactNN;
398
399 // check on input variables
400
401 if (neighpix >3 || neighpix < 2)
402 return(-1);
403
404 if (cell<0 || cell> gsNCells-1)
405 return(-2);
406
407 if (lut<0 || lut> gsNLutInCell-1)
408 return(-3);
409
410
411 //
412 // Warning!!! Following configurations are valid only for the standard MAGIC
413 // trigger geometry; FIXME! these should be coded somewhere else....
414 //
415
416 // LUT 1 and 2 are similar; LUT 3 differs
417
418 for(int j=0; j< gsNPixInLut; j++)
419 {
420 lutpix = gsPixelsInLut[lut][j]-1;
421 // *fLog << "j=" <<j<<" lutpix="<<lutpix<<" Cell="<<cell<<endl;
422 a[j] = fPixels[lutpix][cell];
423 }
424
425 //
426 // Look for compact pixels 2NN
427 //
428 if (neighpix==2)
429 {
430 // case Lut 1 and 2
431 if (lut == 0 || lut == 1)
432 {
433 if (a[0] && a[1] && a[4])
434 size ++;
435 if (a[1] && ((a[0] && a[4]) ||
436 (a[4] && a[5]) ||
437 (a[5] && a[2]) ))
438 size ++;
439 if (a[2] && ((a[1] && a[5]) ||
440 (a[5] && a[6]) ||
441 (a[6] && a[3]) ))
442 size ++;
443 if (a[3] && ((a[2] && a[6]) ||
444 (a[6] && a[7]) ))
445 size ++;
446 if (a[4] && ((a[0] && a[1]) ||
447 (a[1] && a[5]) ||
448 (a[5] && a[8]) ))
449 size ++;
450 if (a[5] && ((a[1] && a[2]) ||
451 (a[2] && a[6]) ||
452 (a[6] && a[9]) ||
453 (a[9] && a[8]) ||
454 (a[8] && a[4]) ||
455 (a[4] && a[1]) ))
456 size ++;
457 if (a[6] && ((a[2] && a[3]) ||
458 (a[3] && a[7]) ||
459 (a[7] && a[10]) ||
460 (a[10] && a[9]) ||
461 (a[9] && a[5]) ||
462 (a[5] && a[2]) ))
463 size ++;
464 if (a[7] && ((a[3] && a[6]) ||
465 (a[6] && a[10]) ||
466 (a[10] && a[11]) ))
467 size ++;
468 if (a[8] && ((a[4] && a[5]) ||
469 (a[5] && a[9]) ))
470 size ++;
471 if (a[9] && ((a[8] && a[5]) ||
472 (a[5] && a[6]) ||
473 (a[6] && a[10]) ))
474 size ++;
475 if (a[10] && ((a[9] && a[6]) ||
476 (a[6] && a[7]) ||
477 (a[7] && a[11]) ))
478 size ++;
479 if (a[11] && a[7] && a[10])
480 size ++;
481 }
482
483 // case Lut 3
484 if (lut==2)
485 {
486 if (a[0] && a[1] && a[5])
487 size ++;
488 if (a[1] && ((a[0] && a[5]) ||
489 (a[5] && a[2]) ))
490 size ++;
491 if (a[2] && ((a[1] && a[5]) ||
492 (a[5] && a[6]) ||
493 (a[3] && a[4]) ||
494 (a[6] && a[3]) ))
495 size ++;
496 if (a[3] && ((a[2] && a[6]) ||
497 (a[6] && a[7]) ||
498 (a[2] && a[4]) ))
499 size ++;
500 if (a[4] && ((a[2] && a[3]) ))
501 size ++;
502 if (a[5] && ((a[1] && a[2]) ||
503 (a[2] && a[6]) ||
504 (a[6] && a[9]) ||
505 (a[9] && a[8]) ))
506 size ++;
507 if (a[6] && ((a[2] && a[3]) ||
508 (a[3] && a[7]) ||
509 (a[7] && a[10]) ||
510 (a[10] && a[9]) ||
511 (a[9] && a[5]) ||
512 (a[5] && a[2]) ))
513 size ++;
514 if (a[7] && ((a[3] && a[6]) ||
515 (a[6] && a[10]) ||
516 (a[10] && a[11]) ))
517 size ++;
518 if (a[8] && a[5] && a[9])
519 size ++;
520 if (a[9] && ((a[8] && a[5]) ||
521 (a[5] && a[6]) ||
522 (a[6] && a[10]) ))
523 size ++;
524 if (a[10] && ((a[9] && a[6]) ||
525 (a[6] && a[7]) ||
526 (a[7] && a[11]) ))
527 size ++;
528 if (a[11] && a[7] && a[10])
529 size ++;
530 }
531 }
532
533 //
534 // Look for compact pixels 3NN
535 //
536 if (neighpix==3)
537 {
538 // case Lut 1 and 2
539 if (lut == 0 || lut == 1)
540 {
541 if (a[0] && a[1] && a[4] && a[5]) // pix 0 is compact if there is a 4NN
542 size ++;
543 if (a[1] && ((a[0] && a[4] && a[5]) ||
544 (a[2] && a[4] && a[5]) ))
545 size ++;
546 if (a[2] && ((a[1] && a[5] && a[6]) ||
547 (a[5] && a[6] && a[3]) ))
548 size ++;
549 if (a[3] && (a[2] && a[6] && a[7] ))
550 size ++;
551 if (a[4] && ((a[0] && a[1] && a[5]) ||
552 (a[1] && a[5] && a[8]) ))
553 size ++;
554 if (a[5] && ((a[1] && a[2] && a[6]) ||
555 (a[2] && a[6] && a[9]) ||
556 (a[6] && a[9] && a[8]) ||
557 (a[9] && a[8] && a[4]) ||
558 (a[8] && a[4] && a[1]) ||
559 (a[4] && a[1] && a[2]) ))
560 size ++;
561 if (a[6] && ((a[2] && a[3] && a[7]) ||
562 (a[3] && a[7] && a[10]) ||
563 (a[7] && a[10] && a[9]) ||
564 (a[10] && a[9] && a[5]) ||
565 (a[9] && a[5] && a[2]) ||
566 (a[5] && a[2] && a[3]) ))
567 size ++;
568 if (a[7] && ((a[3] && a[6] && a[10]) ||
569 (a[6] && a[10] && a[11]) ))
570 size ++;
571 if (a[8] && (a[4] && a[5] && a[9] ))
572 size ++;
573 if (a[9] && ((a[8] && a[5] && a[6]) ||
574 (a[5] && a[6] && a[10]) ))
575 size ++;
576 if (a[10] && ((a[9] && a[6] && a[7]) ||
577 (a[6] && a[7] && a[11]) ))
578 size ++;
579 if (a[11] && a[7] && a[10] && a[6]) //pix 0 is compact if there is a 4NN
580 size ++;
581 }
582
583 // case Lut 3
584 if (lut==2)
585 {
586 if (a[0] && a[1] && a[5] && a[8]) // pix 0 is compact if there is a 4NN
587 size ++;
588 if (a[1] && (a[0] && a[5] && a[2])) //pix 0 is compact if there is a 4NN
589 size ++;
590 if (a[2] && ((a[1] && a[5] && a[6]) ||
591 (a[3] && a[5] && a[6]) ||
592 (a[6] && a[3] && a[4]) ))
593 size ++;
594 if (a[3] && ((a[2] && a[4] && a[6]) ||
595 (a[2] && a[6] && a[7]) ))
596 size ++;
597 if (a[4] && (a[2] && a[3] && a[6] ))
598 size ++;
599 if (a[5] && ((a[1] && a[2] && a[6]) ||
600 (a[2] && a[6] && a[9]) ||
601 (a[6] && a[9] && a[8]) ))
602 size ++;
603 if (a[6] && ((a[2] && a[3] && a[7]) ||
604 (a[3] && a[7] && a[10]) ||
605 (a[7] && a[10] && a[9]) ||
606 (a[10] && a[9] && a[5]) ||
607 (a[9] && a[5] && a[2]) ||
608 (a[5] && a[2] && a[3]) ))
609 size ++;
610 if (a[7] && ((a[3] && a[6] && a[10]) ||
611 (a[6] && a[10] && a[11]) ))
612 size ++;
613 if (a[8] && a[5] && a[9] && a[6]) //pix 0 is compact if there is a 4NN
614 size ++;
615 if (a[9] && ((a[8] && a[5] && a[6]) ||
616 (a[5] && a[6] && a[10]) ))
617 size ++;
618 if (a[10] && ((a[9] && a[6] && a[7]) ||
619 (a[6] && a[7] && a[11]) ))
620 size ++;
621 if (a[11] && a[7] && a[10] && a[6]) //pix 0 is compact if there is a 4NN
622 size ++;
623 }
624 }
625
626
627 if(size<0 ||size>gsNPixInLut)
628 *fLog << "*" << size <<"-";
629
630 return size;
631}
632
633
634// --------------------------------------------------------------------------
635//
636// Look for clusters and calculate the PseudoSize of one event,
637// defined as the multiplicity of the bigger cluster.
638//
639//
640//
641void MMcTriggerLvl2::CalcPseudoSize()
642{
643 // Fill the fCompactPixel array with the compact pixels
644 CalcCompactPixels(fGeomCam);
645
646 // seek the LUT with higher number of compact pixels
647 //
648 //int fMaxCell = CalcBiggerCellPseudoSize();
649 int sizetemp=0;
650 int maxlut=0;
651
652 for (int i=0;i<gsNLutInCell;i++)
653 if (GetLutCompactPixel(fMaxCell,i) > sizetemp)
654 {
655 maxlut=i;
656 sizetemp = GetLutCompactPixel(fMaxCell,i);
657 }
658
659 int startpix;
660 //
661 // seek a starting pixel for the iteration inside the lut
662 //
663 int check=1;
664 for (int pixlut=0;pixlut<gsNPixInLut;pixlut++)
665 {
666 int pixcell =gsPixelsInLut[maxlut][pixlut]-1;
667 startpix = gsPixelsInCell[pixcell][fMaxCell]-1;
668 //*fLog << "pix, compact:" << startpix << "@"<<fCompactPixel[startpix];
669 if (fCompactPixel[startpix]) // a starting pixel was found
670 break;
671 check++;
672 }
673
674 //*fLog << "check = " << check << endl;
675 // A LUT contains 12 pixels
676 if (check > gsNPixInLut)
677 {
678 *fLog <<"Warning: a starting pixel was not found! - PseudoSize = "<< fPseudoSize << endl;
679 fPseudoSize=-1;
680 return;
681 }
682
683 //
684 // Bulding cluster
685 //
686 Int_t cluster[gsNTrigPixels];
687 int pnt=0;
688 int pnt2=0; //pointer in the array fCluster_pix, needed for filling
689
690 memset (cluster,0,gsNTrigPixels*sizeof(Int_t));
691 memset (fCluster_pix,-1,gsNTrigPixels*sizeof(Int_t));
692
693 cluster[startpix]=1;
694 fCluster_pix[0]=startpix; //the starting pix is the first in cluster
695
696 // Look at neighbour pixs if they are compact (iterative search)
697 // until the array (fCluster_pix) has no more compact pixels.
698 // pnt points to the pixel in the array fCluster_pix whose neighbors are
699 // under study; pnt2 points to the last element of this array.
700 //
701 while (fCluster_pix[pnt] != -1)
702 {
703 const MGeomPix &pix=(*fGeomCam)[fCluster_pix[pnt]];
704
705 for (int i=0;i<pix.GetNumNeighbors();i++)
706 {
707 int pix_neigh = pix.GetNeighbor(i);
708 // check if pixel is fired and doesn't belong to cluster
709 if (fCompactPixel[pix_neigh] && !cluster[pix_neigh])
710 {
711 cluster[pix_neigh] = 1;
712 fCluster_pix[++pnt2] = pix_neigh;
713 }
714 }
715 pnt++;
716 }
717
718 fPseudoSize = pnt;
719
720 //if (fPseudoSize < 4)
721 // *fLog << "fPseudoSize = " << fPseudoSize << endl;
722
723 // *fLog << "ClusterID:" <<(*clust).GetClusterId() << " Mult:" << (*clust).GetMultiplicity()<<endl;
724 /*
725 *fLog <<"PSize: "<< fPseudoSize << " in cell:" << fMaxCell << " lut:" <<maxlut << endl << " Pixels: ";
726 for (int i =0;i<fPseudoSize; i++)
727 *fLog << fCluster_pix[i]+1 <<"; ";
728 *fLog << endl;
729 */
730 return;
731}
732
733// --------------------------------------------------------------------------
734//
735// Fill the fCompactPixels array with the pixels which are compact
736//
737// neighpix is the number of Next-Neighbors which defines the compact pixel
738// topology (it can be 2-NN or 3-NN)
739//
740// Note: it is a *global* method; it looks in the all camera as a whole
741//
742//
743void MMcTriggerLvl2::CalcCompactPixels(MGeomCam *geom)
744{
745 memset (fCompactPixel,0,gsNTrigPixels*sizeof(Int_t));
746 // *fLog << endl << "NEW Event!";
747 for(Int_t pixid=0; pixid<gsNTrigPixels; pixid++)
748 {
749 // Look if the pixel is fired, otherwise continue
750 if (!fFiredPixel[pixid])
751 continue;
752
753 const MGeomPix &pix=(*geom)[pixid];
754
755 // Reshuffle pixneighbour order, to arrange them (anti)clockwise
756 // around the current pixel (Bubble sorting)
757 // The NeighPixOrdered array has a doubledd size so that
758 // the content in the first index is repeated in the last one
759 // to have a closed loop
760
761 Int_t NeighPixOrdered[2*pix.GetNumNeighbors()];
762
763 for (Int_t j=0; j<pix.GetNumNeighbors(); j++)
764 NeighPixOrdered[j] = pix.GetNeighbor(j);
765
766 for (Int_t j=0; j<pix.GetNumNeighbors()-1; j++)
767 {
768 for (Int_t jk=pix.GetNumNeighbors()-1;jk>j;--jk)
769 {
770 UInt_t tmp = NeighPixOrdered[j+1];
771 const MGeomPix &pixneigh = (*geom)[NeighPixOrdered[j]];
772 for (int k=0; k<pix.GetNumNeighbors(); k++)
773 if (NeighPixOrdered[jk] == pixneigh.GetNeighbor(k))
774 {
775 NeighPixOrdered[j+1] = NeighPixOrdered[jk];
776 NeighPixOrdered[jk] = tmp;
777 }
778 }
779 }
780
781 // Duplicate second half to take into account configurations
782 // containing tha last and first pixel
783 for (Int_t j=pix.GetNumNeighbors(); j<2*pix.GetNumNeighbors(); j++)
784 NeighPixOrdered[j] = NeighPixOrdered[j-pix.GetNumNeighbors()];
785
786 // Look for compact pixels.
787 // A compact pixel must have at least fCompactNN adjacent fired neighbors
788 // It checks the 6 different configurations of neighbors pixels; if
789 // one has fCompacNN adjacent fired pixels than the pixel pixid is
790 // promoted to compact pixel and saved into the fCompactPixel array.
791 //
792 for (int i=0;i<pix.GetNumNeighbors();i++)
793 {
794 int j=0; // j counts the adjacent fired pixels
795 //*fLog << pixid <<"->"<< pix.GetNeighbor(i+j) << endl;
796 while ((i+j < 2*pix.GetNumNeighbors()) && (fFiredPixel[NeighPixOrdered[i+j]]==1) && (j < fCompactNN))
797 j++;
798 if (j>=fCompactNN) // configuration satisfies the compact condition
799 {
800 fCompactPixel[pixid]=1; // pixel is compact
801 // *fLog << ","<<pixid;
802 break;
803 }
804 }
805 }
806
807}
808
809// --------------------------------------------------------------------------
810//
811// The Energy has to be given by this class to the Energy-PSSize correlation
812// histogram (MHMcTriggerLvl2)
813//
814void MMcTriggerLvl2::GetEnergy(MMcEvt *fMcEvt)
815{
816 const MMcEvt &h = *(MMcEvt *)fMcEvt;
817 fEnergy = h.GetEnergy();
818 return;
819}
820
821
822// --------------------------------------------------------------------------
823//
824// Looks for a x-NN compact pattern in the whole camera
825// We have x-NN compact pattern when a triggered pix has
826// x-1 triggered neighbor pixels.
827// The variable fTriggerPattern = x is computed
828// (x= 2,3,4,5,6,7)
829//
830// x=2 * *
831//
832// x=3 * *
833// * *
834// x=4 * *
835// *
836// * *
837// x=5 * * *
838// *
839// x=6 * * *
840// * *
841// * *
842// x=7 * * *
843// * *
844//
845void MMcTriggerLvl2::CalcTriggerPattern(MGeomCam *geom)
846{
847 fTriggerPattern=0; //initialize
848
849 for(Int_t pixid=0;pixid<gsNTrigPixels;pixid++)
850 {
851 // Look if the pixel is fired, otherwise continue
852 if (!fFiredPixel[pixid])
853 continue;
854
855 const MGeomPix &pix=(*geom)[pixid];
856
857 // Look for x-NN compact pattern
858 // If a x-NN pattern exists then a pixel must have
859 // at least x-1 adjacent neighbors (look at patterns)
860 // For each triggered pixel the number of adjacent triggered pixels
861 // is counted.
862 //
863 int j=1;
864 for (int i=0;i<pix.GetNumNeighbors();i++)
865 if (fFiredPixel[pix.GetNeighbor(i)]==1) j++;
866
867 if (j > fTriggerPattern)
868 fTriggerPattern=j;
869
870 if (fTriggerPattern==7)
871 break; // the 7-NN (max) pattern was found: exit
872
873 } // next pixel
874}
875
876
877// --------------------------------------------------------------------------
878// Look for the cell with higher number of Compact pixels
879// Return the cell number (starting from 0)
880//
881Int_t MMcTriggerLvl2::CalcBiggerCellPseudoSize()
882{
883 Int_t maxcell=-1;
884
885 fCellPseudoSize=0;
886
887 for (Int_t i=0;i<gsNCells;i++)
888 {
889 int size = GetCellCompactPixel(i,fGeomCam);
890 if (size > fCellPseudoSize)
891 {
892 fCellPseudoSize = size;
893 maxcell = i;
894 }
895 }
896
897 //*fLog << "fCellPseudoSize = " << fCellPseudoSize << " in cell N. " << fMaxCell+1 << endl;
898
899 return maxcell;
900}
901
902// --------------------------------------------------------------------------
903// Compute the number of compact pixels in one cell
904//
905Int_t MMcTriggerLvl2::GetCellCompactPixel(int cell, MGeomCam *geom)
906{
907 int size=0;
908
909 // check on input variables
910
911 if (cell<0 || cell>gsNCells-1)
912 return(-2);
913
914 //*fLog << " CNN:" << fCompactNN;
915 //*fLog << "Cell: " << cell+1 << " Compat Pixels:";
916
917 for(Int_t id=0; id<gsNPixInCell; id++)
918 {
919 UInt_t pixid = gsPixelsInCell[id][cell]-1;
920
921 // Look if the pixel is fired, otherwise continue
922 if (!fFiredPixel[pixid])
923 continue;
924
925 //*fLog << "Fired pix:"<<pixid+1 << " ";
926
927 const MGeomPix &pix=(*geom)[pixid];
928
929 // Reshuffle pixneighbour order, to arrange them (anti)clockwise
930 // around the current pixel (Bubble sorting)
931 // The NeighPixOrdered has one index more so that
932 // the content in the first index is repeated in the last one
933 // to have a closed loop
934 Int_t NeighPixOrdered[2*pix.GetNumNeighbors()];
935 for (Int_t j=0; j<pix.GetNumNeighbors(); j++)
936 NeighPixOrdered[j] = pix.GetNeighbor(j);
937
938 for (Int_t j=0; j<pix.GetNumNeighbors()-1; j++)
939 {
940 for (Int_t jk=pix.GetNumNeighbors()-1;jk>j;--jk)
941 {
942 UInt_t tmp = NeighPixOrdered[j+1];
943 const MGeomPix &pixneigh = (*geom)[NeighPixOrdered[j]];
944 for (int k=0; k<pix.GetNumNeighbors(); k++)
945 if (NeighPixOrdered[jk] == pixneigh.GetNeighbor(k))
946 {
947 NeighPixOrdered[j+1] = NeighPixOrdered[jk];
948 NeighPixOrdered[jk] = tmp;
949 }
950 }
951 }
952
953 // Duplicate second half of the array to take into account configurations
954 // containing tha last and first pixel
955 for (Int_t j=pix.GetNumNeighbors(); j<2*pix.GetNumNeighbors(); j++)
956 NeighPixOrdered[j] = NeighPixOrdered[j-pix.GetNumNeighbors()];
957
958 // Look for compact pixels.
959
960 // A compact pixel must have at least fCompactNN adjacent neighbors
961 // It checks the 6 different configurations of neighbors pixels.
962 // The neighbour pixels must belong to the cell
963
964 // *fLog << "cell:"<< cell << " ordered pixels:";
965 /*
966 for (int i=0;i<2*pix.GetNumNeighbors();i++)
967 {
968 if (fFiredPixel[NeighPixOrdered[i]])
969 *fLog << NeighPixOrdered[i]+1 << "*;";
970 else
971 *fLog << NeighPixOrdered[i]+1 << ";";
972 }
973 */
974 //*fLog <<endl;
975 //*fLog << pixid <<"->"<< pix.GetNumNeighbors() << " CNN="<< fCompactNN <<endl;
976
977
978 for (int i=0;i<pix.GetNumNeighbors();i++)
979 {
980 int j=0;
981 while ((i+j < 2*pix.GetNumNeighbors()) && (fFiredPixel[NeighPixOrdered[i+j]]==1) && (j < fCompactNN) && IsPixelInCell(NeighPixOrdered[i+j],cell) )
982 j++;
983
984 if (j>=fCompactNN) //configuration satisfies the compact condition
985 {
986 size++; // pixel is compact
987 //*fLog << "->" << pixid+1;
988 break; // (new pixel)
989 }
990 }
991 }
992
993 //*fLog <<" - size:" << size << endl<<endl;
994
995 return size;
996
997}
998
999//---------------------------------------------------------------------
1000// Check if a given pixel belongs to a given cell
1001//
1002Bool_t MMcTriggerLvl2::IsPixelInCell(Int_t pixel, Int_t cell)
1003{
1004 for (int i=0; i<gsNPixInCell; i++)
1005 if ((gsPixelsInCell[i][cell]-1) == pixel)
1006 return kTRUE;
1007
1008 return kFALSE;
1009}
1010
1011//---------------------------------------------------------------------
1012// Check if a given pixel is in the trigger region
1013//
1014Bool_t MMcTriggerLvl2::IsPixelInTrigger(Int_t pixel) const
1015{
1016 for (int cell=0; cell<gsNCells; cell++)
1017 for (int i=0; i<gsNPixInCell; i++)
1018 if ((gsPixelsInCell[i][cell]-1) == pixel)
1019 return kTRUE;
1020
1021 return kFALSE;
1022}
1023
1024// --------------------------------------------------------------------------
1025//
1026// Returns, depending on the type flag:
1027//
1028// 0: 1,0 if the pixel is triggered (1) or not (0)
1029//
1030//
1031Bool_t MMcTriggerLvl2::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
1032{
1033 // Pixel in no-trigger region are set to 0
1034 val = this->IsPixelInTrigger(idx) ? fFiredPixel[idx]+fCompactPixel[idx] : 0;
1035
1036 return kTRUE;
1037}
1038
1039void MMcTriggerLvl2::DrawPixelContent(Int_t num) const
1040{
1041 *fLog << "MMcTriggerLvl2::DrawPixelContent - not available." << endl;
1042}
Note: See TracBrowser for help on using the repository browser.