33 double tscattmax=0, tscatttot=0, tscatt;
39 dir = (traj2-traj1).
Unit();
40 if((dir-olddir).Mag()<1
e-7)
continue;
41 if(
p>0&&dir.Mag()<=1.0&&olddir.Mag()<=1.0&&(traj1-vert).
Mag()>dist&&(traj2-vert).
Mag()>
dist) {
42 tscatt =
std::acos(dir.Dot(olddir)) * 180 / TMath::Pi();
43 if(tscatt > tscattmax) tscattmax = tscatt;
56 TVector3 traj1, traj2;
62 for (
unsigned int trackHit = 0; trackHit < trackHits.
size(); ++trackHit){
64 double xyz[3] = {0.,0.,0.}, junk = 0;
67 if(trackHits[trackHit]->
View()==
geo::kX) dist2d =
sqrt(
pow(xyz[0]-vert.X(),2) +
pow(xyz[2]-vert.Z(),2));
68 else if(trackHits[trackHit]->
View()==
geo::kY) dist2d =
sqrt(
pow(xyz[1]-vert.Y(),2) +
pow(xyz[2]-vert.Z(),2));
69 if(dist2d>dist) survivingCluster.
Add(trackHits[trackHit]);
72 if(survivingCluster.
NCell()>0) gev = survivingCluster.
TotalGeV();
82 if((traj1-vert).
Mag() > dist && (traj2-vert).
Mag() > dist)
length += (traj1-traj2).
Mag();
83 else if((traj1-vert).Mag() < dist && (traj2-vert).
Mag() <
dist)
continue;
86 if((traj1-vert).Mag() <
dist) {
88 dir = (traj2-traj1).
Unit();
92 dir = (traj1-traj2).
Unit();
94 TVector3 point =
start;
95 bool keepgoing =
true;
98 if( (point - vert).Mag() > (dist - 0.5)) {
102 if((point-start).Mag() > (traj1-traj2).
Mag() + 1) {
104 std::cout <<
"SOMETHING WENT WRONG IN PION RECO DEDX CALC, TOO FAR DISTANCE" <<
std::endl;
113 else if(
length == 0)
return 0;
120 TVector3 trackend = track->
Stop();
126 for (
unsigned int sliceHit = 0; sliceHit < sliceHits.
size(); ++sliceHit){
128 for (
unsigned int trackHit = 0; trackHit < trackHits.
size(); ++trackHit){
129 match = (sliceHits[sliceHit] == trackHits[trackHit]);
134 double xyz[3] = {0.,0.,0.}, junk = 0;
136 float disttovert = 0, disttoend = 0;
138 disttovert =
sqrt(
pow(xyz[0]-vert.X(),2) +
pow(xyz[2]-vert.Z(),2));
139 disttoend =
sqrt(
pow(xyz[0]-trackend.X(),2) +
pow(xyz[2]-trackend.Z(),2));
142 disttovert =
sqrt(
pow(xyz[1]-vert.Y(),2) +
pow(xyz[2]-vert.Z(),2));
143 disttoend =
sqrt(
pow(xyz[1]-trackend.Y(),2) +
pow(xyz[2]-trackend.Z(),2));
145 if (!match && disttovert>dist && disttoend<activitydist) activityCluster.
Add(sliceHits[sliceHit]);
147 if(activityCluster.
NCell()>0)
return activityCluster.
TotalGeV();
157 TVector3 prox,
min = {9999, 9999, 9999}, start1, start2, end1, end2;
159 for(
uint i=0;
i<sliceTracks.size()-1;
i++) {
161 if((
int)
i==trknum)
continue;
165 start1 = first->
Start();
166 start2 = second->
Start();
168 end1 = first->
Stop();
169 end2 = second->
Stop();
174 prox = start1-start2;
175 if(prox.Mag() < min.Mag() && (start1-vert).
Mag() > dist && (start2-vert).
Mag() >
dist) min = prox;
179 if(prox.Mag() < min.Mag() && (start1-vert).
Mag() > dist && (end2-vert).
Mag() >
dist) min = prox;
183 if(prox.Mag() < min.Mag() && (end1-vert).
Mag() > dist && (start2-vert).
Mag() >
dist) min = prox;
187 if(prox.Mag() < min.Mag() && (end1-vert).
Mag() > dist && (end2-vert).
Mag() >
dist) min = prox;
189 if(min.Mag()>17000)
return -1;
197 float gev = 0,
length = 0;
198 TVector3 traj1, traj2;
204 for (
unsigned int prongHit = 0; prongHit < prongHits.
size(); ++prongHit){
206 double xyz[3] = {0.,0.,0.}, junk = 0;
208 TVector3 cellcenter(xyz);
210 if(prongHits[prongHit]->
View()==
geo::kX) dist2d =
sqrt(
pow(xyz[0]-vert.X(),2) +
pow(xyz[2]-vert.Z(),2));
211 else if(prongHits[prongHit]->
View()==
geo::kY) dist2d =
sqrt(
pow(xyz[1]-vert.Y(),2) +
pow(xyz[2]-vert.Z(),2));
212 if(dist2d>dist) survivingCluster.
Add(prongHits[prongHit]);
215 if(survivingCluster.
NCell()>0) gev = survivingCluster.
TotalGeV();
218 traj1 = prong->
Start();
221 if((traj1-vert).
Mag() > dist && (traj2-vert).
Mag() > dist)
length += (traj1-traj2).
Mag();
222 else if((traj1-vert).Mag() < dist && (traj2-vert).
Mag() <
dist)
length = 0;
225 if((traj1-vert).Mag() <
dist) {
227 dir = (traj2-traj1).
Unit();
231 dir = (traj1-traj2).
Unit();
233 TVector3 point =
start;
234 bool keepgoing =
true;
237 if( (point - vert).Mag() > (dist - 0.5)) {
241 if((point-start).Mag() > (traj1-traj2).
Mag() + 1) {
243 std::cout <<
"SOMETHING WENT WRONG IN PION RECO DEDX CALC, TOO FAR DISTANCE" <<
std::endl;
251 else if(
length == 0)
return 0;
262 for (
unsigned int sliceHit = 0; sliceHit < sliceHits.
size(); ++sliceHit){
264 for (
unsigned int prongHit = 0; prongHit < prongHits.
size(); ++prongHit){
265 match = (sliceHits[sliceHit] == prongHits[prongHit]);
270 double xyz[3] = {0.,0.,0.}, junk = 0;
272 float disttovert = 0, disttoend = 0;
274 disttovert =
sqrt(
pow(xyz[0]-vert.X(),2) +
pow(xyz[2]-vert.Z(),2));
275 disttoend =
sqrt(
pow(xyz[0]-prongend.X(),2) +
pow(xyz[2]-prongend.Z(),2));
278 disttovert =
sqrt(
pow(xyz[1]-vert.Y(),2) +
pow(xyz[2]-vert.Z(),2));
279 disttoend =
sqrt(
pow(xyz[1]-prongend.Y(),2) +
pow(xyz[2]-prongend.Z(),2));
281 if (!match && disttovert>dist && disttoend<activitydist) activityCluster.
Add(sliceHits[sliceHit]);
283 if(activityCluster.
NCell()>0)
return activityCluster.
TotalGeV();
293 TVector3 prox,
min = {9999, 9999, 9999}, start1, start2, end1, end2;
295 for(
uint i=0;
i<sliceProngs.size()-1;
i++) {
297 if((
int)
i==pngnum)
continue;
301 start1 = first->
Start();
302 start2 = second->
Start();
310 prox = start1-start2;
311 if(prox.Mag() < min.Mag() && (start1-vert).
Mag() > dist && (start2-vert).
Mag() >
dist) min = prox;
315 if(prox.Mag() < min.Mag() && (start1-vert).
Mag() > dist && (end2-vert).
Mag() >
dist) min = prox;
319 if(prox.Mag() < min.Mag() && (end1-vert).
Mag() > dist && (start2-vert).
Mag() >
dist) min = prox;
323 if(prox.Mag() < min.Mag() && (end1-vert).
Mag() > dist && (end2-vert).
Mag() >
dist) min = prox;
325 if(min.Mag()>17000)
return -1;
size_t NTrajectoryPoints() const
unsigned int NCell(geo::View_t view) const
Number of cells in view view.
void GetCenter(double *xyz, double localz=0.0) const
TVector3 TrajectoryPoint(unsigned int i) const
The ith point on the trajectory, a 3-vector in cm.
const CellGeo * Cell(int icell) const
float getProngActivity(const art::Ptr< rb::Prong > prong, TVector3 vert, float dist, art::PtrVector< rb::CellHit > const &sliceHits, float activitydist)
float getTrackProximity(const std::vector< art::Ptr< rb::Track > > sliceTracks, int trknum, TVector3 vertex, float dist)
Vertical planes which measure X.
A collection of associated CellHits.
float getProngProximity(const std::vector< art::Ptr< rb::Prong > > sliceProngs, int pngnum, TVector3 vertex, float dist)
const PlaneGeo * Plane(unsigned int i) const
virtual TVector3 Start() const
float getTrackScatt(const art::Ptr< rb::Track > track, TVector3 vert, float dist)
float getTrackDedx(const art::Ptr< rb::Track > track, TVector3 vert, float dist)
Horizontal planes which measure Y.
art::PtrVector< rb::CellHit > AllCells() const
Get all cells from both views.
virtual void Add(const art::Ptr< rb::CellHit > &cell, double weight=1)
virtual double TotalLength() const
Distance along prong to reach last cell hit.
olddir
Are we streaming the files via xrootd? txtcmd="cat %s | xargs -n1 samweb2xrootd > xrootd_inFile...
float getProngDedx(const art::Ptr< rb::Prong > prong, TVector3 vert, float dist)
virtual TVector3 Dir() const
Unit vector describing prong direction.
Geometry information for a single readout plane.
Collect Geo headers and supply basic geometry functions.
static float min(const float a, const float b, const float c)
float getTrackActivity(const art::Ptr< rb::Track > track, TVector3 vert, float dist, art::PtrVector< rb::CellHit > const &sliceHits, float activitydist)
double TotalGeV(EEnergyCalcScheme escheme=kRecomputeEnergy) const
Simple sum of the estimated GeV of all the hits.
TVector3 Stop() const
Position of the final trajectory point.
produce variables useful to reconstruction pion energy
Encapsulate the geometry of one entire detector (near, far, ndos)