source: tags/Mars-V0.9/manalysis/MMcTriggerLvl2.cc

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