Universita' di GenovaINFN Sezione di Genova  
AIRWATCH / EUSO Genova

LblChipTrackingTrgEngine - source file

// ESAF : Euso Simulation and Analysis Framework
// $Id: LblChipTrackingTrgEngine.cc,v 1.3 2005/10/10 16:28:03 ejudd Exp $

#include "LblChipTrackingTrgEngine.hh"
#include "MacroCellData.hh"
#include "MacroCell.hh"
#include "ElementaryCell.hh"
#include "FrontEndChip.hh"
#include "BoolAlgebra.hh"
#include "LblTrackSegment.hh"
#include "EEvent.hh"
#include "ELblTrackTriggerDataAdder.hh"

ClassImp(LblChipTrackingTrgEngine)

//______________________________________________________________________________
 LblChipTrackingTrgEngine::LblChipTrackingTrgEngine() : 
TriggerEngine(string("LblChipTrackingTrgEngine"), kLblChipTrackingTrigger ) {

    // constructor

    // Read all the control flags from the LblChipTrackingTrgEngine.cfg file

    // First the flag to control which pxiels get used to make track segments
    fPxlthr = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Pxlthr");

    // Now the flags to control how track segments are found
    fColumn = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Column");
    fColthr = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Colthr");
    fFrmsinseg = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Frmsinseg");
    fSegs = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Segs");
    fTrksumthr = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Trksumthr");
    fWantsumsub = (bool)Conf()->GetNum("LblChipTrackingTrgEngine.Wantsumsub");

    // Next the flags to control how segments are joined to form tracks
    // and tracks are tested for persistency
    fFrmsinpsum = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Frmsinpsum");
    fMinpersum = (Int_t)Conf()->GetNum("LblChipTrackingTrgEngine.Minpersum");

    // Finally the flags to control the level of output messages
    fVerbose = (bool)Conf()->GetNum("LblChipTrackingTrgEngine.Verbose");
    fDebug = (bool)Conf()->GetNum("LblChipTrackingTrgEngine.Debug");

    if (fVerbose) {
	Msg(EsafMsg::Info) << "fPxlthr = " << fPxlthr << MsgDispatch;

	Msg(EsafMsg::Info) << "fColumn = " << fColumn << MsgDispatch;
	Msg(EsafMsg::Info) << "fColthr = " << fColthr << MsgDispatch;
	Msg(EsafMsg::Info) << "fFrmsinseg = " << fFrmsinseg << MsgDispatch;
	Msg(EsafMsg::Info) << "fSegs = " << fSegs << MsgDispatch;
	Msg(EsafMsg::Info) << "fTrksumthr = " << fTrksumthr << MsgDispatch;
	Msg(EsafMsg::Info) << "fWantsumsub = " << fWantsumsub << MsgDispatch;

	Msg(EsafMsg::Info) << "fFrmsinpsum = " << fFrmsinpsum << MsgDispatch;
	Msg(EsafMsg::Info) << "fMinpersum = " << fMinpersum << MsgDispatch;

	Msg(EsafMsg::Info) << "fDebug = " << fDebug << MsgDispatch;
    }

}


//______________________________________________________________________________
 LblChipTrackingTrgEngine::~LblChipTrackingTrgEngine() {
    
    // dtor

    // FIXME - hardwired array sizes are BAD
    Clear(199,12);

}


