Skip to content

Commit

Permalink
Merge pull request #219 from JeffersonLab/dirc_tree_fix
Browse files Browse the repository at this point in the history 
Many DIRC changes to correctly add North OB to monitoring and several analysis updates
  • Loading branch information
@sdobbs
sdobbs authored Oct 29, 2019
2 parents 1b69309 + 8a20e89 commit 72c8b90
Show file tree
Hide file tree
Showing 26 changed files with 403 additions and 163 deletions.
29 changes: 23 additions & 6 deletions src/libraries/ANALYSIS/DEventWriterROOT.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -691,6 +691,8 @@ void DEventWriterROOT::Create_Branches_ChargedHypotheses(DTreeBranchRegister& lo
//DIRC:
if(DIRC_OUTPUT) {
locBranchRegister.Register_FundamentalArray<Int_t>(Build_BranchName(locParticleBranchName, "NumPhotons_DIRC"), locArraySizeString, dInitNumTrackArraySize);
locBranchRegister.Register_FundamentalArray<Float_t>(Build_BranchName(locParticleBranchName, "ExtrapolatedX_DIRC"), locArraySizeString, dInitNumTrackArraySize);
locBranchRegister.Register_FundamentalArray<Float_t>(Build_BranchName(locParticleBranchName, "ExtrapolatedY_DIRC"), locArraySizeString, dInitNumTrackArraySize);
locBranchRegister.Register_FundamentalArray<Float_t>(Build_BranchName(locParticleBranchName, "ThetaC_DIRC"), locArraySizeString, dInitNumTrackArraySize);
locBranchRegister.Register_FundamentalArray<Float_t>(Build_BranchName(locParticleBranchName, "Lele_DIRC"), locArraySizeString, dInitNumTrackArraySize);
locBranchRegister.Register_FundamentalArray<Float_t>(Build_BranchName(locParticleBranchName, "Lpi_DIRC"), locArraySizeString, dInitNumTrackArraySize);
Expand Down Expand Up @@ -1646,17 +1648,32 @@ void DEventWriterROOT::Fill_ChargedHypo(DTreeFillData* locTreeFillData, unsigned

// DIRC
if(DIRC_OUTPUT) {
int locDIRCNumPhotons = (locChargedTrackHypothesis->Get_DIRCMatchParams() != NULL) ? locChargedTrackHypothesis->Get_DIRCMatchParams()->dNPhotons : 0;
int locDIRCNumPhotons = 0;
double locDIRCExtrapolatedX = 999;
double locDIRCExtrapolatedY = 999;
double locDIRCThetaC = 999.;
double locDIRCLele = 999.;
double locDIRCLpi = 999.;
double locDIRCLk = 999.;
double locDIRCLp = 999.;
auto locDIRCMatchParams = locChargedTrackHypothesis->Get_DIRCMatchParams();
if(locDIRCMatchParams != NULL) {
locDIRCExtrapolatedX = locDIRCMatchParams->dExtrapolatedPos.X();
locDIRCExtrapolatedY = locDIRCMatchParams->dExtrapolatedPos.Y();
locDIRCNumPhotons = locDIRCMatchParams->dNPhotons;
locDIRCThetaC = locDIRCMatchParams->dThetaC;
locDIRCLele = locDIRCMatchParams->dLikelihoodElectron;
locDIRCLpi = locDIRCMatchParams->dLikelihoodPion;
locDIRCLk = locDIRCMatchParams->dLikelihoodKaon;
locDIRCLp = locDIRCMatchParams->dLikelihoodProton;
}
locTreeFillData->Fill_Array<Int_t>(Build_BranchName(locParticleBranchName, "NumPhotons_DIRC"), locDIRCNumPhotons, locArrayIndex);
double locDIRCThetaC = (locChargedTrackHypothesis->Get_DIRCMatchParams() != NULL) ? locChargedTrackHypothesis->Get_DIRCMatchParams()->dThetaC : 0.0;
locTreeFillData->Fill_Array<Float_t>(Build_BranchName(locParticleBranchName, "ExtrapolatedX_DIRC"), locDIRCExtrapolatedX, locArrayIndex);
locTreeFillData->Fill_Array<Float_t>(Build_BranchName(locParticleBranchName, "ExtrapolatedY_DIRC"), locDIRCExtrapolatedY, locArrayIndex);
locTreeFillData->Fill_Array<Float_t>(Build_BranchName(locParticleBranchName, "ThetaC_DIRC"), locDIRCThetaC, locArrayIndex);
double locDIRCLele = (locChargedTrackHypothesis->Get_DIRCMatchParams() != NULL) ? locChargedTrackHypothesis->Get_DIRCMatchParams()->dLikelihoodElectron : 0.0;
locTreeFillData->Fill_Array<Float_t>(Build_BranchName(locParticleBranchName, "Lele_DIRC"), locDIRCLele, locArrayIndex);
double locDIRCLpi = (locChargedTrackHypothesis->Get_DIRCMatchParams() != NULL) ? locChargedTrackHypothesis->Get_DIRCMatchParams()->dLikelihoodPion : 0.0;
locTreeFillData->Fill_Array<Float_t>(Build_BranchName(locParticleBranchName, "Lpi_DIRC"), locDIRCLpi, locArrayIndex);
double locDIRCLk = (locChargedTrackHypothesis->Get_DIRCMatchParams() != NULL) ? locChargedTrackHypothesis->Get_DIRCMatchParams()->dLikelihoodKaon : 0.0;
locTreeFillData->Fill_Array<Float_t>(Build_BranchName(locParticleBranchName, "Lk_DIRC"), locDIRCLk, locArrayIndex);
double locDIRCLp = (locChargedTrackHypothesis->Get_DIRCMatchParams() != NULL) ? locChargedTrackHypothesis->Get_DIRCMatchParams()->dLikelihoodProton : 0.0;
locTreeFillData->Fill_Array<Float_t>(Build_BranchName(locParticleBranchName, "Lp_DIRC"), locDIRCLp, locArrayIndex);
}
}
Expand Down
65 changes: 47 additions & 18 deletions src/libraries/DIRC/DDIRCLut.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -26,19 +26,23 @@ DDIRCLut::DDIRCLut()
gPARMS->SetDefaultParameter("DIRC:TRUTH_PIXELTIME",DIRC_TRUTH_PIXELTIME);

// timing cuts for photons
DIRC_CUT_TDIFFD = 2; // direct cut in ns
DIRC_CUT_TDIFFD = 3; // direct cut in ns
gPARMS->SetDefaultParameter("DIRC:CUT_TDIFFD",DIRC_CUT_TDIFFD);
DIRC_CUT_TDIFFR = 3; // reflected cut in ns
DIRC_CUT_TDIFFR = 4; // reflected cut in ns
gPARMS->SetDefaultParameter("DIRC:CUT_TDIFFR",DIRC_CUT_TDIFFR);

// Gives DeltaT = 0, but shouldn't it be v=20.3767 [cm/ns] for 1.47125
DIRC_LIGHT_V = 19.80; // v=19.8 [cm/ns] for 1.5141
gPARMS->SetDefaultParameter("DIRC:LIGHT_V",DIRC_LIGHT_V);

// sigma (thetaC for single photon) in radiams
DIRC_SIGMA_THETAC = 0.0085;
DIRC_SIGMA_THETAC = 0.01;
gPARMS->SetDefaultParameter("DIRC:SIGMA_THETAC",DIRC_SIGMA_THETAC);

// Rotate tracks angle based on bar box survey data
DIRC_ROTATE_TRACK = false;
gPARMS->SetDefaultParameter("DIRC:ROTATE_TRACK",DIRC_ROTATE_TRACK);

// set PID for different passes in debuging histograms
dFinalStatePIDs.push_back(Positron);
dFinalStatePIDs.push_back(PiPlus);
Expand All @@ -64,6 +68,18 @@ bool DDIRCLut::brun(JEventLoop *loop) {
loop->Get(locDIRCGeometry);
dDIRCGeometry = locDIRCGeometry[0];

// rotation angles for bar boxes (reverse sign from survey)
for(int i=0; i<12; i++) {
dRotationX[i] = -0.2134*TMath::DegToRad(); // +pitch // -0.00219;
dRotationY[i] = 0.0963*TMath::DegToRad(); // -yaw // 0.00181
dRotationZ[i] = -0.0355*TMath::DegToRad(); // -roll // -0.00062;
}
for(int i=12; i<24; i++) {
dRotationX[i] = -0.0791*TMath::DegToRad(); // +pitch // -0.00138;
dRotationY[i] = 0.1003*TMath::DegToRad(); // -yaw // 0.00175;
dRotationZ[i] = -0.0381*TMath::DegToRad(); // -roll // -0.00066;
}

return true;
}

Expand Down Expand Up @@ -119,10 +135,8 @@ bool DDIRCLut::CreateDebugHistograms() {
return true;
}

bool DDIRCLut::CalcLUT(TVector3 locProjPos, TVector3 locProjMom, const vector<const DDIRCPmtHit*> locDIRCHits, double locFlightTime, Particle_t locPID, shared_ptr<DDIRCMatchParams>& locDIRCMatchParams, const vector<const DDIRCTruthBarHit*> locDIRCBarHits, map<shared_ptr<const DDIRCMatchParams>, vector<const DDIRCPmtHit*> >& locDIRCTrackMatchParams) const
bool DDIRCLut::CalcLUT(TVector3 locProjPos, TVector3 locProjMom, const vector<const DDIRCPmtHit*> locDIRCHits, double locFlightTime, double locMass, shared_ptr<DDIRCMatchParams>& locDIRCMatchParams, const vector<const DDIRCTruthBarHit*> locDIRCBarHits, map<shared_ptr<const DDIRCMatchParams>, vector<const DDIRCPmtHit*> >& locDIRCTrackMatchParams) const
{
double locMass = ParticleMass(locPID);

// get bar and track position/momentum from extrapolation
TVector3 momInBar = locProjMom;
TVector3 posInBar = locProjPos;
Expand Down Expand Up @@ -156,13 +170,13 @@ bool DDIRCLut::CalcLUT(TVector3 locProjPos, TVector3 locProjMom, const vector<co
locDIRCMatchParams = std::make_shared<DDIRCMatchParams>();

// initialize variables for LUT summary
double locAngle = CalcAngle(momInBar, locMass);
map<Particle_t, double> locExpectedAngle = CalcExpectedAngles(momInBar);
double locAngle = CalcAngle(momInBar.Mag(), locMass);
map<Particle_t, double> locExpectedAngle = CalcExpectedAngles(momInBar.Mag());
map<Particle_t, double> logLikelihoodSum;
Particle_t locHypothesisPID = PiPlus;
for(uint loc_i = 0; loc_i<dFinalStatePIDs.size(); loc_i++) {
logLikelihoodSum[dFinalStatePIDs[loc_i]]=0;
if(fabs(ParticleMass(dFinalStatePIDs[loc_i])-ParticleMass(locPID)) < 0.01) locHypothesisPID = dFinalStatePIDs[loc_i];
if(fabs(ParticleMass(dFinalStatePIDs[loc_i])-locMass) < 0.01) locHypothesisPID = dFinalStatePIDs[loc_i];
}

// loop over DIRC hits
Expand Down Expand Up @@ -234,7 +248,16 @@ vector<pair<double,double>> DDIRCLut::CalcPhoton(const DDIRCPmtHit *locDIRCHit,
int box = (bar < 24) ? 1 : 0;
if((box == 0 && channel < dMaxChannels) || (box == 1 && channel >= dMaxChannels))
return locDIRCPhotons;


int pmt = channel/64;
double rotX = 0, rotY = 0, rotZ=0;
int xbin = (int)(posInBar.X()/5.0) + 19;
if(bar>=0 && bar<=23 && xbin>=0 && xbin<=39) {
rotX = dRotationX[bar];
rotY = dRotationY[bar];
rotZ = dRotationZ[bar];
}

// use hit time to determine if reflected or not
double hitTime = locDIRCHit->t - locFlightTime;

Expand All @@ -247,7 +270,7 @@ vector<pair<double,double>> DDIRCLut::CalcPhoton(const DDIRCPmtHit *locDIRCHit,
}

// needs to be X dependent choice for reflection cut (from CCDB?)
bool reflected = hitTime>0; // try all photons as reflected for now
bool reflected = hitTime>35; // try only some photons as reflected for now

// get position along bar for calculated time
double radiatorL = dDIRCGeometry->GetBarLength(bar);
Expand Down Expand Up @@ -290,10 +313,16 @@ vector<pair<double,double>> DDIRCLut::CalcPhoton(const DDIRCPmtHit *locDIRCHit,
if(r) dir.SetXYZ( -dir.X(), dir.Y(), dir.Z());
if(dir.Angle(fnY1) < dCriticalAngle || dir.Angle(fnZ1) < dCriticalAngle) continue;

TVector3 trackMom = momInBar;
if(DIRC_ROTATE_TRACK) { // rotate tracks to bar plane from survey data
trackMom.RotateX(rotX);
trackMom.RotateY(rotY);
trackMom.RotateZ(rotZ);
}
tangle = trackMom.Angle(dir);

luttheta = dir.Angle(TVector3(-1,0,0));
if(luttheta > TMath::PiOver2()) luttheta = TMath::Pi()-luttheta;
tangle = momInBar.Angle(dir);//-0.002; //correction

double bartime = lenz/cos(luttheta)/DIRC_LIGHT_V;
double totalTime = bartime+evtime;

Expand Down Expand Up @@ -336,7 +365,7 @@ vector<pair<double,double>> DDIRCLut::CalcPhoton(const DDIRCPmtHit *locDIRCHit,
}

// remove photon candidates not used in likelihood
if(fabs(tangle-0.5*(locExpectedAngle[PiPlus]+locExpectedAngle[KPlus]))>0.05) continue;
if(fabs(tangle-0.5*(locExpectedAngle[PiPlus]+locExpectedAngle[KPlus]))>0.03) continue;

// save good photons to matched list
isGood = true;
Expand Down Expand Up @@ -376,15 +405,15 @@ double DDIRCLut::CalcLikelihood(double locExpectedThetaC, double locThetaC) cons
return locLikelihood;
}

double DDIRCLut::CalcAngle(TVector3 momInBar, double locMass) const {
return acos(sqrt(momInBar.Mag()*momInBar.Mag() + locMass*locMass)/momInBar.Mag()/dIndex);
double DDIRCLut::CalcAngle(double locP, double locMass) const {
return acos(sqrt(locP*locP + locMass*locMass)/locP/dIndex);
}

map<Particle_t, double> DDIRCLut::CalcExpectedAngles(TVector3 momInBar) const {
map<Particle_t, double> DDIRCLut::CalcExpectedAngles(double locP) const {

map<Particle_t, double> locExpectedAngles;
for(uint loc_i = 0; loc_i<dFinalStatePIDs.size(); loc_i++) {
locExpectedAngles[dFinalStatePIDs[loc_i]] = acos(sqrt(momInBar.Mag()*momInBar.Mag() + ParticleMass(dFinalStatePIDs[loc_i])*ParticleMass(dFinalStatePIDs[loc_i]))/momInBar.Mag()/dIndex);
locExpectedAngles[dFinalStatePIDs[loc_i]] = acos(sqrt(locP*locP + ParticleMass(dFinalStatePIDs[loc_i])*ParticleMass(dFinalStatePIDs[loc_i]))/locP/dIndex);
}
return locExpectedAngles;
}
9 changes: 6 additions & 3 deletions src/libraries/DIRC/DDIRCLut.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,12 +36,12 @@ class DDIRCLut: public JObject {

bool brun(JEventLoop *loop);
bool CreateDebugHistograms();
bool CalcLUT(TVector3 locProjPos, TVector3 locProjMom, const vector<const DDIRCPmtHit*> locDIRCHits, double locFlightTime, Particle_t locPID, shared_ptr<DDIRCMatchParams>& locDIRCMatchParams, const vector<const DDIRCTruthBarHit*> locDIRCBarHits, map<shared_ptr<const DDIRCMatchParams>, vector<const DDIRCPmtHit*> >& locDIRCTrackMatchParams) const;
bool CalcLUT(TVector3 locProjPos, TVector3 locProjMom, const vector<const DDIRCPmtHit*> locDIRCHits, double locFlightTime, double locMass, shared_ptr<DDIRCMatchParams>& locDIRCMatchParams, const vector<const DDIRCTruthBarHit*> locDIRCBarHits, map<shared_ptr<const DDIRCMatchParams>, vector<const DDIRCPmtHit*> >& locDIRCTrackMatchParams) const;
vector<pair<double,double>> CalcPhoton(const DDIRCPmtHit *locDIRCHit, double locFlightTime, TVector3 posInBar, TVector3 momInBar, map<Particle_t, double> locExpectedAngle, double locAngle, Particle_t locPID, map<Particle_t, double> &logLikelihoodSum, int &nPhotonsThetaC, double &meanThetaC, double &meanDeltaT, bool &isGood) const;
vector<pair<double,double>> CalcPhoton(const DDIRCPmtHit *locDIRCHit, double locFlightTime, TVector3 posInBar, TVector3 momInBar, map<Particle_t, double> locExpectedAngle, double locAngle, Particle_t locPID, map<Particle_t, double> &logLikelihoodSum) const;
double CalcLikelihood(double locExpectedThetaC, double locThetaC) const;
double CalcAngle(TVector3 momInBar, double locMass) const;
map<Particle_t, double> CalcExpectedAngles(TVector3 momInBar) const;
double CalcAngle(double locP, double locMass) const;
map<Particle_t, double> CalcExpectedAngles(double locP) const;

private:
DApplication *dapp;
Expand All @@ -51,12 +51,15 @@ class DDIRCLut: public JObject {
bool DIRC_DEBUG_HISTS;
bool DIRC_TRUTH_BARHIT;
bool DIRC_TRUTH_PIXELTIME;
bool DIRC_ROTATE_TRACK;

double DIRC_CUT_TDIFFD;
double DIRC_CUT_TDIFFR;
double DIRC_SIGMA_THETAC;
double DIRC_LIGHT_V;

double dRotationX[48], dRotationY[48], dRotationZ[48];

int dMaxChannels;
double dCriticalAngle, dIndex;

Expand Down
2 changes: 1 addition & 1 deletion src/libraries/DIRC/DDIRCLutReader.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,7 @@ DDIRCLutReader::DDIRCLutReader(JApplication *japp, unsigned int run_number)
jcalib = japp->GetJCalibration(run_number);
if(jcalib->GetCalib("/DIRC/LUT/lut_map", lut_map_name))
jout << "Can't find requested /DIRC/LUT/lut_map in CCDB for this run!" << endl;
else if(lut_map_name.find("map_name") != lut_map_name.end() && lut_map_name["map_name"] != "None") {
else if(lut_map_name.find("map_name") != lut_map_name.end() && lut_map_name["map_name"] != "None" && lut_file.empty()) {
jresman = japp->GetJResourceManager(run_number);
lut_file = jresman->GetResource(lut_map_name["map_name"]);
}
Expand Down
9 changes: 8 additions & 1 deletion src/libraries/DIRC/DDIRCPmtHit_factory.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,9 @@ using namespace jana;
//------------------
jerror_t DDIRCPmtHit_factory::init(void)
{
DIRC_SKIP = false;
gPARMS->SetDefaultParameter("DIRC:SKIP",DIRC_SKIP);

// initialize calibration tables
vector<double> new_t0s(DIRC_MAX_CHANNELS);
vector<int> new_status(DIRC_MAX_CHANNELS);
Expand Down Expand Up @@ -84,6 +87,9 @@ jerror_t DDIRCPmtHit_factory::evnt(JEventLoop *loop, uint64_t eventnumber)
/// data in HDDM format. The HDDM event source will copy
/// the precalibrated values directly into the _data vector.

if(DIRC_SKIP)
return NOERROR;

// check that SSP board timestamps match for all modules
vector<const DDIRCTriggerTime*> timestamps;
loop->Get(timestamps);
Expand Down Expand Up @@ -125,6 +131,7 @@ jerror_t DDIRCPmtHit_factory::evnt(JEventLoop *loop, uint64_t eventnumber)
if(digihit_lead->channel != digihit_trail->channel) continue;
int channel = digihit_lead->channel;
int box = (channel < DIRC_MAX_CHANNELS) ? 1 : 0; // North=0 and South=1
if(box == 0) channel -= DIRC_MAX_CHANNELS; // box-local channel to index CCDB tables

// get time-over-threshold
timeOverThreshold = (double)digihit_trail->time - (double)digihit_lead->time;
Expand Down Expand Up @@ -152,7 +159,7 @@ jerror_t DDIRCPmtHit_factory::evnt(JEventLoop *loop, uint64_t eventnumber)
double slope = 0.3;
double timeOverThresholdPeak = 50;
if(applyTimeOffset) {
hit->t = hit->t - time_offsets[box][channel] + t_base[box];
hit->t = hit->t + t_base[box] - time_offsets[box][channel];
}
if(applyTimewalk) {
hit->t += slope*(timeOverThreshold - timeOverThresholdPeak);
Expand Down
1 change: 1 addition & 0 deletions src/libraries/DIRC/DDIRCPmtHit_factory.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -44,6 +44,7 @@ class DDIRCPmtHit_factory:public jana::JFactory<DDIRCPmtHit>{

double t_base[2];
bool applyTimeOffset, applyTimewalk;
bool DIRC_SKIP;
};

#endif // _DDIRCPmtHit_factory_
Expand Down
6 changes: 5 additions & 1 deletion src/libraries/HDDM/DEventSourceREST.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1307,6 +1307,10 @@ jerror_t DEventSourceREST::Extract_DDetectorMatches(JEventLoop* locEventLoop, hd
for(; dircIter != dircList.end(); ++dircIter)
{
size_t locTrackIndex = dircIter->getTrack();
if(locTrackIndex > locTrackTimeBasedVector.size()) continue;

auto locTrackTimeBased = locTrackTimeBasedVector[locTrackIndex];
if( !locTrackTimeBased ) continue;

auto locDIRCMatchParams = std::make_shared<DDIRCMatchParams>();
map<shared_ptr<const DDIRCMatchParams> ,vector<const DDIRCPmtHit*> > locDIRCTrackMatchParams;
Expand All @@ -1317,7 +1321,7 @@ jerror_t DEventSourceREST::Extract_DDetectorMatches(JEventLoop* locEventLoop, hd
TVector3 locProjMom(dircIter->getPx(),dircIter->getPy(),dircIter->getPz());
double locFlightTime = dircIter->getT();

if( locParticleID->Get_DIRCLut()->CalcLUT(locProjPos, locProjMom, locDIRCHits, locFlightTime, locTrackTimeBasedVector[locTrackIndex]->PID(), locDIRCMatchParams, locDIRCBarHits, locDIRCTrackMatchParams) )
if( locParticleID->Get_DIRCLut()->CalcLUT(locProjPos, locProjMom, locDIRCHits, locFlightTime, locTrackTimeBased->mass(), locDIRCMatchParams, locDIRCBarHits, locDIRCTrackMatchParams) )
locDetectorMatches->Add_Match(locTrackTimeBasedVector[locTrackIndex], std::const_pointer_cast<const DDIRCMatchParams>(locDIRCMatchParams));
}
else {
Expand Down
1 change: 1 addition & 0 deletions src/libraries/PID/DChargedTrackHypothesis.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -338,6 +338,7 @@ inline void DChargedTrackHypothesis::DTrackingInfo::Reset(void)
dTOFHitMatchParams = nullptr;
dBCALShowerMatchParams = nullptr;
dFCALShowerMatchParams = nullptr;
dDIRCMatchParams = nullptr;
}

#endif // _DChargedTrackHypothesis_
2 changes: 1 addition & 1 deletion src/libraries/PID/DParticleID.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1781,7 +1781,7 @@ bool DParticleID::Cut_MatchDIRC(const vector<DTrackFitter::Extrapolation_t> &ext
}

// Calculate DIRC LUT
return dDIRCLut->CalcLUT(locProjPos, locProjMom, locDIRCHits, locFlightTime, locPID, locDIRCMatchParams, locDIRCBarHits, locDIRCTrackMatchParams);
return dDIRCLut->CalcLUT(locProjPos, locProjMom, locDIRCHits, locFlightTime, ParticleMass(locPID), locDIRCMatchParams, locDIRCBarHits, locDIRCTrackMatchParams);

}

Expand Down
8 changes: 4 additions & 4 deletions src/plugins/Analysis/dirc_hists/DEventProcessor_dirc_hists.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -178,8 +178,8 @@ jerror_t DEventProcessor_dirc_hists::evnt(JEventLoop *loop, uint64_t eventnumber
hExtrapolatedBarHitXY[locPID]->Fill(posInBar.X(), posInBar.Y());
japp->RootFillUnLock(this); //RELEASE ROOT FILL LOCK

double locAngle = dDIRCLut->CalcAngle(momInBar, locMass);
map<Particle_t, double> locExpectedAngle = dDIRCLut->CalcExpectedAngles(momInBar);
double locAngle = dDIRCLut->CalcAngle(locP, locMass);
map<Particle_t, double> locExpectedAngle = dDIRCLut->CalcExpectedAngles(locP);

// get map of DIRCMatches to PMT hits
map<shared_ptr<const DDIRCMatchParams>, vector<const DDIRCPmtHit*> > locDIRCTrackMatchParamsMap;
Expand All @@ -193,7 +193,7 @@ jerror_t DEventProcessor_dirc_hists::evnt(JEventLoop *loop, uint64_t eventnumber
vector<pair<double, double>> locDIRCPhotons = dDIRCLut->CalcPhoton(locDIRCPmtHits[loc_i], locExtrapolatedTime, posInBar, momInBar, locExpectedAngle, locAngle, locPID, logLikelihoodSum);
double locHitTime = locDIRCPmtHits[loc_i]->t - locExtrapolatedTime;
int locChannel = locDIRCPmtHits[loc_i]->ch%dMaxChannels;
if(locHitTime > 0 && locHitTime < 100) locPhotonInclusive++;
if(locHitTime > 0 && locHitTime < 150) locPhotonInclusive++;

int pixel_row = locDIRCGeometry->GetPixelRow(locChannel);
int pixel_col = locDIRCGeometry->GetPixelColumn(locChannel);
Expand Down Expand Up @@ -231,7 +231,7 @@ jerror_t DEventProcessor_dirc_hists::evnt(JEventLoop *loop, uint64_t eventnumber
}

// fill histograms for candidate photons in timing cut
if(fabs(locDeltaT) < 100.0) {
if(fabs(locDeltaT) < 5.0) {
hThetaC[locPID]->Fill(locThetaC);
hDeltaThetaC[locPID]->Fill(locThetaC-locExpectedThetaC);
hDeltaThetaCVsP[locPID]->Fill(momInBar.Mag(), locThetaC-locExpectedThetaC);
Expand Down
Loading

0 comments on commit 72c8b90

Please sign in to comment.