source: trunk/MagicSoft/Mars/manalysis/MMcTriggerLvl2.cc@ 2216

Last change on this file since 2216 was 2210, checked in by tbretz, 21 years ago
*** empty log message ***
File size: 24.0 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!
21! Copyright: MAGIC Software Development, 2000-2003
22!
23!
24\* ======================================================================== */
25
26
27/////////////////////////////////////////////////////////////////////////////
28// //
29// MMcTriggerLvl2 //
30// Storage container for the 2nd level trigger selection parameters //
31// as part of the 2nd level trigger simulation //
32// //
33// input parameter: //
34// fCompactNN number of next neighboors that define a compact pixel //
35// //
36// //
37// Basic L2T Selection Parameters: //
38// //
39// fLutPseudoSize number of compact pixels in one LUT //
40// fPseudoSize Multiplicity of the bigger cluster //
41// fSizeBiggerCell higher number of fired pixel in cell //
42// //
43/////////////////////////////////////////////////////////////////////////////
44
45#include "MMcTriggerLvl2.h"
46
47#include "MGeomCam.h"
48#include "MGeomPix.h"
49#include "MGeomCamMagic.h"
50#include "MHCamera.h"
51
52#include "MMcTrig.hxx"
53
54#include "MLog.h"
55
56#include <TCanvas.h>
57
58ClassImp(MMcTriggerLvl2);
59
60using namespace std;
61
62// ---------------------------
63// Initialization of static const data members pixel_in_cell and pixel_in_lut
64
65//
66// Correspondence TABLE between pixel numbering in the trigger cells and
67// the standard spiral counting
68// (Note: Pixels start to count from 1 instead of 0)
69//
70const Int_t MMcTriggerLvl2::gsPixelsInCell[36][19] = {
71 {26, 91, 66, 71, 76, 81, 86, 269, 224, 233, 242, 251, 260, 391, 336, 347, 358, 369, 380},
72 {25, 61, 41, 45, 49, 53, 57, 215, 175, 183, 191, 199, 207, 325, 275, 285, 295, 305, 315},
73 {24, 37, 22, 25, 28, 31, 34, 167, 132, 139, 146, 153, 160, 265, 220, 229, 238, 247, 256},
74 {23, 19, 9, 11, 13, 15, 17, 125, 95, 101, 107, 113, 119, 211, 171, 179, 187, 195, 203},
75 {27, 62, 67, 72, 77, 82, 87, 270, 225, 234, 243, 252, 261, 392, 337, 348, 359, 370, 381},
76 {12, 38, 42, 46, 50, 54, 58, 216, 176, 184, 192, 200, 208, 326, 276, 286, 296, 306, 316},
77 {11, 20, 23, 26, 29, 32, 35, 168, 133, 140, 147, 154, 161, 266, 221, 230, 239, 248, 257},
78 {10, 8, 10, 12, 14, 16, 18, 126, 96, 102, 108, 114, 120, 212, 172, 180, 188, 196, 204},
79 {22, 2, 3, 4, 5, 6, 7, 90, 65, 70, 75, 80, 85, 164, 129, 136, 143, 150, 157},
80 {28, 93, 99, 105, 111, 117, 123, 271, 226, 235, 244, 253, 262, 393, 338, 349, 360, 371, 382},
81 {13, 63, 68, 73, 78, 83, 88, 217, 177, 185, 193, 201, 209, 327, 277, 287, 297, 307, 317},
82 { 4, 39, 43, 47, 51, 55, 59, 169, 134, 141, 148, 155, 162, 267, 222, 231, 240, 249, 258},
83 { 3, 21, 24, 27, 30, 33, 36, 127, 97, 103, 109, 115, 121, 213, 173, 181, 189, 197, 205},
84 { 9, 9, 11, 13, 15, 17, 19, 91, 66, 71, 76, 81, 86, 165, 130, 137, 144, 151, 158},
85 {21, 3, 4, 5, 6, 7, 2, 61, 41, 45, 49, 53, 57, 123, 93, 99, 105, 111, 117},
86 {29, 130, 137, 144, 151, 158, 165, 218, 227, 236, 245, 254, 263, 394, 339, 350, 361, 372, 383},
87 {14, 94, 100, 106, 112, 118, 124, 170, 178, 186, 194, 202, 210, 328, 278, 288, 298, 308, 318},
88 { 5, 64, 69, 74, 79, 84, 89, 128, 135, 142, 149, 156, 163, 268, 223, 232, 241, 250, 259},
89 { 1, 40, 44, 48, 52, 56, 60, 92, 98, 104, 110, 116, 122, 214, 174, 182, 190, 198, 206},
90 { 2, 22, 25, 28, 31, 34, 37, 62, 67, 72, 77, 82, 87, 166, 131, 138, 145, 152, 159},
91 { 8, 10, 12, 14, 16, 18, 8, 38, 42, 46, 50, 54, 58, 124, 94, 100, 106, 112, 118},
92 {20, 11, 13, 15, 17, 19, 9, 20, 23, 26, 29, 32, 35, 88, 63, 68, 73, 78, 83},
93 {30, 131, 138, 145, 152, 159, 166, 219, 228, 237, 246, 255, 264, 329, 279, 289, 299, 309, 319},
94 {15, 95, 101, 107, 113, 119, 125, 171, 179, 187, 195, 203, 211, 269, 224, 233, 242, 251, 260},
95 { 6, 65, 70, 75, 80, 85, 90, 129, 136, 143, 150, 157, 164, 215, 175, 183, 191, 199, 207},
96 { 7, 41, 45, 49, 53, 57, 61, 93, 99, 105, 111, 117, 123, 167, 132, 139, 146, 153, 160},
97 {19, 23, 26, 29, 32, 35, 20, 63, 68, 73, 78, 83, 88, 125, 95, 101, 107, 113, 119},
98 {37, 24, 27, 30, 33, 36, 21, 39, 43, 47, 51, 55, 59, 89, 64, 69, 74, 79, 84},
99 {31, 132, 139, 146, 153, 160, 167, 220, 229, 238, 247, 256, 265, 270, 225, 234, 243, 252, 261},
100 {16, 96, 102, 108, 114, 120, 126, 172, 180, 188, 196, 204, 212, 216, 176, 184, 192, 200, 208},
101 {17, 66, 71, 76, 81, 86, 91, 130, 137, 144, 151, 158, 165, 168, 133, 140, 147, 154, 161},
102 {18, 42, 46, 50, 54, 58, 38, 94, 100, 106, 112, 118, 124, 126, 96, 102, 108, 114, 120},
103 {36, 43, 47, 51, 55, 59, 39, 64, 69, 74, 79, 84, 89, 90, 65, 70, 75, 80, 85},
104 {32, 133, 140, 147, 154, 161, 168, 221, 230, 239, 248, 257, 266, 217, 177, 185, 193, 201, 209},
105 {33, 97, 103, 109, 115, 121, 127, 173, 181, 189, 197, 205, 213, 169, 134, 141, 148, 155, 162},
106 {35, 68, 73, 78, 83, 88, 63, 95, 101, 107, 113, 119, 125, 91, 66, 71, 76, 81, 86}
107 };
108
109//
110// corrispondence between pixels in cell and lut (pix numbering starts from 1)
111//
112const Int_t MMcTriggerLvl2::gsPixelsInLut[3][12] = {
113 { 1, 2, 3, 4, 6, 7, 8, 9, 12, 13, 14, 15},
114 {34, 29, 23, 16, 30, 24, 17, 10, 25, 18, 11, 5},
115 {35, 31, 26, 20, 19, 32, 27, 21, 36, 33, 28, 22},
116 };
117
118
119// --------------------------------------------------------------------------
120//
121// Default constructor
122//
123MMcTriggerLvl2::MMcTriggerLvl2(const char *name, const char *title)
124{
125 fName = name ? name : ClassName();
126 fTitle = title;
127
128 *fLog << "created MMcTriggerLvl2" << endl;
129
130 //
131 // Initialization of the fPixels array to zero
132 //
133 memset (fPixels,0,36*19*sizeof(Int_t));
134 memset (fFiredPixel,0,397*sizeof(Int_t));
135
136 // create a new camera
137 SetNewCamera(new MGeomCamMagic);
138}
139
140// --------------------------------------------------------------------------
141//
142// Destructor
143//
144MMcTriggerLvl2::~MMcTriggerLvl2()
145{
146 delete fCam;
147 delete fGeomCam;
148}
149
150// --------------------------------------------------------------------------
151//
152// Print functions for debugging and testing
153//
154// Options available:
155//
156// "cell":
157// Print out the pixel number of a given trigger cell
158//
159// "status":
160// Print a table with the status (pixel is fired or not after Lvl1) for
161// all the pixels in all the trigger cells
162//
163// no option:
164// Print the value of the selection parameters
165//
166//
167void MMcTriggerLvl2::Print(Option_t *opt) const
168{
169 TString str(opt);
170
171 if (str.Contains("status", TString::kIgnoreCase))
172 {
173 *fLog << " Status of cells 1-9" <<endl;
174 for(int i=0; i<36; i++)
175 {
176 for(int j=0; j<9; j++)
177 {
178 // *fLog.width(3);
179 *fLog <<gsPixelsInCell[i][j]-1 << ":" << fPixels[i][j] << "\t ";
180 }
181 *fLog << endl;
182 }
183 }
184 else if (str.Contains("cell", TString::kIgnoreCase))
185 {
186 *fLog << " Pixel numbering in cells 1-9" <<endl;
187 for (int i=0;i<36;i++)
188 {
189 for(int j=0; j<9; j++)
190 {
191 *fLog << gsPixelsInCell[i][j]-1 << "\t";
192 }
193 *fLog << endl;
194 }
195 }
196 else
197 {
198 *fLog << " L2T selection parameters:" << endl;
199 *fLog << " - LutPseudoSize = " << fLutPseudoSize << endl;
200 *fLog << " - PseudoSize = " << fPseudoSize << endl;
201 *fLog << " - BiggerCellSize = " << fSizeBiggerCell << endl;
202 }
203
204}
205
206// --------------------------------------------------------------------------
207//
208// Take the information supplied by the Lvl1 (it reads from the object MMcTrig
209// the status of all pixels after Lvl1) and pass it to the Lvl2 (fill the
210// array fPixels)
211//
212//
213void MMcTriggerLvl2::SetLv1(MMcTrig *trig)
214{
215 if (!trig)
216 return;
217
218 fMcTrig = trig;
219
220 //fMcTrig->PrintPixelsFirstLevel();
221
222 for(int i=0; i<36; i++)
223 {
224 for(int j=0; j<19; j++)
225 {
226 int pixel = gsPixelsInCell[i][j]-1;
227 fPixels[i][j] = (fMcTrig->IsPixelFired(pixel,0)) ? 1 : 0;
228 fFiredPixel[pixel]=(fMcTrig->IsPixelFired(pixel,0)) ? 1 : 0;
229 }
230 }
231}
232
233
234// --------------------------------------------------------------------------
235//
236// Set the trigger status ( pixel fired(=1) or not(=0) ) manually for a given
237// pixel of a given cell
238//
239//
240void MMcTriggerLvl2::SetPixelFired(Int_t pixel, Int_t fired)
241{
242 for(int i=0; i<36; i++)
243 {
244 for(int j=0; j<19; j++)
245 {
246 if(gsPixelsInCell[i][j]-1==pixel) fPixels[i][j]=fired;
247 }
248 }
249
250}
251
252// --------------------------------------------------------------------------
253//
254// Calculate the Level 2 Trigger (L2T) parameters. They can be used
255// to select the events that have passed the L2T selection rule.
256//
257//
258void MMcTriggerLvl2::Calc()
259{
260
261 // Find the Lut and cell with the higher LutPseudoSize.
262 int lutcell = CalcBiggerLutPseudoSize();
263 int maxcell = lutcell/10;
264 int maxlut = lutcell-maxcell*10;
265 fLutPseudoSize = GetLutCompactPixel(maxcell,maxlut);
266
267 CalcPseudoSize();
268
269 fSizeBiggerCell = GetCellNumberFired(CalcBiggerFiredCell());
270
271}
272
273
274// --------------------------------------------------------------------------
275//
276// Display the MAGIC Camera a draw the pixels corresponding to a given cell
277// (preliminary) <to be fixed>
278//
279void MMcTriggerLvl2::DrawCell(Int_t cell)
280{
281
282 if(cell>18) return;
283
284 //
285 // Use MHCamera class variable for avoiding to create a MHCamera each
286 // time this function is called. Also, now all the hexagons are drawn in
287 // the same camera pad
288 //
289 if (!fGeomCam)
290 fGeomCam = new MGeomCamMagic;
291 if (!fCam)
292 {
293 fCam = new MHCamera(*fGeomCam);
294 fCam->Draw();
295 fCam->DrawPixelNumbers();
296 }
297
298 //fCam->Draw();
299 fCam->Reset();
300
301 int color=0;
302
303 for(int i=0; i<36; i++)
304 {
305 color = (fPixels[i][cell]) ? 5 : 3;
306 fCam->SetPix( gsPixelsInCell[i][cell]-1, color, 1, 5 );
307 }
308
309 //
310 // Update the display (paint the camera with the new colors)
311 //
312 gPad->Modified();
313 gPad->Update();
314}
315
316
317
318// --------------------------------------------------------------------------
319//
320// Display the MAGIC camera and draw all the pixel fired after Lvl1
321// (preliminary) <to be fixed>
322//
323//
324void MMcTriggerLvl2::DrawLv1()
325{
326 //
327 // Use MHCamera class variable for avoiding to create a MHCamera each
328 // time this function is called. Also, now all the hexagons are drawn in
329 // the same camera pad
330 //
331 // !FixMe! Delete this method when the trigger display is available!
332 //
333 if (!fGeomCam)
334 fGeomCam = new MGeomCamMagic;
335 if (!fCam)
336 {
337 fCam = new MHCamera(*fGeomCam);
338 fCam->Draw();
339 fCam->DrawPixelNumbers();
340 }
341
342 //
343 // Set the array of colors for each pixel (that will be painted after
344 // updating the dispaly)
345 //
346 int color=0;
347
348 for(int i=0; i<577; i++)
349 {
350 color = (fMcTrig->IsPixelFired(i,0)) ? 1 : 2;
351 fCam->SetPix( i, color, 1, 5 );
352 }
353
354 //
355 // Update the display (paint the camera with the new colors)
356 //
357 gPad->Modified();
358 gPad->Update();
359
360}
361
362// --------------------------------------------------------------------------
363//
364// For a given cell, just count how many pixels have been fired after Lvl1
365//
366Int_t MMcTriggerLvl2::GetCellNumberFired(int cell)
367{
368 int size=0;
369
370 for(int i=0; i<36; i++)
371 {
372 size += fPixels[i][cell];
373 }
374
375 return size;
376
377}
378
379
380// --------------------------------------------------------------------------
381//
382// Find the cell which the bigger number of fired pixels
383//
384//
385Int_t MMcTriggerLvl2::CalcBiggerFiredCell()
386{
387 int size=-1;
388 int cell=-1;
389
390 for(int j=0; j<19; j++)
391 {
392 if (GetCellNumberFired(j) > size)
393 {
394 size = GetCellNumberFired(j);
395 cell = j;
396 }
397 }
398
399 // *fLog << "size " <<size <<" in cell " << cell << endl;
400
401 return cell;
402
403}
404
405// --------------------------------------------------------------------------
406//
407// Calculate the higher LutPseudoSize of one event defined as the higher number
408// of compact pixel in the LUTs of the trigger cells.
409// neighpix is the number of Next-Neighbors which defines the compact pixel
410// topology (it can be 2-NN or 3-NN)
411//
412// It returns the cell and the LUT with the bigger LutPseudoSize, coded
413// accordingly to: cell*10 + LUT
414//
415Int_t MMcTriggerLvl2::CalcBiggerLutPseudoSize()
416{
417 int size=0;
418 int cell=-1;
419 int lut=-1;
420
421 for(int j=0; j<19; j++)
422 {
423 for(int i=0; i<3; i++)
424 {
425 if (GetLutCompactPixel(j,i) >= size)
426 {
427 size = GetLutCompactPixel(j,i);
428 cell = j;
429 lut = i;
430 }
431 }
432 }
433
434 // *fLog <<"Max cell: " << cell+1 << " Max Lut: " << lut+1 << " PseudoSize: " << size <<endl;
435
436 return cell*10+lut;
437}
438
439
440// --------------------------------------------------------------------------
441//
442// Identify compact pixels in one LUT and calculate the LutPseudoSize
443// (=number of compact pixels in one LUT)
444// neighpix is the number of Next-Neighbors which defines the compact pixel
445// topology.
446// Up to now only 2NN or 3NN are considered
447// <!if changed: change check made by Preprocess in MMcTriggerLvl2Calc>
448//
449// Returns:
450// -1 wrong neighpix
451// -2 wrong cell (<1 or >19)
452// -3 wrong lut (<1 or >3)
453// else:
454// the number of compact pixels in one LUT.
455//
456Int_t MMcTriggerLvl2::GetLutCompactPixel(int cell, int lut)
457{
458 int size=0;
459 int lutpix, a[12];
460 int neighpix= (*this).fCompactNN;
461
462 // check on input variables
463
464 if (neighpix >3 || neighpix < 2)
465 return(-1);
466
467 if (cell<0 || cell>18)
468 return(-2);
469
470 if (lut<0 || lut>2)
471 return(-3);
472
473
474 // LUT 1 and 2 are similar; LUT 3 differs
475
476 for(int j=0; j<12; j++)
477 {
478 lutpix = gsPixelsInLut[lut][j]-1;
479 // *fLog << "j=" <<j<<" lutpix="<<lutpix<<" Cell="<<cell<<endl;
480 a[j] = fPixels[lutpix][cell];
481 }
482
483 //
484 // Look for compact pixels 2NN
485 //
486 if (neighpix==2)
487 {
488 // case Lut 1 and 2
489 if (lut == 0 || lut == 1)
490 {
491 if (a[0] && a[1] && a[4])
492 size ++;
493 if (a[1] && ((a[0] && a[4]) ||
494 (a[4] && a[5]) ||
495 (a[5] && a[2]) ))
496 size ++;
497 if (a[2] && ((a[1] && a[5]) ||
498 (a[5] && a[6]) ||
499 (a[6] && a[3]) ))
500 size ++;
501 if (a[3] && ((a[2] && a[6]) ||
502 (a[6] && a[7]) ))
503 size ++;
504 if (a[4] && ((a[0] && a[1]) ||
505 (a[1] && a[5]) ||
506 (a[5] && a[8]) ))
507 size ++;
508 if (a[5] && ((a[1] && a[2]) ||
509 (a[2] && a[6]) ||
510 (a[6] && a[9]) ||
511 (a[9] && a[8]) ||
512 (a[8] && a[4]) ||
513 (a[4] && a[1]) ))
514 size ++;
515 if (a[6] && ((a[2] && a[3]) ||
516 (a[3] && a[7]) ||
517 (a[7] && a[10]) ||
518 (a[10] && a[9]) ||
519 (a[9] && a[5]) ||
520 (a[5] && a[2]) ))
521 size ++;
522 if (a[7] && ((a[3] && a[6]) ||
523 (a[6] && a[10]) ||
524 (a[10] && a[11]) ))
525 size ++;
526 if (a[8] && ((a[4] && a[5]) ||
527 (a[5] && a[9]) ))
528 size ++;
529 if (a[9] && ((a[8] && a[5]) ||
530 (a[5] && a[6]) ||
531 (a[6] && a[10]) ))
532 size ++;
533 if (a[10] && ((a[9] && a[6]) ||
534 (a[6] && a[7]) ||
535 (a[7] && a[11]) ))
536 size ++;
537 if (a[11] && a[7] && a[10])
538 size ++;
539 }
540
541 // case Lut 3
542 if (lut==2)
543 {
544 if (a[0] && a[1] && a[5])
545 size ++;
546 if (a[1] && ((a[0] && a[5]) ||
547 (a[5] && a[2]) ))
548 size ++;
549 if (a[2] && ((a[1] && a[5]) ||
550 (a[5] && a[6]) ||
551 (a[3] && a[4]) ||
552 (a[6] && a[3]) ))
553 size ++;
554 if (a[3] && ((a[2] && a[6]) ||
555 (a[6] && a[7]) ||
556 (a[2] && a[4]) ))
557 size ++;
558 if (a[4] && ((a[2] && a[3]) ))
559 size ++;
560 if (a[5] && ((a[1] && a[2]) ||
561 (a[2] && a[6]) ||
562 (a[6] && a[9]) ||
563 (a[9] && a[8]) ))
564 size ++;
565 if (a[6] && ((a[2] && a[3]) ||
566 (a[3] && a[7]) ||
567 (a[7] && a[10]) ||
568 (a[10] && a[9]) ||
569 (a[9] && a[5]) ||
570 (a[5] && a[2]) ))
571 size ++;
572 if (a[7] && ((a[3] && a[6]) ||
573 (a[6] && a[10]) ||
574 (a[10] && a[11]) ))
575 size ++;
576 if (a[8] && a[5] && a[9])
577 size ++;
578 if (a[9] && ((a[8] && a[5]) ||
579 (a[5] && a[6]) ||
580 (a[6] && a[10]) ))
581 size ++;
582 if (a[10] && ((a[9] && a[6]) ||
583 (a[6] && a[7]) ||
584 (a[7] && a[11]) ))
585 size ++;
586 if (a[11] && a[7] && a[10])
587 size ++;
588 }
589 }
590
591 //
592 // Look for compact pixels 3NN
593 //
594 if (neighpix==3)
595 {
596 // case Lut 1 and 2
597 if (lut == 0 || lut == 1)
598 {
599 if (a[0] && a[1] && a[4] && a[5]) // pix 0 is compact if there is a 4NN
600 size ++;
601 if (a[1] && ((a[0] && a[4] && a[5]) ||
602 (a[2] && a[4] && a[5]) ))
603 size ++;
604 if (a[2] && ((a[1] && a[5] && a[6]) ||
605 (a[5] && a[6] && a[3]) ))
606 size ++;
607 if (a[3] && (a[2] && a[6] && a[7] ))
608 size ++;
609 if (a[4] && ((a[0] && a[1] && a[5]) ||
610 (a[1] && a[5] && a[8]) ))
611 size ++;
612 if (a[5] && ((a[1] && a[2] && a[6]) ||
613 (a[2] && a[6] && a[9]) ||
614 (a[6] && a[9] && a[8]) ||
615 (a[9] && a[8] && a[4]) ||
616 (a[8] && a[4] && a[1]) ||
617 (a[4] && a[1] && a[2]) ))
618 size ++;
619 if (a[6] && ((a[2] && a[3] && a[7]) ||
620 (a[3] && a[7] && a[10]) ||
621 (a[7] && a[10] && a[9]) ||
622 (a[10] && a[9] && a[5]) ||
623 (a[9] && a[5] && a[2]) ||
624 (a[5] && a[2] && a[3]) ))
625 size ++;
626 if (a[7] && ((a[3] && a[6] && a[10]) ||
627 (a[6] && a[10] && a[11]) ))
628 size ++;
629 if (a[8] && (a[4] && a[5] && a[9] ))
630 size ++;
631 if (a[9] && ((a[8] && a[5] && a[6]) ||
632 (a[5] && a[6] && a[10]) ))
633 size ++;
634 if (a[10] && ((a[9] && a[6] && a[7]) ||
635 (a[6] && a[7] && a[11]) ))
636 size ++;
637 if (a[11] && a[7] && a[10] && a[6]) //pix 0 is compact if there is a 4NN
638 size ++;
639 }
640
641 // case Lut 3
642 if (lut==2)
643 {
644 if (a[0] && a[1] && a[5] && a[8]) // pix 0 is compact if there is a 4NN
645 size ++;
646 if (a[1] && (a[0] && a[5] && a[2])) //pix 0 is compact if there is a 4NN
647 size ++;
648 if (a[2] && ((a[1] && a[5] && a[6]) ||
649 (a[3] && a[5] && a[6]) ||
650 (a[6] && a[3] && a[4]) ))
651 size ++;
652 if (a[3] && ((a[2] && a[4] && a[6]) ||
653 (a[2] && a[6] && a[7]) ))
654 size ++;
655 if (a[4] && (a[2] && a[3] && a[6] ))
656 size ++;
657 if (a[5] && ((a[1] && a[2] && a[6]) ||
658 (a[2] && a[6] && a[9]) ||
659 (a[6] && a[9] && a[8]) ))
660 size ++;
661 if (a[6] && ((a[2] && a[3] && a[7]) ||
662 (a[3] && a[7] && a[10]) ||
663 (a[7] && a[10] && a[9]) ||
664 (a[10] && a[9] && a[5]) ||
665 (a[9] && a[5] && a[2]) ||
666 (a[5] && a[2] && a[3]) ))
667 size ++;
668 if (a[7] && ((a[3] && a[6] && a[10]) ||
669 (a[6] && a[10] && a[11]) ))
670 size ++;
671 if (a[8] && a[5] && a[9] && a[6]) //pix 0 is compact if there is a 4NN
672 size ++;
673 if (a[9] && ((a[8] && a[5] && a[6]) ||
674 (a[5] && a[6] && a[10]) ))
675 size ++;
676 if (a[10] && ((a[9] && a[6] && a[7]) ||
677 (a[6] && a[7] && a[11]) ))
678 size ++;
679 if (a[11] && a[7] && a[10] && a[6]) //pix 0 is compact if there is a 4NN
680 size ++;
681 }
682 }
683
684
685 if(size<0 ||size>12)
686 *fLog << "*" << size <<"-";
687
688 return size;
689}
690
691
692// --------------------------------------------------------------------------
693//
694// Look for clusters and calculate the PseudoSize of one event,
695// defined as the multiplicity of the bigger cluster.
696//
697//
698//
699void MMcTriggerLvl2::CalcPseudoSize()
700{
701 // Fill the fCompactPixel array with the compact pixels
702 CalcCompactPixels(fGeomCam);
703
704 // seek the LUT with higher number of compact pixels
705 //
706 int cellut = CalcBiggerLutPseudoSize();
707 int maxcell = cellut/10;
708 int maxlut = cellut - maxcell*10;
709 int startpix;
710
711 //
712 // seek a starting pixel in the lut for the iteration
713 //
714 int check=1;
715 for (int pixlut=0;pixlut<12;pixlut++)
716 {
717 int pixcell =gsPixelsInLut[maxlut][pixlut];
718 startpix = gsPixelsInCell[pixcell][maxcell]-1;
719 *fLog << "pix, compact:" << startpix << "@"<<fCompactPixel[startpix];
720 if (fCompactPixel[startpix]) // a starting pixel was found
721 break;
722 check++;
723 }
724
725 // A LUT contains 12 pixels
726 if (check >= 12)
727 {
728 *fLog <<"Error: a starting pixels was not found!"<<endl;
729 return;
730 }
731 //
732 // Bulding cluster
733 //
734 Int_t cluster[397];
735 int pnt=0;
736 int pnt2=0; //pointer in the array fCluster_pix, needed for filling
737
738 memset (cluster,0,397*sizeof(Int_t));
739 memset (fCluster_pix,-1,397*sizeof(Int_t));
740
741 cluster[startpix]=1;
742 fCluster_pix[0]=startpix; //the starting pix is the first in cluster
743
744 // Look at neighbour pixs if they are compact (iterative search)
745 // until the array (fCluster_pix) has no more compact pixels.
746 // pnt points to the pixel in the array fCluster_pix whose neighbors are
747 // under study; pnt2 points to the last element of this array.
748 //
749 while (fCluster_pix[pnt] != -1)
750 {
751 const MGeomPix &pix=(*fGeomCam)[fCluster_pix[pnt]];
752
753 for (int i=0;i<pix.GetNumNeighbors();i++)
754 {
755 int pix_neigh = pix.GetNeighbor(i);
756 // check if pixel is fired and doesn't belong to cluster
757 if (fFiredPixel[pix_neigh] && !cluster[pix_neigh])
758 {
759 cluster[pix_neigh] = 1;
760 fCluster_pix[++pnt2] = pix_neigh;
761 }
762 }
763 pnt++;
764 }
765
766 fPseudoSize = pnt;
767
768 // *fLog << "ClusterID:" <<(*clust).GetClusterId() << " Mult:" << (*clust).GetMultiplicity()<<endl;
769
770 *fLog <<"PSize: "<< fPseudoSize << " in cell:" << maxcell << " lut:" <<maxlut <<endl;
771
772 return;
773}
774
775// --------------------------------------------------------------------------
776//
777// Fill the fCompactPixels array with the pixels which are compact
778//
779// neighpix is the number of Next-Neighbors which defines the compact pixel
780// topology (it can be 2-NN or 3-NN)
781//
782// Note: it is a *global* method; it looks in the all camera as a whole
783//
784//
785void MMcTriggerLvl2::CalcCompactPixels(MGeomCam *geom)
786{
787 memset (fCompactPixel,0,397*sizeof(Int_t));
788
789 for(UInt_t pixid=0;pixid<397;pixid++)
790 {
791 // Look if the pixel is fired, otherwise continue
792 if (!fFiredPixel[pixid])
793 continue;
794
795 const MGeomPix &pix=(*geom)[pixid];
796
797 // Look for compact pixels.
798 // A compact pixel must have at least fCompactNN adjacent neighbors
799 // It checks the 6 different configurations of neighbors pixels
800 for (int i=0;i<pix.GetNumNeighbors();i++)
801 {
802 int j=0;
803 while (fFiredPixel[pix.GetNeighbor(i+j)]==1 && j < fCompactNN)
804 j++;
805 if (j!=fCompactNN) continue; // configuration doesn't satisfy compact condition
806
807 fCompactPixel[pixid]=1; // pixel is compact
808
809 *fLog << ","<<pixid;
810
811 break;
812 }
813 }
814}
815
816
817
Note: See TracBrowser for help on using the repository browser.