//______________________________________________________________________________
 void LblChipTrackingTrgEngine::Simulate( MacroCellData* pData) {    

    Int_t tot_tracks;
    Int_t ec_tracks;
    Int_t tot_trigs;
    Int_t ec_trigs;

    // simulate chip tracking algorithm
    if (fVerbose)
	Msg(EsafMsg::Info) << "LBL CHIP TRACKING STARTED MACROCELL " << pData->Cell()->Id() << MsgDispatch;

    // Get maps of ChipGtuData for this MacroCell
    // m1 is a reference (i.e. a short alias) for pData->fChipData2
    // fChipData2 is a map: first element is the Id of this Front End
    //                      second element is a pointer to a second map
    map<Int_t, map<Int_t,ChipGtuData*>* > &m1 = pData->fChipData2;
    map<Int_t, map<Int_t,ChipGtuData*>* >::const_iterator it1;
    Int_t ncells = m1.size();
    if (fVerbose)
	Msg(EsafMsg::Info) << "MacroCell contains " << ncells << " Elementary Cells" << MsgDispatch;

    // Quit if there is no data in this map element
    if (!ncells) return;

    // Now do the main loop over the Elementary Cells, find segments and tracks
    // then look for persistency.
    tot_tracks = 0;
    tot_trigs = 0;
    for (it1=m1.begin(); it1!=m1.end(); it1++) {

	if (fVerbose)
	    Msg(EsafMsg::Info) << "Procesing Elementary Cell " << it1->first << MsgDispatch;

	// Get maps of ChipGtuData for this specific Elementary Cell
	// m2 is a reference to the map object that it1->second points to
	// This map has 2 elements: first element is GTU number
	//                          second element is a pointer to ChipGtuData for this GTU
	map<Int_t,ChipGtuData*> &m2 = *(it1->second);
	map<Int_t,ChipGtuData*>::const_iterator it2;

	// Check to see if this Elementary Cell has enough GTUs
	it2 = m2.begin();
	Int_t first_gtu = it2->first;
	it2 = m2.end();
	it2--;
	Int_t last_gtu = it2->first;
	Int_t gtot = last_gtu - first_gtu + 1;
	if (fVerbose) Msg(EsafMsg::Info) << "First GTU = " << first_gtu << " Last GTU = " << last_gtu << " Elementary Cell contains " << gtot << " GTUs" << MsgDispatch;
	if (last_gtu >= 200) {
	    if (fVerbose) Msg(EsafMsg::Info) << "Last GTU is too late" << MsgDispatch;
	    continue;
	}
	if (gtot < (fSegs*fFrmsinseg)) {
	    if (fVerbose) Msg(EsafMsg::Info) << "Not enough GTUs in this Elementary Cell" << MsgDispatch;
	    continue;
	}

	// Get the number of channels in this chip from the first GTU
	it2 = m2.begin();
	Int_t chan = (it2->second)->FrontEnd()->Channels();
	Int_t NPIX = (it2->second)->FrontEnd()->NumSide();
	if (fDebug)
	    Msg(EsafMsg::Debug) << "Number of channels is " << chan << " Number per side is " << NPIX << MsgDispatch;

	// Zero the arrays
	Clear(last_gtu,NPIX);

	// Find the tracks and segments
	ec_tracks = Find_Tracks(m2,last_gtu,chan,NPIX);

	// Check to see if there is a good track
	ec_trigs = Check_Persistency(last_gtu,NPIX);

	if (fVerbose)
	    Msg(EsafMsg::Info) << "EC " << it1->first << " generated " << ec_trigs << " triggers from " << ec_tracks << " tracks" << MsgDispatch;

	tot_tracks += ec_tracks;
	tot_trigs += ec_trigs;

	// Save the tracks in the root output
	FillELblTrackTrigger(pData->Cell()->Id(),it1->first);

    } // End of it1 loop over Elementary Cells

    // Print total track count
    if (fVerbose)
	Msg(EsafMsg::Info) << tot_tracks << " tracks were found in cell " << pData->Cell()->Id() << MsgDispatch;

    // Set trigger flags and save count if any triggers were issued
    if (tot_trigs > 0) {
	SetTrigger(tot_trigs);
	Msg(EsafMsg::Info) << tot_trigs << " triggers have been generated in cell " << pData->Cell()->Id() << MsgDispatch;
    }

}

//______________________________________________________________________________
 void LblChipTrackingTrgEngine::Clear(Int_t last_gtu,Int_t NPIX) {
    
    // Clear the fRawSums, fSums, fPers and fTrks arrays
    for (Int_t i=0; i<(last_gtu+1); i++) {
	for (Int_t j=0; j<NPIX; j++) {
	    for (Int_t k=0; k<NPIX; k++) {
		fRawSums[i][j][k] = 0;
		fSums[i][j][k] = 0;
		fPers[i][j][k] = 0;
		fTrks[i][j][k] = -1;
	    }
	}
    }

    // Clear the track vectors
    fTracks.clear();
}


