Newer
Older
\ingroup Detectors
\brief Class for an Aerogel detector consisting of pairs of PMT's
attached to a diffuser box
#include "THcAerogel.h"
#include "THcHodoscope.h"
#include "TClonesArray.h"
#include "THcSignalHit.h"
#include "THaEvData.h"
#include "THaDetMap.h"
#include "THcDetectorMap.h"
#include "THcGlobals.h"
#include "THaCutList.h"
#include "THcParmList.h"
#include "THaApparatus.h"
#include "VarDef.h"
#include "VarType.h"
#include "THaTrack.h"
#include "TClonesArray.h"
#include "TMath.h"
#include "THaTrackProj.h"
#include "THcRawAdcHit.h"
#include "THcHallCSpectrometer.h"
#include <cstring>
#include <cstdio>
#include <cstdlib>
#include <iostream>
using namespace std;
//_____________________________________________________________________________
THcAerogel::THcAerogel( const char* name, const char* description,
THaApparatus* apparatus ) :
fAdcPosTimeWindowMin(0), fAdcPosTimeWindowMax(0), fAdcNegTimeWindowMin(0),
fAdcNegTimeWindowMax(0), fRegionValue(0), fPosGain(0), fNegGain(0),
frPosAdcPedRaw(0), frPosAdcPulseIntRaw(0), frPosAdcPulseAmpRaw(0),
frPosAdcPulseTimeRaw(0), frPosAdcPed(0), frPosAdcPulseInt(0),
frPosAdcPulseAmp(0), frPosAdcPulseTime(0), frNegAdcPedRaw(0),
frNegAdcPulseIntRaw(0), frNegAdcPulseAmpRaw(0), frNegAdcPulseTimeRaw(0),
frNegAdcPed(0), frNegAdcPulseInt(0), frNegAdcPulseAmp(0),
frNegAdcPulseTime(0), fPosAdcErrorFlag(0),
fNegAdcErrorFlag(0), fPosPedSum(0), fPosPedSum2(0), fPosPedLimit(0),
fPosPedCount(0), fNegPedSum(0), fNegPedSum2(0), fNegPedLimit(0), fNegPedCount(0),
fA_Pos(0), fA_Neg(0), fA_Pos_p(0), fA_Neg_p(0), fT_Pos(0), fT_Neg(0),
fPosPed(0), fPosSig(0), fPosThresh(0), fNegPed(0), fNegSig(0),
fNegThresh(0), fPosPedMean(0), fNegPedMean(0),
fPosTDCHits(0), fNegTDCHits(0), fPosADCHits(0), fNegADCHits(0)
{
}
//_____________________________________________________________________________
THcAerogel::THcAerogel( ) :
fAdcPosTimeWindowMin(0), fAdcPosTimeWindowMax(0), fAdcNegTimeWindowMin(0),
fAdcNegTimeWindowMax(0), fRegionValue(0), fPosGain(0), fNegGain(0),
frPosAdcPedRaw(0), frPosAdcPulseIntRaw(0), frPosAdcPulseAmpRaw(0),
frPosAdcPulseTimeRaw(0), frPosAdcPed(0), frPosAdcPulseInt(0),
frPosAdcPulseAmp(0), frPosAdcPulseTime(0), frNegAdcPedRaw(0),
frNegAdcPulseIntRaw(0), frNegAdcPulseAmpRaw(0), frNegAdcPulseTimeRaw(0),
frNegAdcPed(0), frNegAdcPulseInt(0), frNegAdcPulseAmp(0),
frNegAdcPulseTime(0), fPosAdcErrorFlag(0),
fNegAdcErrorFlag(0), fPosPedSum(0), fPosPedSum2(0), fPosPedLimit(0),
fPosPedCount(0), fNegPedSum(0), fNegPedSum2(0), fNegPedLimit(0), fNegPedCount(0),
fA_Pos(0), fA_Neg(0), fA_Pos_p(0), fA_Neg_p(0), fT_Pos(0), fT_Neg(0),
fPosPed(0), fPosSig(0), fPosThresh(0), fNegPed(0), fNegSig(0),
fNegThresh(0), fPosPedMean(0), fNegPedMean(0),
fPosTDCHits(0), fNegTDCHits(0), fPosADCHits(0), fNegADCHits(0)
//_____________________________________________________________________________
THcAerogel::~THcAerogel()
{
// Destructor
DeleteArrays();
}
//_____________________________________________________________________________
void THcAerogel::DeleteArrays()
{
// Delete all dynamically allocated memory
delete frPosAdcPedRaw; frPosAdcPedRaw = NULL;
delete frPosAdcPulseIntRaw; frPosAdcPulseIntRaw = NULL;
delete frPosAdcPulseAmpRaw; frPosAdcPulseAmpRaw = NULL;
delete frPosAdcPulseTimeRaw; frPosAdcPulseTimeRaw = NULL;
delete frPosAdcPed; frPosAdcPed = NULL;
delete frPosAdcPulseInt; frPosAdcPulseInt = NULL;
delete frPosAdcPulseAmp; frPosAdcPulseAmp = NULL;
delete frPosAdcPulseTime; frPosAdcPulseTime = NULL;
delete frNegAdcPedRaw; frNegAdcPedRaw = NULL;
delete frNegAdcPulseIntRaw; frNegAdcPulseIntRaw = NULL;
delete frNegAdcPulseAmpRaw; frNegAdcPulseAmpRaw = NULL;
delete frNegAdcPulseTimeRaw; frNegAdcPulseTimeRaw = NULL;
delete frNegAdcPed; frNegAdcPed = NULL;
delete frNegAdcPulseInt; frNegAdcPulseInt = NULL;
delete frNegAdcPulseAmp; frNegAdcPulseAmp = NULL;
delete frNegAdcPulseTime; frNegAdcPulseTime = NULL;
delete fPosAdcErrorFlag; fPosAdcErrorFlag = NULL;
delete fNegAdcErrorFlag; fNegAdcErrorFlag = NULL;
delete [] fRegionValue; fRegionValue = 0;
delete [] fAdcPosTimeWindowMin; fAdcPosTimeWindowMin = 0;
delete [] fAdcPosTimeWindowMax; fAdcPosTimeWindowMax = 0;
delete [] fAdcNegTimeWindowMin; fAdcNegTimeWindowMin = 0;
delete [] fAdcNegTimeWindowMax; fAdcNegTimeWindowMax = 0;
// 6 GeV variables
delete fPosTDCHits; fPosTDCHits = NULL;
delete fNegTDCHits; fNegTDCHits = NULL;
delete fPosADCHits; fPosADCHits = NULL;
delete fNegADCHits; fNegADCHits = NULL;
delete [] fPosGain; fPosGain = NULL;
delete [] fNegGain; fNegGain = NULL;
delete [] fA_Pos; fA_Pos = NULL;
delete [] fA_Neg; fA_Neg = NULL;
delete [] fA_Pos_p; fA_Pos_p = NULL;
delete [] fA_Neg_p; fA_Neg_p = NULL;
delete [] fT_Pos; fT_Pos = NULL;
delete [] fT_Neg; fT_Neg = NULL;
if (fSixGevData)
DeletePedestalArrays();
}
//_____________________________________________________________________________
void THcAerogel::DeletePedestalArrays()
{
// Delete all dynamically allocated memory for pedestal processing
delete [] fPosPedLimit; fPosPedLimit = NULL;
delete [] fNegPedLimit; fNegPedLimit = NULL;
delete [] fPosPedMean; fPosPedMean = NULL;
delete [] fNegPedMean; fNegPedMean = NULL;
delete [] fPosPedSum; fPosPedSum = NULL;
delete [] fPosPedSum2; fPosPedSum2 = NULL;
delete [] fPosPedCount; fPosPedCount = NULL;
delete [] fNegPedSum; fNegPedSum = NULL;
delete [] fNegPedSum2; fNegPedSum2 = NULL;
delete [] fNegPedCount; fNegPedCount = NULL;
delete [] fPosPed; fPosPed = NULL;
delete [] fPosSig; fPosSig = NULL;
delete [] fPosThresh; fPosThresh = NULL;
delete [] fNegPed; fNegPed = NULL;
delete [] fNegSig; fNegSig = NULL;
delete [] fNegThresh; fNegThresh = NULL;
//_____________________________________________________________________________
THaAnalysisObject::EStatus THcAerogel::Init( const TDatime& date )
{
// cout << "THcAerogel::Init for: " << GetName() << endl;
char EngineDID[] = "xAERO";
EngineDID[0] = toupper(GetApparatus()->GetName()[0]);
if( gHcDetectorMap->FillMap(fDetMap, EngineDID) < 0 ) {
static const char* const here = "Init()";
_det_logger->error( "{} Error filling detectormap for {}.", "THcAerogel::Init", EngineDID );
return kInitError;
EStatus status;
if( (status = THaNonTrackingDetector::Init( date )) )
return fStatus=status;
// Should probably put this in ReadDatabase as we will know the
// maximum number of hits after setting up the detector map
InitHitList(fDetMap, "THcAerogelHit", fDetMap->GetTotNumChan()+1,
0, fADC_RefTimeCut);
THcHallCSpectrometer *app=dynamic_cast<THcHallCSpectrometer*>(GetApparatus());
if( !app ||
!(fglHod = dynamic_cast<THcHodoscope*>(app->GetDetector("hod"))) ) {
static const char* const here = "ReadDatabase()";
_det_logger->warn("{} Hodoscope \"{}\" not found. ",Here(here),"hod");
fPresentP = 0;
THaVar* vpresent = gHaVars->Find(Form("%s.present",GetApparatus()->GetName()));
if(vpresent) {
fPresentP = (Bool_t *) vpresent->GetValuePointer();
}
return fStatus = kOK;
}
//_____________________________________________________________________________
Int_t THcAerogel::ReadDatabase( const TDatime& date )
{
// This function is called by THaDetectorBase::Init() once at the beginning
// of the analysis.
// cout << "THcAerogel::ReadDatabase for: " << GetName() << endl;
char prefix[2];
prefix[0]=tolower(GetApparatus()->GetName()[0]);
prefix[1]='\0';
CreateMissReportParms(Form("%saero",prefix));
fNRegions = 1; // Default if not set in parameter file
{"aero_num_pairs", &fNelem, kInt},
gHcParms->LoadParmValues((DBRequest*)&listextra, prefix);
Bool_t optional = true;
_det_logger->info("Created aerogel detector {}.{} with {} PMT pairs.", GetApparatus()->GetName(), GetName(), fNelem);
//cout << "Created aerogel detector " << GetApparatus()->GetName() << "."
// << GetName() << " with " << fNelem << " PMT pairs" << endl;
DeleteArrays(); // avoid memory leak when reinitializing
fPosGain = new Double_t[fNelem];
fNegGain = new Double_t[fNelem];
// 6 GeV variables
fTdcOffset = 0; // Offset to make reference time subtracted times positve
fA_Pos = new Float_t[fNelem];
fA_Neg = new Float_t[fNelem];
fA_Pos_p = new Float_t[fNelem];
fA_Neg_p = new Float_t[fNelem];
fT_Pos = new Float_t[fNelem];
fT_Neg = new Float_t[fNelem];
// Normal constructor with name and description
frPosAdcPedRaw = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frPosAdcPulseIntRaw = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frPosAdcPulseAmpRaw = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frPosAdcPulseTimeRaw = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frPosAdcPed = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frPosAdcPulseInt = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frPosAdcPulseAmp = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frPosAdcPulseTime = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frNegAdcPedRaw = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frNegAdcPulseIntRaw = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frNegAdcPulseAmpRaw = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frNegAdcPulseTimeRaw = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frNegAdcPed = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frNegAdcPulseInt = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frNegAdcPulseAmp = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
frNegAdcPulseTime = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
fPosAdcErrorFlag = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
fNegAdcErrorFlag = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
fNumPosAdcHits.assign(fNelem, 0);
fNumGoodPosAdcHits.assign(fNelem, 0);
fNumNegAdcHits.assign(fNelem, 0);
fNumGoodNegAdcHits.assign(fNelem, 0);
fNumTracksMatched.assign(fNelem, 0);
fNumTracksFired.assign(fNelem, 0);
fPosNpe.assign(fNelem, 0.0);
fNegNpe.assign(fNelem, 0.0);
fGoodPosAdcPed.assign(fNelem, 0.0);
fGoodPosAdcMult.assign(fNelem, 0.0);
fGoodPosAdcPulseInt.assign(fNelem, 0.0);
fGoodPosAdcPulseIntRaw.assign(fNelem, 0.0);
fGoodPosAdcPulseAmp.assign(fNelem, 0.0);
fGoodPosAdcPulseTime.assign(fNelem, 0.0);
fGoodPosAdcTdcDiffTime.assign(fNelem, 0.0);
fGoodNegAdcPed.assign(fNelem, 0.0);
fGoodNegAdcMult.assign(fNelem, 0.0);
fGoodNegAdcPulseInt.assign(fNelem, 0.0);
fGoodNegAdcPulseIntRaw.assign(fNelem, 0.0);
fGoodNegAdcPulseAmp.assign(fNelem, 0.0);
fGoodNegAdcPulseTime.assign(fNelem, 0.0);
fGoodNegAdcTdcDiffTime.assign(fNelem, 0.0);
// 6 GeV variables
fPosTDCHits = new TClonesArray("THcSignalHit", fNelem*16);
fNegTDCHits = new TClonesArray("THcSignalHit", fNelem*16);
fPosADCHits = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
fNegADCHits = new TClonesArray("THcSignalHit", fNelem*MaxNumAdcPulse);
fPosNpeSixGev.assign(fNelem, 0.0);
fNegNpeSixGev.assign(fNelem, 0.0);
// Region parameters
fRegionsValueMax = fNRegions * 8;
fRegionValue = new Double_t[fRegionsValueMax];
fAdcPosTimeWindowMin = new Double_t [fNelem];
fAdcPosTimeWindowMax = new Double_t [fNelem];
fAdcNegTimeWindowMin = new Double_t [fNelem];
fAdcNegTimeWindowMax = new Double_t [fNelem];
{"aero_num_regions", &fNRegions, kInt},
{"aero_red_chi2_min", &fRedChi2Min, kDouble},
{"aero_red_chi2_max", &fRedChi2Max, kDouble},
{"aero_beta_min", &fBetaMin, kDouble},
{"aero_beta_max", &fBetaMax, kDouble},
{"aero_enorm_min", &fENormMin, kDouble},
{"aero_enorm_max", &fENormMax, kDouble},
{"aero_dp_min", &fDpMin, kDouble},
{"aero_dp_max", &fDpMax, kDouble},
{"aero_diff_box_zpos", &fDiffBoxZPos, kDouble},
{"aero_npe_thresh", &fNpeThresh, kDouble},
//// {"aero_adcTimeWindowMin", &fAdcTimeWindowMin, kDouble},
//// {"aero_adcTimeWindowMax", &fAdcTimeWindowMax, kDouble},
{"aero_adcPosTimeWindowMin", fAdcPosTimeWindowMin, kDouble, static_cast<UInt_t>(fNelem), 1},
{"aero_adcPosTimeWindowMax", fAdcPosTimeWindowMax, kDouble, static_cast<UInt_t>(fNelem), 1},
{"aero_adcNegTimeWindowMin", fAdcNegTimeWindowMin, kDouble, static_cast<UInt_t>(fNelem), 1},
{"aero_adcNegTimeWindowMax", fAdcNegTimeWindowMax, kDouble, static_cast<UInt_t>(fNelem), 1},
{"aero_adc_tdc_offset", &fAdcTdcOffset, kDouble, 0, 1},
{"aero_debug_adc", &fDebugAdc, kInt, 0, 1},
{"aero_six_gev_data", &fSixGevData, kInt, 0, 1},
{"aero_pos_gain", fPosGain, kDouble, (UInt_t) fNelem},
{"aero_neg_gain", fNegGain, kDouble, (UInt_t) fNelem},
{"aero_tdc_offset", &fTdcOffset, kInt, 0, optional},
{"aero_region", &fRegionValue[0], kDouble, (UInt_t) fRegionsValueMax},
{"aero_adcrefcut", &fADC_RefTimeCut, kInt, 0, 1},
for(Int_t ip=0;ip<fNelem;ip++) {
fAdcPosTimeWindowMin[ip] = -1000.;
fAdcNegTimeWindowMin[ip] = -1000.;
fAdcPosTimeWindowMax[ip] = 1000.;
fAdcNegTimeWindowMax[ip] = 1000.;
}
fSixGevData = 0; // Set 6 GeV data parameter to false unless set in parameter file
fDebugAdc = 0; // Set ADC debug parameter to false unless set in parameter file
fAdcTdcOffset = 0.0;
gHcParms->LoadParmValues((DBRequest*)&list, prefix);
if (fSixGevData) {
// Create arrays to hold pedestal results
InitializePedestals();
DBRequest list2[]={
{"aero_pos_ped_limit", fPosPedLimit, kInt, (UInt_t) fNelem, optional},
{"aero_neg_ped_limit", fNegPedLimit, kInt, (UInt_t) fNelem, optional},
{"aero_pos_ped_mean", fPosPedMean, kDouble, (UInt_t) fNelem, optional},
{"aero_neg_ped_mean", fNegPedMean, kDouble, (UInt_t) fNelem, optional},
{"aero_min_peds", &fMinPeds, kInt, 0, optional},
{0}
};
gHcParms->LoadParmValues((DBRequest*)&list2, prefix);
cout << "6 GeV Data Analysis Flag Set To TRUE" << endl;
}
// cout << "Track Matching Parameters for: " << GetApparatus()->GetName()
// << "." << GetName() << endl;
// for (Int_t iregion = 0; iregion < fNRegions; iregion++) {
// cout << "Region = " << iregion + 1 << endl;
// for (Int_t ivalue = 0; ivalue < 8; ivalue++)
// cout << fRegionValue[GetIndex(iregion, ivalue)] << " ";
// cout << endl;
// }
return kOK;
}
//_____________________________________________________________________________
Int_t THcAerogel::DefineVariables( EMode mode )
{
// Initialize global variables for histogramming and tree
// cout << "THcAerogel::DefineVariables called for: " << GetName() << endl;
if( mode == kDefine && fIsSetup ) return kOK;
fIsSetup = ( mode == kDefine );
// Register variables in global list
// Do we need to put the number of pos/neg TDC/ADC hits into the variables?
// No. They show up in tree as Ndata.H.aero.postdchits for example
if (fDebugAdc) {
RVarDef vars[] = {
{"posGain", "Positive PMT gains", "fPosGain"},
{"negGain", "Negative PMT gains", "fNegGain"},
{"numPosAdcHits", "Number of Positive ADC Hits Per PMT", "fNumPosAdcHits"}, // Aerogel occupancy
{"totNumPosAdcHits", "Total Number of Positive ADC Hits", "fTotNumPosAdcHits"}, // Aerogel multiplicity
{"numNegAdcHits", "Number of Negative ADC Hits Per PMT", "fNumNegAdcHits"}, // Aerogel occupancy
{"totNumNegAdcHits", "Total Number of Negative ADC Hits", "fTotNumNegAdcHits"}, // Aerogel multiplicity
{"totnumAdcHits", "Total Number of ADC Hits Per PMT", "fTotNumAdcHits"}, // Aerogel multiplicity
{"posAdcPedRaw", "Positive Raw ADC pedestals", "frPosAdcPedRaw.THcSignalHit.GetData()"},
{"posAdcPulseIntRaw", "Positive Raw ADC pulse integrals", "frPosAdcPulseIntRaw.THcSignalHit.GetData()"},
{"posAdcPulseAmpRaw", "Positive Raw ADC pulse amplitudes", "frPosAdcPulseAmpRaw.THcSignalHit.GetData()"},
{"posAdcPulseTimeRaw", "Positive Raw ADC pulse times", "frPosAdcPulseTimeRaw.THcSignalHit.GetData()"},
{"posAdcPed", "Positive ADC pedestals", "frPosAdcPed.THcSignalHit.GetData()"},
{"posAdcPulseInt", "Positive ADC pulse integrals", "frPosAdcPulseInt.THcSignalHit.GetData()"},
{"posAdcPulseAmp", "Positive ADC pulse amplitudes", "frPosAdcPulseAmp.THcSignalHit.GetData()"},
{"posAdcPulseTime", "Positive ADC pulse times", "frPosAdcPulseTime.THcSignalHit.GetData()"},
{"negAdcPedRaw", "Negative Raw ADC pedestals", "frNegAdcPedRaw.THcSignalHit.GetData()"},
{"negAdcPulseIntRaw", "Negative Raw ADC pulse integrals", "frNegAdcPulseIntRaw.THcSignalHit.GetData()"},
{"negAdcPulseAmpRaw", "Negative Raw ADC pulse amplitudes", "frNegAdcPulseAmpRaw.THcSignalHit.GetData()"},
{"negAdcPulseTimeRaw", "Negative Raw ADC pulse times", "frNegAdcPulseTimeRaw.THcSignalHit.GetData()"},
{"negAdcPed", "Negative ADC pedestals", "frNegAdcPed.THcSignalHit.GetData()"},
{"negAdcPulseInt", "Negative ADC pulse integrals", "frNegAdcPulseInt.THcSignalHit.GetData()"},
{"negAdcPulseAmp", "Negative ADC pulse amplitudes", "frNegAdcPulseAmp.THcSignalHit.GetData()"},
{"negAdcPulseTime", "Negative ADC pulse times", "frNegAdcPulseTime.THcSignalHit.GetData()"},
{ 0 }
};
DefineVarsFromList( vars, mode);
} //end debug statement
if (fSixGevData) {
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
RVarDef vars[] = {
{"apos", "Positive Raw ADC Amplitudes", "fA_Pos"},
{"aneg", "Negative Raw ADC Amplitudes", "fA_Neg"},
{"apos_p", "Positive Ped-subtracted ADC Amplitudes", "fA_Pos_p"},
{"aneg_p", "Negative Ped-subtracted ADC Amplitudes", "fA_Neg_p"},
{"tpos", "Positive Raw TDC", "fT_Pos"},
{"tneg", "Negative Raw TDC", "fT_Neg"},
{"ntdc_pos_hits", "Number of Positive Tube Hits", "fNTDCPosHits"},
{"ntdc_neg_hits", "Number of Negative Tube Hits", "fNTDCNegHits"},
{"posadchits", "Positive ADC hits", "fPosADCHits.THcSignalHit.GetPaddleNumber()"},
{"negadchits", "Negative ADC hits", "fNegADCHits.THcSignalHit.GetPaddleNumber()"},
{"postdchits", "Positive TDC hits", "fPosTDCHits.THcSignalHit.GetPaddleNumber()"},
{"negtdchits", "Negative TDC hits", "fNegTDCHits.THcSignalHit.GetPaddleNumber()"},
{"nGoodHits", "Total number of good hits", "fNGoodHits"},
{"posNpeSixGev", "Number of Positive PEs", "fPosNpeSixGev"},
{"negNpeSixGev", "Number of Negative PEs", "fNegNpeSixGev"},
{"posNpeSumSixGev", "Total Number of Positive PEs", "fPosNpeSumSixGev"},
{"negNpeSumSixGev", "Total Number of Negative PEs", "fNegNpeSumSixGev"},
{"npeSumSixGev", "Total Number of PEs", "fNpeSumSixGev"},
{ 0 }
};
DefineVarsFromList( vars, mode);
} //end fSixGevData statement
RVarDef vars[] = {
{"posAdcCounter", "Positive ADC counter numbers", "frPosAdcPulseIntRaw.THcSignalHit.GetPaddleNumber()"},
{"negAdcCounter", "Negative ADC counter numbers", "frNegAdcPulseIntRaw.THcSignalHit.GetPaddleNumber()"},
{"posAdcErrorFlag", "Error Flag for When FPGA Fails", "fPosAdcErrorFlag.THcSignalHit.GetData()"},
{"negAdcErrorFlag", "Error Flag for When FPGA Fails", "fNegAdcErrorFlag.THcSignalHit.GetData()"},
{"numGoodPosAdcHits", "Number of Good Positive ADC Hits Per PMT", "fNumGoodPosAdcHits"}, // Aerogel occupancy
{"numGoodNegAdcHits", "Number of Good Negative ADC Hits Per PMT", "fNumGoodNegAdcHits"}, // Aerogel occupancy
{"totNumGoodPosAdcHits", "Total Number of Good Positive ADC Hits", "fTotNumGoodPosAdcHits"}, // Aerogel multiplicity
{"totNumGoodNegAdcHits", "Total Number of Good Negative ADC Hits", "fTotNumGoodNegAdcHits"}, // Aerogel multiplicity
{"totnumGoodAdcHits", "TotalNumber of Good ADC Hits Per PMT", "fTotNumGoodAdcHits"}, // Aerogel multiplicity
{"numTracksMatched", "Number of Tracks Matched Per Region", "fNumTracksMatched"},
{"numTracksFired", "Number of Tracks that Fired Per Region", "fNumTracksFired"},
{"totNumTracksMatched", "Total Number of Tracks Matched Per Region", "fTotNumTracksMatched"},
{"totNumTracksFired", "Total Number of Tracks that Fired", "fTotNumTracksFired"},
{"xAtAero", "Track X at Aero diffusion box", "fXAtAero"},
{"yAtAero", "Track Y at Aero diffusion box", "fYAtAero"},
{"posNpe", "Number of Positive PEs", "fPosNpe"},
{"negNpe", "Number of Negative PEs", "fNegNpe"},
{"posNpeSum", "Total Number of Positive PEs", "fPosNpeSum"},
{"negNpeSum", "Total Number of Negative PEs", "fNegNpeSum"},
{"npeSum", "Total Number of PEs", "fNpeSum"},
{"goodPosAdcPed", "Good Negative ADC pedestals", "fGoodPosAdcPed"},
{"goodPosAdcMult", "Good Positive ADC mult", "fGoodPosAdcMult"},
{"goodPosAdcPulseInt", "Good Negative ADC pulse integrals", "fGoodPosAdcPulseInt"},
{"goodPosAdcPulseIntRaw", "Good Negative ADC raw pulse integrals", "fGoodPosAdcPulseIntRaw"},
{"goodPosAdcPulseAmp", "Good Negative ADC pulse amplitudes", "fGoodPosAdcPulseAmp"},
{"goodPosAdcPulseTime", "Good Negative ADC pulse times", "fGoodPosAdcPulseTime"},
{"goodPosAdcTdcDiffTime", "Good Positive hodo Start - ADC pulse times", "fGoodPosAdcTdcDiffTime"},
{"goodNegAdcPed", "Good Negative ADC pedestals", "fGoodNegAdcPed"},
{"goodNegAdcMult", "Good Negative ADC Mult", "fGoodNegAdcMult"},
{"goodNegAdcPulseInt", "Good Negative ADC pulse integrals", "fGoodNegAdcPulseInt"},
{"goodNegAdcPulseIntRaw", "Good Negative ADC raw pulse integrals", "fGoodNegAdcPulseIntRaw"},
{"goodNegAdcPulseAmp", "Good Negative ADC pulse amplitudes", "fGoodNegAdcPulseAmp"},
{"goodNegAdcPulseTime", "Good Negative ADC pulse times", "fGoodNegAdcPulseTime"},
{"goodNegAdcTdcDiffTime", "Good Negative hodo Start - ADC pulse times", "fGoodNegAdcTdcDiffTime"},
{ 0 }
};
return DefineVarsFromList(vars, mode);
//_____________________________________________________________________________
void THcAerogel::Clear(Option_t* opt)
{
fNhits = 0;
fTotNumAdcHits = 0;
fTotNumGoodAdcHits = 0;
fTotNumPosAdcHits = 0;
fTotNumNegAdcHits = 0;
fTotNumGoodPosAdcHits = 0;
fTotNumGoodNegAdcHits = 0;
fTotNumTracksMatched = 0;
fTotNumTracksFired = 0;
fXAtAero = 0.0;
fYAtAero = 0.0;
fPosNpeSum = 0.0;
fNegNpeSum = 0.0;
frPosAdcPedRaw->Clear();
frPosAdcPulseIntRaw->Clear();
frPosAdcPulseAmpRaw->Clear();
frPosAdcPulseTimeRaw->Clear();
frPosAdcPed->Clear();
frPosAdcPulseInt->Clear();
frPosAdcPulseAmp->Clear();
frNegAdcPedRaw->Clear();
frNegAdcPulseIntRaw->Clear();
frNegAdcPulseAmpRaw->Clear();
frNegAdcPulseTimeRaw->Clear();
frNegAdcPed->Clear();
frNegAdcPulseInt->Clear();
frNegAdcPulseAmp->Clear();
fPosAdcErrorFlag->Clear();
fNegAdcErrorFlag->Clear();
for (UInt_t ielem = 0; ielem < fNumPosAdcHits.size(); ielem++)
fNumPosAdcHits.at(ielem) = 0;
for (UInt_t ielem = 0; ielem < fNumNegAdcHits.size(); ielem++)
fNumNegAdcHits.at(ielem) = 0;
for (UInt_t ielem = 0; ielem < fNumGoodPosAdcHits.size(); ielem++)
fNumGoodPosAdcHits.at(ielem) = 0;
for (UInt_t ielem = 0; ielem < fNumGoodNegAdcHits.size(); ielem++)
fNumGoodNegAdcHits.at(ielem) = 0;
for (UInt_t ielem = 0; ielem < fNumTracksMatched.size(); ielem++)
fNumTracksMatched.at(ielem) = 0;
for (UInt_t ielem = 0; ielem < fNumTracksFired.size(); ielem++)
fNumTracksFired.at(ielem) = 0;
for (UInt_t ielem = 0; ielem < fGoodPosAdcPed.size(); ielem++) {
fGoodPosAdcPed.at(ielem) = 0.0;
fGoodPosAdcPulseInt.at(ielem) = 0.0;
fGoodPosAdcPulseIntRaw.at(ielem) = 0.0;
fGoodPosAdcPulseAmp.at(ielem) = 0.0;
fGoodPosAdcTdcDiffTime.at(ielem) = kBig;
fPosNpe.at(ielem) = 0.0;
}
for (UInt_t ielem = 0; ielem < fGoodNegAdcPed.size(); ielem++) {
fGoodNegAdcPed.at(ielem) = 0.0;
fGoodNegAdcPulseInt.at(ielem) = 0.0;
fGoodNegAdcPulseIntRaw.at(ielem) = 0.0;
fGoodNegAdcPulseAmp.at(ielem) = 0.0;
fGoodNegAdcTdcDiffTime.at(ielem) = kBig;
fNegNpe.at(ielem) = 0.0;
}
// 6 GeV variables
fNGoodHits = 0;
fNADCPosHits = 0;
fNADCNegHits = 0;
fNTDCPosHits = 0;
fNTDCNegHits = 0;
fNpeSumSixGev = 0.0;
fPosNpeSumSixGev = 0.0;
fNegNpeSumSixGev = 0.0;
fPosTDCHits->Clear();
fNegTDCHits->Clear();
fPosADCHits->Clear();
fNegADCHits->Clear();
for (UInt_t ielem = 0; ielem < fPosNpeSixGev.size(); ielem++)
fPosNpeSixGev.at(ielem) = 0.0;
for (UInt_t ielem = 0; ielem < fNegNpeSixGev.size(); ielem++)
fNegNpeSixGev.at(ielem) = 0.0;
for(Int_t itube = 0;itube < fNelem;itube++) {
fA_Pos[itube] = 0;
fA_Neg[itube] = 0;
fA_Pos_p[itube] = 0;
fA_Neg_p[itube] = 0;
fT_Pos[itube] = 0;
fT_Neg[itube] = 0;
//_____________________________________________________________________________
Int_t THcAerogel::Decode( const THaEvData& evdata )
{
// Get the Hall C style hitlist (fRawHitList) for this event
Bool_t present = kTRUE; // Suppress reference time warnings
if(fPresentP) { // if this spectrometer not part of trigger
present = *fPresentP;
}
fNhits = DecodeToHitList(evdata, !present);
if (fSixGevData) {
if(gHaCuts->Result("Pedestal_event")) {
AccumulatePedestals(fRawHitList);
fAnalyzePedestals = 1; // Analyze pedestals first normal events
return(0);
}
if(fAnalyzePedestals) {
CalculatePedestals();
Print("");
fAnalyzePedestals = 0; // Don't analyze pedestals next event
}
UInt_t nrPosAdcHits = 0;
UInt_t nrNegAdcHits = 0;
THcAerogelHit* hit = (THcAerogelHit*) fRawHitList->At(ihit);
Int_t npmt = hit->fCounter;
THcRawAdcHit& rawPosAdcHit = hit->GetRawAdcHitPos();
THcRawAdcHit& rawNegAdcHit = hit->GetRawAdcHitNeg();
for (UInt_t thit=0; thit<rawPosAdcHit.GetNPulses(); ++thit) {
((THcSignalHit*) frPosAdcPedRaw->ConstructedAt(nrPosAdcHits))->Set(npmt, rawPosAdcHit.GetPedRaw());
((THcSignalHit*) frPosAdcPed->ConstructedAt(nrPosAdcHits))->Set(npmt, rawPosAdcHit.GetPed());
((THcSignalHit*) frPosAdcPulseIntRaw->ConstructedAt(nrPosAdcHits))->Set(npmt, rawPosAdcHit.GetPulseIntRaw(thit));
((THcSignalHit*) frPosAdcPulseInt->ConstructedAt(nrPosAdcHits))->Set(npmt, rawPosAdcHit.GetPulseInt(thit));
((THcSignalHit*) frPosAdcPulseAmpRaw->ConstructedAt(nrPosAdcHits))->Set(npmt, rawPosAdcHit.GetPulseAmpRaw(thit));
((THcSignalHit*) frPosAdcPulseAmp->ConstructedAt(nrPosAdcHits))->Set(npmt, rawPosAdcHit.GetPulseAmp(thit));
((THcSignalHit*) frPosAdcPulseTimeRaw->ConstructedAt(nrPosAdcHits))->Set(npmt, rawPosAdcHit.GetPulseTimeRaw(thit));
((THcSignalHit*) frPosAdcPulseTime->ConstructedAt(nrPosAdcHits))->Set(npmt, rawPosAdcHit.GetPulseTime(thit)+fAdcTdcOffset);
if (rawPosAdcHit.GetPulseAmpRaw(thit) > 0) ((THcSignalHit*) fPosAdcErrorFlag->ConstructedAt(nrPosAdcHits))->Set(npmt, 0);
if (rawPosAdcHit.GetPulseAmpRaw(thit) <= 0) ((THcSignalHit*) fPosAdcErrorFlag->ConstructedAt(nrPosAdcHits))->Set(npmt, 1);
++nrPosAdcHits;
fTotNumAdcHits++;
fTotNumPosAdcHits++;
fNumPosAdcHits.at(npmt-1) = npmt;
for (UInt_t thit=0; thit<rawNegAdcHit.GetNPulses(); ++thit) {
((THcSignalHit*) frNegAdcPedRaw->ConstructedAt(nrNegAdcHits))->Set(npmt, rawNegAdcHit.GetPedRaw());
((THcSignalHit*) frNegAdcPed->ConstructedAt(nrNegAdcHits))->Set(npmt, rawNegAdcHit.GetPed());
((THcSignalHit*) frNegAdcPulseIntRaw->ConstructedAt(nrNegAdcHits))->Set(npmt, rawNegAdcHit.GetPulseIntRaw(thit));
((THcSignalHit*) frNegAdcPulseInt->ConstructedAt(nrNegAdcHits))->Set(npmt, rawNegAdcHit.GetPulseInt(thit));
((THcSignalHit*) frNegAdcPulseAmpRaw->ConstructedAt(nrNegAdcHits))->Set(npmt, rawNegAdcHit.GetPulseAmpRaw(thit));
((THcSignalHit*) frNegAdcPulseAmp->ConstructedAt(nrNegAdcHits))->Set(npmt, rawNegAdcHit.GetPulseAmp(thit));
((THcSignalHit*) frNegAdcPulseTimeRaw->ConstructedAt(nrNegAdcHits))->Set(npmt, rawNegAdcHit.GetPulseTimeRaw(thit));
((THcSignalHit*) frNegAdcPulseTime->ConstructedAt(nrNegAdcHits))->Set(npmt, rawNegAdcHit.GetPulseTime(thit));
if (rawNegAdcHit.GetPulseAmpRaw(thit) > 0) ((THcSignalHit*) fNegAdcErrorFlag->ConstructedAt(nrNegAdcHits))->Set(npmt, 0);
if (rawNegAdcHit.GetPulseAmpRaw(thit) <= 0) ((THcSignalHit*) fNegAdcErrorFlag->ConstructedAt(nrNegAdcHits))->Set(npmt, 1);
++nrNegAdcHits;
fTotNumAdcHits++;
fTotNumNegAdcHits++;
fNumNegAdcHits.at(npmt-1) = npmt;
return ihit;
}
//_____________________________________________________________________________
Int_t THcAerogel::ApplyCorrections( void )
{
return(0);
}
//_____________________________________________________________________________
Int_t THcAerogel::CoarseProcess( TClonesArray& ) //tracks
{
Double_t StartTime = 0.0;
if( fglHod ) StartTime = fglHod->GetStartTime();
//cout << " starttime = " << StartTime << endl;
// Loop over the elements in the TClonesArray
for(Int_t ielem = 0; ielem < frPosAdcPulseInt->GetEntries(); ielem++) {
Int_t npmt = ((THcSignalHit*) frPosAdcPulseInt->ConstructedAt(ielem))->GetPaddleNumber() - 1;
Double_t pulsePed = ((THcSignalHit*) frPosAdcPed->ConstructedAt(ielem))->GetData();
Double_t pulseInt = ((THcSignalHit*) frPosAdcPulseInt->ConstructedAt(ielem))->GetData();
Double_t pulseIntRaw = ((THcSignalHit*) frPosAdcPulseIntRaw->ConstructedAt(ielem))->GetData();
Double_t pulseAmp = ((THcSignalHit*) frPosAdcPulseAmp->ConstructedAt(ielem))->GetData();
Eric Pooser
committed
Double_t pulseTime = ((THcSignalHit*) frPosAdcPulseTime->ConstructedAt(ielem))->GetData();
Double_t adctdcdiffTime = StartTime-pulseTime;
Bool_t errorFlag = ((THcSignalHit*) fPosAdcErrorFlag->ConstructedAt(ielem))->GetData();
//// Bool_t pulseTimeCut = adctdcdiffTime > fAdcTimeWindowMin && adctdcdiffTime < fAdcTimeWindowMax;
Bool_t pulseTimeCut = adctdcdiffTime > fAdcPosTimeWindowMin[npmt] && adctdcdiffTime < fAdcPosTimeWindowMax[npmt];
// By default, the last hit within the timing cut will be considered "good"
if (!errorFlag)
{
fGoodPosAdcMult.at(npmt) += 1;
}
if (!errorFlag && pulseTimeCut) {
fGoodPosAdcPed.at(npmt) = pulsePed;
// cout << " out = " << npmt << " " << frPosAdcPulseInt->GetEntries() << " " <<fGoodPosAdcMult.at(npmt);
fGoodPosAdcPulseInt.at(npmt) = pulseInt;
fGoodPosAdcPulseIntRaw.at(npmt) = pulseIntRaw;
fGoodPosAdcPulseAmp.at(npmt) = pulseAmp;
fGoodPosAdcPulseTime.at(npmt) = pulseTime;
fGoodPosAdcTdcDiffTime.at(npmt) = adctdcdiffTime;
fPosNpe.at(npmt) = fPosGain[npmt]*fGoodPosAdcPulseInt.at(npmt);
fPosNpeSum += fPosNpe.at(npmt);
fTotNumGoodAdcHits++;
fTotNumGoodPosAdcHits++;
fNumGoodPosAdcHits.at(npmt) = npmt + 1;
}
}
// Loop over the elements in the TClonesArray
for(Int_t ielem = 0; ielem < frNegAdcPulseInt->GetEntries(); ielem++) {
Int_t npmt = ((THcSignalHit*) frNegAdcPulseInt->ConstructedAt(ielem))->GetPaddleNumber() - 1;
Double_t pulsePed = ((THcSignalHit*) frNegAdcPed->ConstructedAt(ielem))->GetData();
Double_t pulseInt = ((THcSignalHit*) frNegAdcPulseInt->ConstructedAt(ielem))->GetData();
Double_t pulseIntRaw = ((THcSignalHit*) frNegAdcPulseIntRaw->ConstructedAt(ielem))->GetData();
Double_t pulseAmp = ((THcSignalHit*) frNegAdcPulseAmp->ConstructedAt(ielem))->GetData();
Eric Pooser
committed
Double_t pulseTime = ((THcSignalHit*) frNegAdcPulseTime->ConstructedAt(ielem))->GetData();
Double_t adctdcdiffTime = StartTime-pulseTime;
Bool_t errorFlag = ((THcSignalHit*) fNegAdcErrorFlag->ConstructedAt(ielem))->GetData();
//// Bool_t pulseTimeCut = adctdcdiffTime > fAdcTimeWindowMin && adctdcdiffTime < fAdcTimeWindowMax;
Bool_t pulseTimeCut = adctdcdiffTime > fAdcNegTimeWindowMin[npmt] && adctdcdiffTime < fAdcNegTimeWindowMax[npmt];
if (!errorFlag)
{
fGoodNegAdcMult.at(npmt) += 1;
}
// By default, the last hit within the timing cut will be considered "good"
if (!errorFlag && pulseTimeCut) {
fGoodNegAdcPed.at(npmt) = pulsePed;
fGoodNegAdcPulseIntRaw.at(npmt) = pulseIntRaw;
fGoodNegAdcPulseAmp.at(npmt) = pulseAmp;
fGoodNegAdcPulseInt.at(npmt) = pulseInt;
fGoodNegAdcPulseTime.at(npmt) = pulseTime;
fGoodNegAdcTdcDiffTime.at(npmt) = adctdcdiffTime;
fNegNpe.at(npmt) = fNegGain[npmt]*fGoodNegAdcPulseInt.at(npmt);
fNegNpeSum += fNegNpe.at(npmt);
fTotNumGoodAdcHits++;
fTotNumGoodNegAdcHits++;
fNumGoodNegAdcHits.at(npmt) = npmt + 1;
}
}
fNpeSum = fNegNpeSum + fPosNpeSum;
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
for(Int_t ihit=0; ihit < fNhits; ihit++) {
Int_t nPosTDCHits = 0;
Int_t nNegTDCHits = 0;
Int_t nPosADCHits = 0;
Int_t nNegADCHits = 0;
// 6 GeV calculations
Int_t adc_pos;
Int_t adc_neg;
Int_t tdc_pos = -1;
Int_t tdc_neg = -1;
if (fSixGevData) {
THcAerogelHit* hit = (THcAerogelHit*) fRawHitList->At(ihit);
Int_t npmt = hit->fCounter - 1;
// Sum positive and negative hits to fill tot_good_hits
if(fPosNpe.at(npmt) > 0.3) {fNADCPosHits++; fNGoodHits++;}
if(fNegNpe.at(npmt) > 0.3) {fNADCNegHits++; fNGoodHits++;}
// ADC positive hit
if((adc_pos = hit->GetRawAdcHitPos().GetPulseInt()) > 0) {
THcSignalHit *sighit = (THcSignalHit*) fPosADCHits->ConstructedAt(nPosADCHits++);
sighit->Set(hit->fCounter, adc_pos);
}
// ADC negative hit
if((adc_neg = hit->GetRawAdcHitNeg().GetPulseInt()) > 0) {
THcSignalHit *sighit = (THcSignalHit*) fNegADCHits->ConstructedAt(nNegADCHits++);
sighit->Set(hit->fCounter, adc_neg);
}
// TDC positive hit
if(hit->GetRawTdcHitPos().GetNHits() > 0) {
THcSignalHit *sighit = (THcSignalHit*) fPosTDCHits->ConstructedAt(nPosTDCHits++);
tdc_pos = hit->GetRawTdcHitPos().GetTime()+fTdcOffset;
sighit->Set(hit->fCounter, tdc_pos);
}
// TDC negative hit
if(hit->GetRawTdcHitNeg().GetNHits() > 0) {
THcSignalHit *sighit = (THcSignalHit*) fNegTDCHits->ConstructedAt(nNegTDCHits++);
tdc_neg = hit->GetRawTdcHitNeg().GetTime()+fTdcOffset;
sighit->Set(hit->fCounter, tdc_neg);
}
// For each TDC, identify the first hit that is positive.
tdc_pos = -1;
tdc_neg = -1;
for(UInt_t thit=0; thit<hit->GetRawTdcHitPos().GetNHits(); thit++) {
Int_t tdc = hit->GetRawTdcHitPos().GetTime(thit);
if(tdc >=0 ) {
tdc_pos = tdc;
break;
}
}
for(UInt_t thit=0; thit<hit->GetRawTdcHitNeg().GetNHits(); thit++) {
Int_t tdc = hit->GetRawTdcHitNeg().GetTime(thit);
if(tdc >= 0) {
tdc_neg = tdc;
break;
}
}
fA_Pos[npmt] = adc_pos;
fA_Neg[npmt] = adc_neg;
fA_Pos_p[npmt] = fA_Pos[npmt] - fPosPedMean[npmt];
fA_Neg_p[npmt] = fA_Neg[npmt] - fNegPedMean[npmt];
fT_Pos[npmt] = tdc_pos;
fT_Neg[npmt] = tdc_neg;
if(fA_Pos[npmt] < 8000) fPosNpeSixGev[npmt] = fPosGain[npmt]*fA_Pos_p[npmt];
else fPosNpeSixGev[npmt] = 100.0;
if(fA_Neg[npmt] < 8000) fNegNpeSixGev[npmt] = fNegGain[npmt]*fA_Neg_p[npmt];
else fNegNpeSixGev[npmt] = 100.0;
fPosNpeSumSixGev += fPosNpeSixGev[npmt];
fNegNpeSumSixGev += fNegNpeSixGev[npmt];
// Sum positive and negative hits to fill tot_good_hits
if(fPosNpeSixGev[npmt] > 0.3) {fNADCPosHits++; fNGoodHits++;}
if(fNegNpeSixGev[npmt] > 0.3) {fNADCNegHits++; fNGoodHits++;}
if(fT_Pos[npmt] > 0 && fT_Pos[npmt] < 8000) fNTDCPosHits++;
if(fT_Neg[npmt] > 0 && fT_Neg[npmt] < 8000) fNTDCNegHits++;
}
if (fPosNpeSumSixGev > 0.5 || fNegNpeSumSixGev > 0.5)
fNpeSumSixGev = fPosNpeSumSixGev + fNegNpeSumSixGev;
else fNpeSumSixGev = 0.0;
// If total hits are 0, then give a noticable ridiculous NPE
if (fNhits < 1) fNpeSumSixGev = 0.0;
}
return 0;
}
//_____________________________________________________________________________
Int_t THcAerogel::FineProcess( TClonesArray& tracks )
{
Int_t nTracks = tracks.GetLast() + 1;
for (Int_t itrack = 0; itrack < nTracks; itrack++) {
THaTrack* track = dynamic_cast<THaTrack*> (tracks[itrack]);
if (track->GetIndex() != 0) continue; // Select the best track
Double_t trackChi2 = track->GetChi2();
Int_t trackNDoF = track->GetNDoF();
Double_t trackRedChi2 = trackChi2/trackNDoF;
Double_t trackBeta = track->GetBeta();
Double_t trackEnergy = track->GetEnergy();
Double_t trackMom = track->GetP();
Double_t trackENorm = trackEnergy/trackMom;
Double_t trackXfp = track->GetX();
Double_t trackYfp = track->GetY();
Double_t trackTheta = track->GetTheta();
Double_t trackPhi = track->GetPhi();
Bool_t trackRedChi2Cut = trackRedChi2 > fRedChi2Min && trackRedChi2 < fRedChi2Max;
Bool_t trackBetaCut = trackBeta > fBetaMin && trackBeta < fBetaMax;
Bool_t trackENormCut = trackENorm > fENormMin && trackENorm < fENormMax;
Bool_t trackDpCut = trackDp > fDpMin && trackDp < fDpMax;
fXAtAero = trackXfp + trackTheta * fDiffBoxZPos;
fYAtAero = trackYfp + trackPhi * fDiffBoxZPos;
if (trackRedChi2Cut && trackBetaCut && trackENormCut && trackDpCut) {
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
// Project the track to the Aerogel diffuser box plane
// cout << "Aerogel Detector: " << GetName() << endl;
// cout << "nTracks = " << nTracks << "\t" << "trackChi2 = " << trackChi2
// << "\t" << "trackNDof = " << trackNDoF << "\t" << "trackRedChi2 = " << trackRedChi2 << endl;
// cout << "trackBeta = " << trackBeta << "\t" << "trackEnergy = " << trackEnergy << "\t"
// << "trackMom = " << trackMom << "\t" << "trackENorm = " << trackENorm << endl;
// cout << "trackXfp = " << trackXfp << "\t" << "trackYfp = " << trackYfp << "\t"
// << "trackTheta = " << trackTheta << "\t" << "trackPhi = " << trackPhi << endl;
// cout << "fDiffBoxZPos = " << fDiffBoxZPos << "\t" << "fXAtAero = " << fXAtAero << "\t" << "fYAtAero = " << fYAtAero << endl;
// cout << "=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:" << endl;
for (Int_t iregion = 0; iregion < fNRegions; iregion++) {
if ((TMath::Abs(fRegionValue[GetIndex(iregion, 0)] - fXAtAero) < fRegionValue[GetIndex(iregion, 4)]) &&
(TMath::Abs(fRegionValue[GetIndex(iregion, 1)] - fYAtAero) < fRegionValue[GetIndex(iregion, 5)]) &&
(TMath::Abs(fRegionValue[GetIndex(iregion, 2)] - trackTheta) < fRegionValue[GetIndex(iregion, 6)]) &&
(TMath::Abs(fRegionValue[GetIndex(iregion, 3)] - trackPhi) < fRegionValue[GetIndex(iregion, 7)])) {
fTotNumTracksMatched++;
fNumTracksMatched.at(iregion) = iregion + 1;
if (fNpeSum > fNpeThresh) {
fTotNumTracksFired++;
fNumTracksFired.at(iregion) = iregion + 1;
} // NPE threshold cut
} // Regional cuts
} // Loop over regions
} // Tracking cuts
} // Track loop
return 0;
}
//_____________________________________________________________________________
// Method for initializing pedestals in the 6 GeV era
void THcAerogel::InitializePedestals()
fMinPeds = 0; // Do not calculate pedestals by default
DeletePedestalArrays();
fPosPedLimit = new Int_t [fNelem];
fNegPedLimit = new Int_t [fNelem];
fPosPedMean = new Double_t[fNelem];
fNegPedMean = new Double_t[fNelem];
fPosPedSum = new Int_t [fNelem];
fPosPedSum2 = new Int_t [fNelem];
fPosPedCount = new Int_t [fNelem];
fNegPedSum = new Int_t [fNelem];
fNegPedSum2 = new Int_t [fNelem];
fNegPedCount = new Int_t [fNelem];
fPosPed = new Double_t [fNelem];
fPosThresh = new Double_t [fNelem];
fNegPed = new Double_t [fNelem];
fNegThresh = new Double_t [fNelem];
for(Int_t i = 0;i < fNelem; i++) {
fPosPedSum[i] = 0;
fPosPedSum2[i] = 0;
fPosPedLimit[i] = 1000; // In engine, this are set in parameter file
fNegPedSum[i] = 0;
fNegPedSum2[i] = 0;
fNegPedLimit[i] = 1000; // In engine, this are set in parameter file
fPosPedMean[i] = 0; // Default pedestal values
fNegPedMean[i] = 0; // Default pedestal values
}
//_____________________________________________________________________________
// Method for accumulating pedestals in the 6 GeV era
void THcAerogel::AccumulatePedestals(TClonesArray* rawhits)
{
// Extract data from the hit list, accumulating into arrays for
// calculating pedestals
Int_t nrawhits = rawhits->GetLast()+1;
while(ihit < nrawhits) {
THcAerogelHit* hit = (THcAerogelHit *) rawhits->At(ihit);
Int_t element = hit->fCounter - 1;
Int_t adcpos = hit->GetRawAdcHitPos().GetPulseInt();
Int_t adcneg = hit->GetRawAdcHitNeg().GetPulseInt();
if(adcpos <= fPosPedLimit[element]) {
fPosPedSum[element] += adcpos;
fPosPedSum2[element] += adcpos*adcpos;
fPosPedCount[element]++;
if(fPosPedCount[element] == fMinPeds/5)
fPosPedLimit[element] = 100 + fPosPedSum[element]/fPosPedCount[element];
}
if(adcneg <= fNegPedLimit[element]) {
fNegPedSum[element] += adcneg;
fNegPedSum2[element] += adcneg*adcneg;
fNegPedCount[element]++;
if(fNegPedCount[element] == fMinPeds/5)
fNegPedLimit[element] = 100 + fNegPedSum[element]/fNegPedCount[element];
}
ihit++;