source: JEM-EUSO/esaf_cc_at_lal/packages/simulation/detector/electronics/src/LblChipTrackingTrgEngine.cc @ 114

Last change on this file since 114 was 114, checked in by moretto, 11 years ago

actual version of ESAF at CCin2p3
File size: 17.0 KB
Line 
1// ESAF : Euso Simulation and Analysis Framework
2// $Id: LblChipTrackingTrgEngine.cc 2168 2005-10-10 16:28:03Z ejudd $
3
4#include "LblChipTrackingTrgEngine.hh"
5#include "MacroCellData.hh"
6#include "MacroCell.hh"
7#include "ElementaryCell.hh"
8#include "FrontEndChip.hh"
9#include "BoolAlgebra.hh"
10#include "LblTrackSegment.hh"
11#include "EEvent.hh"
12#include "ELblTrackTriggerDataAdder.hh"
13
14ClassImp(LblChipTrackingTrgEngine)
15
16//______________________________________________________________________________
17LblChipTrackingTrgEngine::LblChipTrackingTrgEngine() : 
18TriggerEngine(string("LblChipTrackingTrgEngine"), kLblChipTrackingTrigger ) {
19
20    // constructor
21
22    // Read all the control flags from the LblChipTrackingTrgEngine.cfg file
23
24    // First the flag to control which pxiels get used to make track segments
25    fPxlthr = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Pxlthr");
26
27    // Now the flags to control how track segments are found
28    fColumn = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Column");
29    fColthr = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Colthr");
30    fFrmsinseg = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Frmsinseg");
31    fSegs = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Segs");
32    fTrksumthr = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Trksumthr");
33    fWantsumsub = (bool)Conf()->GetNum("LblChipTrackingTrgEngine.Wantsumsub");
34
35    // Next the flags to control how segments are joined to form tracks
36    // and tracks are tested for persistency
37    fFrmsinpsum = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Frmsinpsum");
38    fMinpersum = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Minpersum");
39
40    // Finally the flags to control the level of output messages
41    fVerbose = (bool)Conf()->GetNum("LblChipTrackingTrgEngine.Verbose");
42    fDebug = (bool)Conf()->GetNum("LblChipTrackingTrgEngine.Debug");
43
44    if (fVerbose) {
45        Msg(EsafMsg::Info) << "fPxlthr = " << fPxlthr << MsgDispatch;
46
47        Msg(EsafMsg::Info) << "fColumn = " << fColumn << MsgDispatch;
48        Msg(EsafMsg::Info) << "fColthr = " << fColthr << MsgDispatch;
49        Msg(EsafMsg::Info) << "fFrmsinseg = " << fFrmsinseg << MsgDispatch;
50        Msg(EsafMsg::Info) << "fSegs = " << fSegs << MsgDispatch;
51        Msg(EsafMsg::Info) << "fTrksumthr = " << fTrksumthr << MsgDispatch;
52        Msg(EsafMsg::Info) << "fWantsumsub = " << fWantsumsub << MsgDispatch;
53
54        Msg(EsafMsg::Info) << "fFrmsinpsum = " << fFrmsinpsum << MsgDispatch;
55        Msg(EsafMsg::Info) << "fMinpersum = " << fMinpersum << MsgDispatch;
56
57        Msg(EsafMsg::Info) << "fDebug = " << fDebug << MsgDispatch;
58    }
59
60}
61
62
63//______________________________________________________________________________
64LblChipTrackingTrgEngine::~LblChipTrackingTrgEngine() {
65   
66    // dtor
67
68    // FIXME - hardwired array sizes are BAD
69    Clear(199,12);
70
71}
72
73
74//______________________________________________________________________________
75void LblChipTrackingTrgEngine::Simulate( MacroCellData* pData) {   
76
77    Int_t tot_tracks;
78    Int_t ec_tracks;
79    Int_t tot_trigs;
80    Int_t ec_trigs;
81
82    // simulate chip tracking algorithm
83    if (fVerbose)
84        Msg(EsafMsg::Info) << "LBL CHIP TRACKING STARTED MACROCELL " << pData->Cell()->Id() << MsgDispatch;
85
86    // Get maps of ChipGtuData for this MacroCell
87    // m1 is a reference (i.e. a short alias) for pData->fChipData2
88    // fChipData2 is a map: first element is the Id of this Front End
89    //                      second element is a pointer to a second map
90    map<Int_t, map<Int_t,ChipGtuData*>* > &m1 = pData->fChipData2;
91    map<Int_t, map<Int_t,ChipGtuData*>* >::const_iterator it1;
92    Int_t ncells = m1.size();
93    if (fVerbose)
94        Msg(EsafMsg::Info) << "MacroCell contains " << ncells << " Elementary Cells" << MsgDispatch;
95
96    // Quit if there is no data in this map element
97    if (!ncells) return;
98
99    // Now do the main loop over the Elementary Cells, find segments and tracks
100    // then look for persistency.
101    tot_tracks = 0;
102    tot_trigs = 0;
103    for (it1=m1.begin(); it1!=m1.end(); it1++) {
104
105        if (fVerbose)
106            Msg(EsafMsg::Info) << "Procesing Elementary Cell " << it1->first << MsgDispatch;
107
108        // Get maps of ChipGtuData for this specific Elementary Cell
109        // m2 is a reference to the map object that it1->second points to
110        // This map has 2 elements: first element is GTU number
111        //                          second element is a pointer to ChipGtuData for this GTU
112        map<Int_t,ChipGtuData*> &m2 = *(it1->second);
113        map<Int_t,ChipGtuData*>::const_iterator it2;
114
115        // Check to see if this Elementary Cell has enough GTUs
116        it2 = m2.begin();
117        Int_t first_gtu = it2->first;
118        it2 = m2.end();
119        it2--;
120        Int_t last_gtu = it2->first;
121        Int_t gtot = last_gtu - first_gtu + 1;
122        if (fVerbose) Msg(EsafMsg::Info) << "First GTU = " << first_gtu << " Last GTU = " << last_gtu << " Elementary Cell contains " << gtot << " GTUs" << MsgDispatch;
123        if (last_gtu >= 200) {
124            if (fVerbose) Msg(EsafMsg::Info) << "Last GTU is too late" << MsgDispatch;
125            continue;
126        }
127        if (gtot < (fSegs*fFrmsinseg)) {
128            if (fVerbose) Msg(EsafMsg::Info) << "Not enough GTUs in this Elementary Cell" << MsgDispatch;
129            continue;
130        }
131
132        // Get the number of channels in this chip from the first GTU
133        it2 = m2.begin();
134        Int_t chan = (it2->second)->FrontEnd()->Channels();
135        Int_t NPIX = (it2->second)->FrontEnd()->NumSide();
136        if (fDebug)
137            Msg(EsafMsg::Debug) << "Number of channels is " << chan << " Number per side is " << NPIX << MsgDispatch;
138
139        // Zero the arrays
140        Clear(last_gtu,NPIX);
141
142        // Find the tracks and segments
143        ec_tracks = Find_Tracks(m2,last_gtu,chan,NPIX);
144
145        // Check to see if there is a good track
146        ec_trigs = Check_Persistency(last_gtu,NPIX);
147
148        if (fVerbose)
149            Msg(EsafMsg::Info) << "EC " << it1->first << " generated " << ec_trigs << " triggers from " << ec_tracks << " tracks" << MsgDispatch;
150
151        tot_tracks += ec_tracks;
152        tot_trigs += ec_trigs;
153
154        // Save the tracks in the root output
155        FillELblTrackTrigger(pData->Cell()->Id(),it1->first);
156
157    } // End of it1 loop over Elementary Cells
158
159    // Print total track count
160    if (fVerbose)
161        Msg(EsafMsg::Info) << tot_tracks << " tracks were found in cell " << pData->Cell()->Id() << MsgDispatch;
162
163    // Set trigger flags and save count if any triggers were issued
164    if (tot_trigs > 0) {
165        SetTrigger(tot_trigs);
166        Msg(EsafMsg::Info) << tot_trigs << " triggers have been generated in cell " << pData->Cell()->Id() << MsgDispatch;
167    }
168
169}
170
171//______________________________________________________________________________
172void LblChipTrackingTrgEngine::Clear(Int_t last_gtu,Int_t NPIX) {
173   
174    // Clear the fRawSums, fSums, fPers and fTrks arrays
175    for (Int_t i=0; i<(last_gtu+1); i++) {
176        for (Int_t j=0; j<NPIX; j++) {
177            for (Int_t k=0; k<NPIX; k++) {
178                fRawSums[i][j][k] = 0;
179                fSums[i][j][k] = 0;
180                fPers[i][j][k] = 0;
181                fTrks[i][j][k] = -1;
182            }
183        }
184    }
185
186    // Clear the track vectors
187    fTracks.clear();
188}
189
190
191//______________________________________________________________________________
192Int_t LblChipTrackingTrgEngine::Find_Tracks(map<Int_t,ChipGtuData*> &m2,Int_t last_gtu,Int_t chan,Int_t NPIX) {
193
194    Int_t gt;
195    Int_t last_track_gtu;
196    Int_t row;
197    Int_t col;
198    Int_t colsum;
199    Int_t sub_row;
200    Int_t sub_col;
201    Int_t sub_chan;
202    ChipGtuData* gtus;
203    Int_t seg1sum;
204    Int_t seg2sum;
205    Int_t seg3sum;
206    Int_t trksum;
207    Int_t row_l;
208    Int_t row_h;
209    Int_t col_l;
210    Int_t col_h;
211    Int_t sub_row2;
212    Int_t sub_col2;
213    Int_t ntracks;
214
215    if (fDebug) {
216        Msg(EsafMsg::Debug) << "Number of channels is " << chan << " Number per side is " << NPIX << MsgDispatch;
217        Msg(EsafMsg::Debug) << "last_gtu = " << last_gtu << MsgDispatch;
218    }
219
220    // Create the iterator for this map
221    // it2 is used to loop over all GTUs
222    map<Int_t,ChipGtuData*>::const_iterator it2;
223
224    // Main loop over GTUs
225    ntracks = 0;
226    for (it2=m2.begin(); it2!=m2.end(); it2++) {
227
228        // Get the data for this specific GTU
229        gt = it2->first;
230        gtus = it2->second;
231        if (fDebug)
232            Msg(EsafMsg::Debug) << "Building Segments starting at GTU Number " << gt << MsgDispatch;
233
234        // Add this data to the sums for each pixel
235        for (Int_t nchan=0; nchan<chan; nchan++) {
236
237            row = gtus->FrontEnd()->Row(nchan);
238            col = gtus->FrontEnd()->Column(nchan);
239
240            // Now calculate the sums
241            switch (fColumn) {
242            case 2:
243                colsum = 0;
244                for (Int_t i=0; i<2; i++) {
245                    sub_row = row+i;
246                    for (Int_t j=0; j<2; j++) {
247                        sub_col = col+j;
248                        if ((sub_row<NPIX)&&(sub_col<NPIX)) {
249                            sub_chan = gtus->FrontEnd()->ChanRowCol(sub_row,sub_col);
250                            if (gtus->GetCounter(sub_chan) > fPxlthr) {
251                                colsum += gtus->GetCounter(sub_chan);
252                            }
253                        }
254                    }
255                }
256                break;
257            case 3:
258                colsum = 0;
259                for (Int_t i=0; i<3; i++) {
260                    sub_row = row+i;
261                    for (Int_t j=0; j<3; j++) {
262                        sub_col = col+j;
263                        if ((sub_row<NPIX)&&(sub_col<NPIX)) {
264                            sub_chan = gtus->FrontEnd()->ChanRowCol(sub_row,sub_col);
265                            if (gtus->GetCounter(sub_chan) > fPxlthr) {
266                                colsum += gtus->GetCounter(sub_chan);
267                            }
268                        }
269                    }
270                }
271                break;
272            default:
273                colsum = 0;
274                // Get this pixel's counts, but only if # counts > threshold
275                if (gtus->GetCounter(nchan) > fPxlthr)
276                    colsum += gtus->GetCounter(nchan);
277                break;
278            } // End of switch (fColumn)
279
280            // Loop over active segments and add in this column
281            Int_t num_active = ((gt+1)<fFrmsinseg) ? (gt+1) : fFrmsinseg;
282            for (Int_t active=0; active<num_active; active++) {
283                fRawSums[gt-active][row][col] += (colsum>fColthr) ? colsum : 0;
284            }
285
286        } // End of nchan loop over Channels (i.e. Pixels)
287
288    } // End of loop over map m2. At this point we have the sums for each pixel for all GTUs.
289
290    // Now either subtract off the background counts, if the fWantsumsub flag is true,
291    // or just copy fRawSums to fSums if fWantsumsub is false.
292    for (gt=0; gt<(last_gtu+1); gt++) {
293        for (row=0; row<NPIX; row++) {
294            for (col=0; col<NPIX; col++) {
295                if (fWantsumsub) {
296                    sub_row = (row + 4)%NPIX;
297                    sub_col = (col + 4)%NPIX;
298                    Int_t sumsub = fRawSums[gt][row][col] - fRawSums[gt][sub_row][sub_col];
299                    fSums[gt][row][col] = (sumsub>0) ? sumsub : 0;
300                }
301                else {
302                    fSums[gt][row][col] = fRawSums[gt][row][col];
303                }
304            }
305        }
306    }
307
308    // Now work with partial sums to find tracksums
309    last_track_gtu = (last_gtu+1) - (fSegs*fFrmsinseg);
310    for (gt=0; gt<(last_track_gtu+1); gt++) {
311        for (row=0; row<NPIX; row++) {
312            for (col=0; col<NPIX; col++) {
313
314                seg1sum = fSums[gt][row][col];
315
316                // Single segment tracks
317                if (fSegs == 1) {
318                    trksum = seg1sum;
319                    if (trksum > fTrksumthr) {
320                        fPers[gt][row][col] += 1;
321                        if (fTrks[gt][row][col] == -1) fTrks[gt][row][col] = ntracks;
322                        ntracks++;
323                        // Save the track information
324                        fTracks.push_back( LblTrackSegment() );
325                        LblTrackSegment &trkseg = fTracks.back();
326                        trkseg.SetGtuStart(gt);
327                        trkseg.SetRowStart(row);
328                        trkseg.SetColStart(col);
329                        trkseg.SetGtuEnd(gt+fFrmsinseg-1);
330                        trkseg.SetRowEnd(row);
331                        trkseg.SetColEnd(col);
332                        trkseg.SetNumSegs(1);
333                        trkseg.SetSum(trksum);
334                        if (fDebug) {
335                            Msg(EsafMsg::Debug) << "trksum " << trksum << " > thr " << fTrksumthr << " for single segment track, fPers = " << fPers[gt][row][col] << " fTrks = " << fTrks[gt][row][col] << MsgDispatch;
336                            Msg(EsafMsg::Debug) << "Single segment start-GTU = " << gt << " start-Row = " << row << " start-Col = " << col << MsgDispatch;
337                        }
338                    }
339                }
340
341                // Two-or-more segment tracks
342                if (fSegs >= 2) {
343                    row_l = (row>0)        ? row-1 : 0;
344                    row_h = (row<(NPIX-1)) ? row+1 : (NPIX-1);
345                    col_l = (col>0)        ? col-1 : 0;
346                    col_h = (col<(NPIX-1)) ? col+1 : (NPIX-1);
347
348                    for (sub_row=row_l; sub_row<=row_h; sub_row++) {
349                        for (sub_col=col_l; sub_col<=col_h; sub_col++) {
350                           
351                            seg2sum = fSums[gt+fFrmsinseg][sub_row][sub_col];
352
353                            // Two segment tracks
354                            if (fSegs == 2) {
355                                trksum = seg1sum + seg2sum;
356                                if (trksum > fTrksumthr) {
357                                    fPers[gt][row][col] += 1;
358                                    if (fTrks[gt][row][col] == -1) fTrks[gt][row][col] = ntracks;
359                                    ntracks++;
360                                    // Save the track information
361                                    fTracks.push_back( LblTrackSegment() );
362                                    LblTrackSegment &trkseg = fTracks.back();
363                                    trkseg.SetGtuStart(gt);
364                                    trkseg.SetRowStart(row);
365                                    trkseg.SetColStart(col);
366                                    trkseg.SetGtuEnd(gt+(2*fFrmsinseg)-1);
367                                    trkseg.SetRowEnd(sub_row);
368                                    trkseg.SetColEnd(sub_col);
369                                    trkseg.SetNumSegs(2);
370                                    trkseg.SetSum(trksum);
371                                    if (fDebug) {
372                                        Msg(EsafMsg::Debug) << "trksum " << trksum << " > thr " << fTrksumthr << " for two segment track, fPers = " << fPers[gt][row][col] << " fTrks = " << fTrks[gt][row][col] << MsgDispatch;
373                                        Msg(EsafMsg::Debug) << "First  segment start-GTU = " << gt << " start-Row = " << row << " start-Col = " << col << MsgDispatch;
374                                        Msg(EsafMsg::Debug) << "Second segment start-GTU = " << gt+fFrmsinseg << " start-Row = " << sub_row << " start-Col = " << sub_col << MsgDispatch;
375                                    }
376                                }
377                            }
378
379                            // Three segment tracks
380                            if (fSegs == 3) {
381                                // Find sub_row2 and sub_col2 of 3rd segment by straight
382                                // line extrapolation from row/col and sub_row/sub_col
383                                sub_row2 = sub_row + (sub_row - row);
384                                if ((sub_row2 < 0) || (sub_row2 >= NPIX)) continue;
385                                sub_col2 = sub_col + (sub_col - col);
386                                if ((sub_col2 < 0) || (sub_col2 >= NPIX)) continue;
387
388                                seg3sum = fSums[gt+(2*fFrmsinseg)][sub_row2][sub_col2];
389
390                                trksum = seg1sum + seg2sum + seg3sum;
391                                if (trksum > fTrksumthr) {
392                                    fPers[gt][row][col] += 1;
393                                    if (fTrks[gt][row][col] == -1) fTrks[gt][row][col] = ntracks;
394                                    ntracks++;
395                                    // Save the track information
396                                    fTracks.push_back( LblTrackSegment() );
397                                    LblTrackSegment &trkseg = fTracks.back();
398                                    trkseg.SetGtuStart(gt);
399                                    trkseg.SetRowStart(row);
400                                    trkseg.SetColStart(col);
401                                    trkseg.SetGtuEnd(gt+(3*fFrmsinseg)-1);
402                                    trkseg.SetRowEnd(sub_row2);
403                                    trkseg.SetColEnd(sub_col2);
404                                    trkseg.SetNumSegs(3);
405                                    trkseg.SetSum(trksum);
406                                    if (fDebug) {
407                                        Msg(EsafMsg::Debug) << "trksum " << trksum << " > thr " << fTrksumthr << " for three segment track, fPers = " << fPers[gt][row][col] << " fTrks = " << fTrks[gt][row][col] << MsgDispatch;
408                                        Msg(EsafMsg::Debug) << "First  segment start-GTU = " << gt << " start-Row = " << row << " start-Col = " << col << MsgDispatch;
409                                        Msg(EsafMsg::Debug) << "Second segment start-GTU = " << gt+fFrmsinseg << " start-Row = " << sub_row << " start-Col = " << sub_col << MsgDispatch;
410                                        Msg(EsafMsg::Debug) << "Third  segment start-GTU = " << (gt+(2*fFrmsinseg)) << " start-Row = " << sub_row2 << " start-Col = " << sub_col2 << MsgDispatch;
411                                    }
412                                }
413                            } // End of fSegs == 3 if statement
414
415                        } // End of sub_col loop
416                    } // End of sub_row loop
417                } // End of fSegs >= 2 if statement
418            } // End of col loop
419        } // End of row loop
420    } // End of gt loop over all GTUs
421   
422    return ntracks;
423
424}
425
426
427//______________________________________________________________________________
428Int_t LblChipTrackingTrgEngine::Check_Persistency(Int_t last_gtu,Int_t NPIX) {
429    Int_t persum;
430    Int_t maxpsum;
431    Int_t gtutrig;
432    bool above;
433    bool lastabove;
434    Int_t ntrigs;
435   
436    if (fDebug)
437        Msg(EsafMsg::Debug) << "last_gtu = " << last_gtu << " NPIX = " << NPIX << MsgDispatch;
438
439    ntrigs = 0;
440    for (Int_t row=0; row<NPIX; row++) {
441        for (Int_t col=0; col<NPIX; col++) {
442            maxpsum = 0;
443            gtutrig = 0;
444            lastabove = false;
445
446            for (Int_t gtu=(fFrmsinpsum-1); gtu<(last_gtu+1); gtu++) {
447
448                persum = 0;
449                for (Int_t i=0; i<fFrmsinpsum; i++) {
450                    if (fPers[gtu-i][row][col] > 0) {
451                        persum++;
452                        if (fDebug) Msg(EsafMsg::Debug) << "GTU " << gtu << " Row " << row << " Column " << col << " persum " << persum << MsgDispatch;
453                    }
454                }
455
456                above = false;
457                if (persum >= fMinpersum) {
458                    above = true;
459                    if (fDebug) {
460                        Msg(EsafMsg::Debug) << "GTU " << gtu << " Row " << row << " Column " << col << MsgDispatch;
461                        Msg(EsafMsg::Debug) << "persum " << persum << " >= fMinpersum " << fMinpersum << MsgDispatch;
462                    }
463                    if (persum > maxpsum) {
464                        maxpsum = persum;
465                        gtutrig = gtu;
466                    }
467                    // a trigger will be generated so flag these tracks
468                    for (Int_t i=0; i<fFrmsinpsum; i++) {
469                        Int_t first_element = fTrks[gtu-i][row][col];
470                        Int_t num_elements = fPers[gtu-i][row][col];
471                        for (Int_t j=0; j<num_elements; j++) {
472                            fTracks[first_element+j].SetTriggerID(ntrigs);
473                            if (fDebug)
474                                Msg(EsafMsg::Debug) << "Flagging track number " << (first_element+j) << MsgDispatch;
475                        }
476                    }
477                }
478
479                if (!above && lastabove) {
480                    ntrigs++;
481                    SetGtuTrigger(gtutrig);
482                    maxpsum = 0;
483                    gtutrig = 0;
484                    if (fDebug) {
485                        Msg(EsafMsg::Debug) << "Track has ended, GENERATING TRIGGER" << MsgDispatch;
486                    }
487                }
488
489                lastabove = above;
490            } // End of gtu loop
491        } // End of col loop
492    } // End of row loop
493
494    return ntrigs;
495}
496
497//_____________________________________________________________________________________
498void LblChipTrackingTrgEngine::FillELblTrackTrigger( Int_t mc_id, Int_t ec_id) {
499
500    // Fill the ELblTrackTrigger branch of the ROOT output with all the tracks
501    vector<LblTrackSegment>::iterator it;
502    for( it=fTracks.begin(); it!=fTracks.end(); it++) {
503        if ( EEvent::GetCurrent() ) {
504            ELblTrackTriggerDataAdder a( &(*it),mc_id,ec_id );
505            EEvent::GetCurrent()->Fill(a);
506        }
507    }
508}
Note: See TracBrowser for help on using the repository browser.