//______________________________________________________________________________
 Int_t LblChipTrackingTrgEngine::Find_Tracks(map<Int_t,ChipGtuData*> &m2,Int_t last_gtu,Int_t chan,Int_t NPIX) {

    Int_t gt;
    Int_t last_track_gtu;
    Int_t row;
    Int_t col;
    Int_t colsum;
    Int_t sub_row;
    Int_t sub_col;
    Int_t sub_chan;
    ChipGtuData* gtus;
    Int_t seg1sum;
    Int_t seg2sum;
    Int_t seg3sum;
    Int_t trksum;
    Int_t row_l;
    Int_t row_h;
    Int_t col_l;
    Int_t col_h;
    Int_t sub_row2;
    Int_t sub_col2;
    Int_t ntracks;

    if (fDebug) {
	Msg(EsafMsg::Debug) << "Number of channels is " << chan << " Number per side is " << NPIX << MsgDispatch;
	Msg(EsafMsg::Debug) << "last_gtu = " << last_gtu << MsgDispatch;
    }

    // Create the iterator for this map
    // it2 is used to loop over all GTUs
    map<Int_t,ChipGtuData*>::const_iterator it2;

    // Main loop over GTUs
    ntracks = 0;
    for (it2=m2.begin(); it2!=m2.end(); it2++) {

	// Get the data for this specific GTU
	gt = it2->first;
	gtus = it2->second;
	if (fDebug)
	    Msg(EsafMsg::Debug) << "Building Segments starting at GTU Number " << gt << MsgDispatch;

	// Add this data to the sums for each pixel
	for (Int_t nchan=0; nchan<chan; nchan++) {

	    row = gtus->FrontEnd()->Row(nchan);
	    col = gtus->FrontEnd()->Column(nchan);

	    // Now calculate the sums
	    switch (fColumn) {
	    case 2:
		colsum = 0;
		for (Int_t i=0; i<2; i++) {
		    sub_row = row+i;
		    for (Int_t j=0; j<2; j++) {
			sub_col = col+j;
			if ((sub_row<NPIX)&&(sub_col<NPIX)) {
			    sub_chan = gtus->FrontEnd()->ChanRowCol(sub_row,sub_col);
			    if (gtus->GetCounter(sub_chan) > fPxlthr) {
				colsum += gtus->GetCounter(sub_chan);
			    }
			}
		    }
		}
		break;
	    case 3:
		colsum = 0;
		for (Int_t i=0; i<3; i++) {
		    sub_row = row+i;
		    for (Int_t j=0; j<3; j++) {
			sub_col = col+j;
			if ((sub_row<NPIX)&&(sub_col<NPIX)) {
			    sub_chan = gtus->FrontEnd()->ChanRowCol(sub_row,sub_col);
			    if (gtus->GetCounter(sub_chan) > fPxlthr) {
				colsum += gtus->GetCounter(sub_chan);
			    }
			}
		    }
		}
		break;
	    default:
		colsum = 0;
		// Get this pixel's counts, but only if # counts > threshold
		if (gtus->GetCounter(nchan) > fPxlthr)
		    colsum += gtus->GetCounter(nchan);
		break;
	    } // End of switch (fColumn)

	    // Loop over active segments and add in this column
	    Int_t num_active = ((gt+1)<fFrmsinseg) ? (gt+1) : fFrmsinseg;
	    for (Int_t active=0; active<num_active; active++) {
		fRawSums[gt-active][row][col] += (colsum>fColthr) ? colsum : 0;
	    }

	} // End of nchan loop over Channels (i.e. Pixels)

    } // End of loop over map m2. At this point we have the sums for each pixel for all GTUs.

    // Now either subtract off the background counts, if the fWantsumsub flag is true,
    // or just copy fRawSums to fSums if fWantsumsub is false.
    for (gt=0; gt<(last_gtu+1); gt++) {
	for (row=0; row<NPIX; row++) {
	    for (col=0; col<NPIX; col++) {
		if (fWantsumsub) {
		    sub_row = (row + 4)%NPIX;
		    sub_col = (col + 4)%NPIX;
		    Int_t sumsub = fRawSums[gt][row][col] - fRawSums[gt][sub_row][sub_col];
		    fSums[gt][row][col] = (sumsub>0) ? sumsub : 0;
		}
		else {
		    fSums[gt][row][col] = fRawSums[gt][row][col];
		}
	    }
	}
    }

    // Now work with partial sums to find tracksums
    last_track_gtu = (last_gtu+1) - (fSegs*fFrmsinseg);
    for (gt=0; gt<(last_track_gtu+1); gt++) {
	for (row=0; row<NPIX; row++) {
	    for (col=0; col<NPIX; col++) {

		seg1sum = fSums[gt][row][col];

		// Single segment tracks
		if (fSegs == 1) {
		    trksum = seg1sum;
		    if (trksum > fTrksumthr) {
			fPers[gt][row][col] += 1;
			if (fTrks[gt][row][col] == -1) fTrks[gt][row][col] = ntracks;
			ntracks++;
			// Save the track information
			fTracks.push_back( LblTrackSegment() );
			LblTrackSegment &trkseg = fTracks.back();
			trkseg.SetGtuStart(gt);
			trkseg.SetRowStart(row);
			trkseg.SetColStart(col);
			trkseg.SetGtuEnd(gt+fFrmsinseg-1);
			trkseg.SetRowEnd(row);
			trkseg.SetColEnd(col);
			trkseg.SetNumSegs(1);
			trkseg.SetSum(trksum);
			if (fDebug) {
			    Msg(EsafMsg::Debug) << "trksum " << trksum << " > thr " << fTrksumthr << " for single segment track, fPers = " << fPers[gt][row][col] << " fTrks = " << fTrks[gt][row][col] << MsgDispatch;
			    Msg(EsafMsg::Debug) << "Single segment start-GTU = " << gt << " start-Row = " << row << " start-Col = " << col << MsgDispatch;
			}
		    }
		}

		// Two-or-more segment tracks
		if (fSegs >= 2) {
		    row_l = (row>0)        ? row-1 : 0;
		    row_h = (row<(NPIX-1)) ? row+1 : (NPIX-1);
		    col_l = (col>0)        ? col-1 : 0;
		    col_h = (col<(NPIX-1)) ? col+1 : (NPIX-1);

		    for (sub_row=row_l; sub_row<=row_h; sub_row++) {
			for (sub_col=col_l; sub_col<=col_h; sub_col++) {
			    
			    seg2sum = fSums[gt+fFrmsinseg][sub_row][sub_col];

			    // Two segment tracks
			    if (fSegs == 2) {
				trksum = seg1sum + seg2sum;
				if (trksum > fTrksumthr) {
				    fPers[gt][row][col] += 1;
				    if (fTrks[gt][row][col] == -1) fTrks[gt][row][col] = ntracks;
				    ntracks++;
				    // Save the track information
				    fTracks.push_back( LblTrackSegment() );
				    LblTrackSegment &trkseg = fTracks.back();
				    trkseg.SetGtuStart(gt);
				    trkseg.SetRowStart(row);
				    trkseg.SetColStart(col);
				    trkseg.SetGtuEnd(gt+(2*fFrmsinseg)-1);
				    trkseg.SetRowEnd(sub_row);
				    trkseg.SetColEnd(sub_col);
				    trkseg.SetNumSegs(2);
				    trkseg.SetSum(trksum);
				    if (fDebug) {
					Msg(EsafMsg::Debug) << "trksum " << trksum << " > thr " << fTrksumthr << " for two segment track, fPers = " << fPers[gt][row][col] << " fTrks = " << fTrks[gt][row][col] << MsgDispatch;
					Msg(EsafMsg::Debug) << "First  segment start-GTU = " << gt << " start-Row = " << row << " start-Col = " << col << MsgDispatch;
					Msg(EsafMsg::Debug) << "Second segment start-GTU = " << gt+fFrmsinseg << " start-Row = " << sub_row << " start-Col = " << sub_col << MsgDispatch;
				    }
				}
			    }

			    // Three segment tracks
			    if (fSegs == 3) {
				// Find sub_row2 and sub_col2 of 3rd segment by straight
				// line extrapolation from row/col and sub_row/sub_col
				sub_row2 = sub_row + (sub_row - row);
				if ((sub_row2 < 0) || (sub_row2 >= NPIX)) continue;
				sub_col2 = sub_col + (sub_col - col);
				if ((sub_col2 < 0) || (sub_col2 >= NPIX)) continue;

				seg3sum = fSums[gt+(2*fFrmsinseg)][sub_row2][sub_col2];

				trksum = seg1sum + seg2sum + seg3sum;
				if (trksum > fTrksumthr) {
				    fPers[gt][row][col] += 1;
				    if (fTrks[gt][row][col] == -1) fTrks[gt][row][col] = ntracks;
				    ntracks++;
				    // Save the track information
				    fTracks.push_back( LblTrackSegment() );
				    LblTrackSegment &trkseg = fTracks.back();
				    trkseg.SetGtuStart(gt);
				    trkseg.SetRowStart(row);
				    trkseg.SetColStart(col);
				    trkseg.SetGtuEnd(gt+(3*fFrmsinseg)-1);
				    trkseg.SetRowEnd(sub_row2);
				    trkseg.SetColEnd(sub_col2);
				    trkseg.SetNumSegs(3);
				    trkseg.SetSum(trksum);
				    if (fDebug) {
					Msg(EsafMsg::Debug) << "trksum " << trksum << " > thr " << fTrksumthr << " for three segment track, fPers = " << fPers[gt][row][col] << " fTrks = " << fTrks[gt][row][col] << MsgDispatch;
					Msg(EsafMsg::Debug) << "First  segment start-GTU = " << gt << " start-Row = " << row << " start-Col = " << col << MsgDispatch;
					Msg(EsafMsg::Debug) << "Second segment start-GTU = " << gt+fFrmsinseg << " start-Row = " << sub_row << " start-Col = " << sub_col << MsgDispatch;
					Msg(EsafMsg::Debug) << "Third  segment start-GTU = " << (gt+(2*fFrmsinseg)) << " start-Row = " << sub_row2 << " start-Col = " << sub_col2 << MsgDispatch;
				    }
				}
			    } // End of fSegs == 3 if statement

			} // End of sub_col loop
		    } // End of sub_row loop
		} // End of fSegs >= 2 if statement
	    } // End of col loop
	} // End of row loop
    } // End of gt loop over all GTUs 
    
    return ntracks;

}


//______________________________________________________________________________
 Int_t LblChipTrackingTrgEngine::Check_Persistency(Int_t last_gtu,Int_t NPIX) {
    Int_t persum;
    Int_t maxpsum;
    Int_t gtutrig;
    bool above;
    bool lastabove;
    Int_t ntrigs;
    
    if (fDebug)
	Msg(EsafMsg::Debug) << "last_gtu = " << last_gtu << " NPIX = " << NPIX << MsgDispatch;

    ntrigs = 0;
    for (Int_t row=0; row<NPIX; row++) {
	for (Int_t col=0; col<NPIX; col++) {
	    maxpsum = 0;
	    gtutrig = 0;
	    lastabove = false;

	    for (Int_t gtu=(fFrmsinpsum-1); gtu<(last_gtu+1); gtu++) {

		persum = 0;
		for (Int_t i=0; i<fFrmsinpsum; i++) {
		    if (fPers[gtu-i][row][col] > 0) {
			persum++;
			if (fDebug) Msg(EsafMsg::Debug) << "GTU " << gtu << " Row " << row << " Column " << col << " persum " << persum << MsgDispatch;
		    }
		}

		above = false;
		if (persum >= fMinpersum) {
		    above = true;
		    if (fDebug) {
			Msg(EsafMsg::Debug) << "GTU " << gtu << " Row " << row << " Column " << col << MsgDispatch;
			Msg(EsafMsg::Debug) << "persum " << persum << " >= fMinpersum " << fMinpersum << MsgDispatch;
		    }
		    if (persum > maxpsum) {
			maxpsum = persum;
			gtutrig = gtu;
		    }
		    // a trigger will be generated so flag these tracks
		    for (Int_t i=0; i<fFrmsinpsum; i++) {
			Int_t first_element = fTrks[gtu-i][row][col];
			Int_t num_elements = fPers[gtu-i][row][col];
			for (Int_t j=0; j<num_elements; j++) {
			    fTracks[first_element+j].SetTriggerID(ntrigs);
			    if (fDebug)
				Msg(EsafMsg::Debug) << "Flagging track number " << (first_element+j) << MsgDispatch;
			}
		    }
		}

		if (!above && lastabove) {
		    ntrigs++;
		    SetGtuTrigger(gtutrig);
		    maxpsum = 0;
		    gtutrig = 0;
		    if (fDebug) {
			Msg(EsafMsg::Debug) << "Track has ended, GENERATING TRIGGER" << MsgDispatch;
		    }
		}

		lastabove = above;
	    } // End of gtu loop
	} // End of col loop
    } // End of row loop

    return ntrigs;
}

//_____________________________________________________________________________________
 void LblChipTrackingTrgEngine::FillELblTrackTrigger( Int_t mc_id, Int_t ec_id) {

    // Fill the ELblTrackTrigger branch of the ROOT output with all the tracks
    vector<LblTrackSegment>::iterator it;
    for( it=fTracks.begin(); it!=fTracks.end(); it++) {
	if ( EEvent::GetCurrent() ) {
	    ELblTrackTriggerDataAdder a( &(*it),mc_id,ec_id );
	    EEvent::GetCurrent()->Fill(a);
	}
    }
}
About Us | EUSO Official Website | Web pages created by Roberto Pesce and Alessandro Thea - Last Update Wed Nov 16 16:57:39 2005 Wed Nov 16 16:29:22 2005