Skip to content
Snippets Groups Projects
Commit c0f3e739 authored by Stephen A. Wood's avatar Stephen A. Wood
Browse files

Replace TrackSP structure with THcDCTrack class

parent 294be162
No related branches found
No related tags found
No related merge requests found
...@@ -20,7 +20,7 @@ SRC = src/THcInterface.cxx src/THcParmList.cxx src/THcAnalyzer.cxx \ ...@@ -20,7 +20,7 @@ SRC = src/THcInterface.cxx src/THcParmList.cxx src/THcAnalyzer.cxx \
src/THcRawDCHit.cxx src/THcDCHit.cxx \ src/THcRawDCHit.cxx src/THcDCHit.cxx \
src/THcDCWire.cxx \ src/THcDCWire.cxx \
src/THcDCLookupTTDConv.cxx src/THcDCTimeToDistConv.cxx \ src/THcDCLookupTTDConv.cxx src/THcDCTimeToDistConv.cxx \
src/THcSpacePoint.cxx \ src/THcSpacePoint.cxx src/THcDCTrack.cxx \
src/THcShower.cxx src/THcShowerPlane.cxx \ src/THcShower.cxx src/THcShowerPlane.cxx \
src/THcShowerHit.cxx \ src/THcShowerHit.cxx \
src/THcAerogel.cxx src/THcAerogelHit.cxx src/THcAerogel.cxx src/THcAerogelHit.cxx
......
...@@ -27,6 +27,7 @@ ...@@ -27,6 +27,7 @@
#pragma link C++ class THcDCLookupTTDConv+; #pragma link C++ class THcDCLookupTTDConv+;
#pragma link C++ class THcDCTimeToDistConv+; #pragma link C++ class THcDCTimeToDistConv+;
#pragma link C++ class THcSpacePoint+; #pragma link C++ class THcSpacePoint+;
#pragma link C++ class THcDCTrack+;
#pragma link C++ class THcShower+; #pragma link C++ class THcShower+;
#pragma link C++ class THcShowerPlane+; #pragma link C++ class THcShowerPlane+;
#pragma link C++ class THcShowerHit+; #pragma link C++ class THcShowerHit+;
......
...@@ -15,6 +15,7 @@ ...@@ -15,6 +15,7 @@
#include "THcDetectorMap.h" #include "THcDetectorMap.h"
#include "THcGlobals.h" #include "THcGlobals.h"
#include "THcParmList.h" #include "THcParmList.h"
#include "THcDCTrack.h"
#include "VarDef.h" #include "VarDef.h"
#include "VarType.h" #include "VarType.h"
#include "THaTrack.h" #include "THaTrack.h"
...@@ -65,6 +66,8 @@ THcDC::THcDC( ...@@ -65,6 +66,8 @@ THcDC::THcDC(
fCentralWire = NULL; fCentralWire = NULL;
fPlaneTimeZero = NULL; fPlaneTimeZero = NULL;
fSigma = NULL; fSigma = NULL;
fDCTracks = new TClonesArray( "THcDCTrack", 20 );
} }
//_____________________________________________________________________________ //_____________________________________________________________________________
...@@ -352,6 +355,7 @@ THcDC::~THcDC() ...@@ -352,6 +355,7 @@ THcDC::~THcDC()
fTrackProj->Clear(); fTrackProj->Clear();
delete fTrackProj; fTrackProj = 0; delete fTrackProj; fTrackProj = 0;
} }
delete fDCTracks;
} }
//_____________________________________________________________________________ //_____________________________________________________________________________
...@@ -460,7 +464,7 @@ Int_t THcDC::CoarseTrack( TClonesArray& /* tracks */ ) ...@@ -460,7 +464,7 @@ Int_t THcDC::CoarseTrack( TClonesArray& /* tracks */ )
} }
// Now link the stubs between chambers // Now link the stubs between chambers
LinkStubs(); LinkStubs();
if(ntracks_fp > 0) TrackFit(); if(fNDCTracks > 0) TrackFit();
// Check for internal TrackFit errors // Check for internal TrackFit errors
// Histogram the focal plane tracks // Histogram the focal plane tracks
// Histograms made in h_fill_dc_fp_hist // Histograms made in h_fill_dc_fp_hist
...@@ -529,7 +533,8 @@ void THcDC::LinkStubs() ...@@ -529,7 +533,8 @@ void THcDC::LinkStubs()
fNSp++; fNSp++;
} }
} }
ntracks_fp=0; // Number of Focal Plane tracks found fNDCTracks=0; // Number of Focal Plane tracks found
fDCTracks->Clear();
Double_t stubminx = 999999; Double_t stubminx = 999999;
Double_t stubminy = 999999; Double_t stubminy = 999999;
Double_t stubminxp = 999999; Double_t stubminxp = 999999;
...@@ -542,9 +547,10 @@ void THcDC::LinkStubs() ...@@ -542,9 +547,10 @@ void THcDC::LinkStubs()
// Now make sure this sp is not already used in a sp. // Now make sure this sp is not already used in a sp.
// Could this be done by having a sp point to the track it is in? // Could this be done by having a sp point to the track it is in?
Int_t tryflag=1; Int_t tryflag=1;
for(Int_t itrack=0;itrack<ntracks_fp;itrack++) { for(Int_t itrack=0;itrack<fNDCTracks;itrack++) {
for(Int_t isp=0;isp<fTrackSP[itrack].nSP;isp++) { THcDCTrack *theDCTrack = static_cast<THcDCTrack*>( fDCTracks->At(itrack));
if(fTrackSP[itrack].spID[isp] == isp1) { for(Int_t isp=0;isp<theDCTrack->GetNSpacePoints();isp++) {
if(theDCTrack->GetSpacePointID(isp) == isp1) {
tryflag=0; tryflag=0;
} }
} }
...@@ -580,22 +586,21 @@ void THcDC::LinkStubs() ...@@ -580,22 +586,21 @@ void THcDC::LinkStubs()
assert(sptracks==0); assert(sptracks==0);
//stubtest=1; Used in h_track_tests.f //stubtest=1; Used in h_track_tests.f
// Make a new track if there are not to many // Make a new track if there are not to many
if(ntracks_fp < MAXTRACKS) { if(fNDCTracks < MAXTRACKS) {
sptracks=0; // Number of tracks with this seed sptracks=0; // Number of tracks with this seed
stub_tracks[sptracks++] = ntracks_fp; stub_tracks[sptracks++] = fNDCTracks;
fTrackSP[ntracks_fp].nSP=2; THcDCTrack *theDCTrack = new( (*fDCTracks)[fNDCTracks++]) THcDCTrack(fNPlanes);
fTrackSP[ntracks_fp].spID[0] = isp1; theDCTrack->AddSpacePoint(isp1);
fTrackSP[ntracks_fp].spID[1] = isp2; theDCTrack->AddSpacePoint(isp2);
// Now save the X, Y and XP for the two stubs // Now save the X, Y and XP for the two stubs
// in arrays hx_sp1, hy_sp1, hy_sp1, ... hxp_sp2 // in arrays hx_sp1, hy_sp1, hy_sp1, ... hxp_sp2
// Why not also YP? // Why not also YP?
// Skip for here. May be a diagnostic thing // Skip for here. May be a diagnostic thing
newtrack = 0; // Make no more tracks in this loop newtrack = 0; // Make no more tracks in this loop
// (But could replace a SP?) // (But could replace a SP?)
ntracks_fp++;
} else { } else {
if (fDebugDC) cout << "EPIC FAIL 1: Too many tracks found in THcDC::LinkStubs" << endl; if (fDebugDC) cout << "EPIC FAIL 1: Too many tracks found in THcDC::LinkStubs" << endl;
ntracks_fp=0; fNDCTracks=0;
// Do something here to fail this event // Do something here to fail this event
return; return;
} }
...@@ -604,14 +609,16 @@ void THcDC::LinkStubs() ...@@ -604,14 +609,16 @@ void THcDC::LinkStubs()
// Check if there is another space point in the same chamber // Check if there is another space point in the same chamber
for(Int_t itrack=0;itrack<sptracks;itrack++) { for(Int_t itrack=0;itrack<sptracks;itrack++) {
Int_t track=stub_tracks[itrack]; Int_t track=stub_tracks[itrack];
THcDCTrack *theDCTrack = static_cast<THcDCTrack*>( fDCTracks->At(track));
Int_t spoint=0; Int_t spoint=0;
Int_t duppoint=0; Int_t duppoint=0;
for(Int_t isp=0;isp<fTrackSP[track].nSP;isp++) { for(Int_t isp=0;isp<theDCTrack->GetNSpacePoints();isp++) {
if(fSp[isp2]->fNChamber == if(fSp[isp2]->fNChamber ==
fSp[fTrackSP[track].spID[isp]]->fNChamber) { fSp[theDCTrack->GetSpacePointID(isp)]->fNChamber) {
spoint=isp; spoint=isp;
} }
if(isp2==fTrackSP[track].spID[isp]) { if(isp2==theDCTrack->GetSpacePointID(isp)) {
duppoint=1; duppoint=1;
} }
} // End loop over sp in tracks with isp1 } // End loop over sp in tracks with isp1
...@@ -619,24 +626,24 @@ void THcDC::LinkStubs() ...@@ -619,24 +626,24 @@ void THcDC::LinkStubs()
// add this space point to current track(2) // add this space point to current track(2)
if(!duppoint) { if(!duppoint) {
if(!spoint) { if(!spoint) {
fTrackSP[track].spID[fTrackSP[track].nSP++] = isp2; theDCTrack->AddSpacePoint(isp2);
} else { } else {
// If there is another point in the same chamber // If there is another point in the same chamber
// in this track create a new track with all the // in this track create a new track with all the
// same space points except spoint // same space points except spoint
if(ntracks_fp < MAXTRACKS) { if(fNDCTracks < MAXTRACKS) {
stub_tracks[sptracks++] = ntracks_fp; stub_tracks[sptracks++] = fNDCTracks;
fTrackSP[ntracks_fp].nSP=fTrackSP[track].nSP; THcDCTrack *newDCTrack = new( (*fDCTracks)[fNDCTracks++]) THcDCTrack(fNPlanes);
for(Int_t isp=0;isp<fTrackSP[track].nSP;isp++) { for(Int_t isp=0;isp<theDCTrack->GetNSpacePoints();isp++) {
if(isp!=spoint) { if(isp!=spoint) {
fTrackSP[ntracks_fp].spID[isp] = fTrackSP[track].spID[isp]; newDCTrack->AddSpacePoint(theDCTrack->GetSpacePointID(isp));
} else { } else {
fTrackSP[ntracks_fp].spID[isp] = isp2; newDCTrack->AddSpacePoint(theDCTrack->GetSpacePointID(isp2));
} // End check for dup on copy } // End check for dup on copy
} // End copy of track } // End copy of track
} else { } else {
if (fDebugDC) cout << "EPIC FAIL 2: Too many tracks found in THcDC::LinkStubs" << endl; if (fDebugDC) cout << "EPIC FAIL 2: Too many tracks found in THcDC::LinkStubs" << endl;
ntracks_fp=0; fNDCTracks=0;
// Do something here to fail this event // Do something here to fail this event
return; // Max # of allowed tracks return; // Max # of allowed tracks
} }
...@@ -649,34 +656,36 @@ void THcDC::LinkStubs() ...@@ -649,34 +656,36 @@ void THcDC::LinkStubs()
} // end isp2 loop over new space points } // end isp2 loop over new space points
} // end test on tryflag } // end test on tryflag
} // end isp1 outer loop over space points } // end isp1 outer loop over space points
} else { // Make track out of each single space point } else { // Make track out of each single space point
for(Int_t isp=0;isp<fNSp;isp++) { for(Int_t isp=0;isp<fNSp;isp++) {
if(ntracks_fp<MAXTRACKS) { if(fNDCTracks<MAXTRACKS) {
fTrackSP[ntracks_fp].nSP=1; // Need some constructed at thingy
fTrackSP[ntracks_fp].spID[0]=isp; THcDCTrack *newDCTrack = new( (*fDCTracks)[fNDCTracks++]) THcDCTrack(fNPlanes);
ntracks_fp++; newDCTrack->AddSpacePoint(isp);
} else { } else {
if (fDebugDC) cout << "EPIC FAIL 3: Too many tracks found in THcDC::LinkStubs" << endl; if (fDebugDC) cout << "EPIC FAIL 3: Too many tracks found in THcDC::LinkStubs" << endl;
ntracks_fp=0; fNDCTracks=0;
// Do something here to fail this event // Do something here to fail this event
return; // Max # of allowed tracks return; // Max # of allowed tracks
} }
} }
} }
// Add the list of hits on the track to the track. // Add the list of hits on the track to the track.
for(Int_t itrack=0;itrack<ntracks_fp;itrack++) { // Looks like it adds all hits for all space points to every track
fTrackSP[itrack].fNHits=0; for(Int_t itrack=0;itrack<fNDCTracks;itrack++) {
fTrackSP[itrack].fHits.clear(); THcDCTrack *theDCTrack = static_cast<THcDCTrack*>( fDCTracks->At(itrack));
for(Int_t isp=0;isp<fNSp;isp++) { theDCTrack->Clear();
for(Int_t ihit=0;ihit<fSp[isp]->GetNHits();ihit++) { // Hit list in the track should have been cleared when created.
fTrackSP[itrack].fHits.push_back(fSp[isp]->GetHit(ihit)); for(Int_t isp=0;isp<theDCTrack->GetNSpacePoints();isp++) {
fTrackSP[itrack].fNHits++; Int_t spid=theDCTrack->GetSpacePointID(isp);
for(Int_t ihit=0;ihit<fSp[spid]->GetNHits();ihit++) {
theDCTrack->AddHit(fSp[spid]->GetHit(ihit));
} }
} }
} }
/// ///
/// ///
if (fDebugDC) cout << "Found " << ntracks_fp << " tracks"<<endl; if (fDebugDC) cout << "Found " << fNDCTracks << " tracks"<<endl;
} }
// Primary track fitting routine // Primary track fitting routine
...@@ -689,31 +698,32 @@ void THcDC::TrackFit() ...@@ -689,31 +698,32 @@ void THcDC::TrackFit()
// Initialize residuals // Initialize residuals
// Need to make these member variables so they can be histogrammed // Need to make these member variables so they can be histogrammed
// Probably an array of vectors. // Probably an array of vectors.
Double_t double_resolution[fNPlanes][ntracks_fp]; Double_t double_resolution[fNPlanes][fNDCTracks];
Double_t single_resolution[fNPlanes][ntracks_fp]; Double_t single_resolution[fNPlanes][fNDCTracks];
Double_t double_res[fNPlanes]; // For the good track Double_t double_res[fNPlanes]; // For the good track
for(Int_t ip=0;ip<fNPlanes;ip++) { for(Int_t ip=0;ip<fNPlanes;ip++) {
double_res[ip] = 1000.0; double_res[ip] = 1000.0;
for(Int_t itrack=0;itrack<ntracks_fp;itrack++) { for(Int_t itrack=0;itrack<fNDCTracks;itrack++) {
double_resolution[ip][itrack] = 1000.0; double_resolution[ip][itrack] = 1000.0;
single_resolution[ip][itrack] = 1000.0; single_resolution[ip][itrack] = 1000.0;
} }
} }
Double_t dummychi2 = 1.0E4; Double_t dummychi2 = 1.0E4;
for(Int_t itrack=0;itrack<ntracks_fp;itrack++) { for(Int_t itrack=0;itrack<fNDCTracks;itrack++) {
// Double_t chi2 = dummychi2; // Double_t chi2 = dummychi2;
// Int_t htrack_fit_num = itrack; // Int_t htrack_fit_num = itrack;
fTrackSP[itrack].fNfree = fTrackSP[itrack].fNHits - NUM_FPRAY; THcDCTrack *theDCTrack = static_cast<THcDCTrack*>( fDCTracks->At(itrack));
theDCTrack->SetNFree(theDCTrack->GetNHits() - NUM_FPRAY);
Double_t chi2 = dummychi2; Double_t chi2 = dummychi2;
if(fTrackSP[itrack].fNfree > 0) { if(theDCTrack->GetNFree() > 0) {
TVectorD TT(NUM_FPRAY); TVectorD TT(NUM_FPRAY);
TMatrixD AA(NUM_FPRAY,NUM_FPRAY); TMatrixD AA(NUM_FPRAY,NUM_FPRAY);
for(Int_t irayp=0;irayp<NUM_FPRAY;irayp++) { for(Int_t irayp=0;irayp<NUM_FPRAY;irayp++) {
TT[irayp] = 0.0; TT[irayp] = 0.0;
for(Int_t ihit=0;ihit < fTrackSP[itrack].fNHits;ihit++) { for(Int_t ihit=0;ihit < theDCTrack->GetNHits();ihit++) {
THcDCHit* hit=fTrackSP[itrack].fHits[ihit]; THcDCHit* hit=theDCTrack->GetHit(ihit);
Int_t plane=hit->GetPlaneNum()-1; Int_t plane=hit->GetPlaneNum()-1;
TT[irayp] += (hit->GetCoord()* TT[irayp] += (hit->GetCoord()*
fPlaneCoeffs[plane][raycoeffmap[irayp]]) fPlaneCoeffs[plane][raycoeffmap[irayp]])
...@@ -726,8 +736,8 @@ void THcDC::TrackFit() ...@@ -726,8 +736,8 @@ void THcDC::TrackFit()
if(jrayp<irayp) { // Symmetric if(jrayp<irayp) { // Symmetric
AA[irayp][jrayp] = AA[jrayp][irayp]; AA[irayp][jrayp] = AA[jrayp][irayp];
} else { } else {
for(Int_t ihit=0;ihit < fTrackSP[itrack].fNHits;ihit++) { for(Int_t ihit=0;ihit < theDCTrack->GetNHits();ihit++) {
THcDCHit* hit=fTrackSP[itrack].fHits[ihit]; THcDCHit* hit=theDCTrack->GetHit(ihit);
Int_t plane=hit->GetPlaneNum()-1; Int_t plane=hit->GetPlaneNum()-1;
AA[irayp][jrayp] += fPlaneCoeffs[plane][raycoeffmap[irayp]]* AA[irayp][jrayp] += fPlaneCoeffs[plane][raycoeffmap[irayp]]*
fPlaneCoeffs[plane][raycoeffmap[jrayp]]/ fPlaneCoeffs[plane][raycoeffmap[jrayp]]/
...@@ -748,42 +758,35 @@ void THcDC::TrackFit() ...@@ -748,42 +758,35 @@ void THcDC::TrackFit()
// } // }
// Calculate hit coordinate for each plane for chi2 and efficiency // Calculate hit coordinate for each plane for chi2 and efficiency
// calculations // calculations
fTrackSP[itrack].fCoords.clear();
fTrackSP[itrack].fResiduals.clear(); // Make sure fCoords, fResiduals, and fDoubleResiduals are clear
fTrackSP[itrack].fDoubleResiduals.clear();
for(Int_t iplane=0;iplane < fNPlanes; iplane++) { for(Int_t iplane=0;iplane < fNPlanes; iplane++) {
fTrackSP[itrack].fCoords.push_back(0.0); Double_t coord=0.0;
fTrackSP[itrack].fResiduals.push_back(1000.0);
fTrackSP[itrack].fDoubleResiduals.push_back(1000.0);
for(Int_t ir=0;ir<NUM_FPRAY;ir++) { for(Int_t ir=0;ir<NUM_FPRAY;ir++) {
fTrackSP[itrack].fCoords[iplane] += fPlaneCoeffs[iplane][raycoeffmap[ir]]*dray[ir]; coord += fPlaneCoeffs[iplane][raycoeffmap[ir]]*dray[ir];
} }
theDCTrack->SetCoord(iplane,coord);
} }
// Compute Chi2 and residuals // Compute Chi2 and residuals
chi2 = 0.0; chi2 = 0.0;
fTrackSP[itrack].fResiduals.clear(); for(Int_t ihit=0;ihit < theDCTrack->GetNHits();ihit++) {
for(Int_t ihit=0;ihit < fTrackSP[itrack].fNHits;ihit++) { THcDCHit* hit=theDCTrack->GetHit(ihit);
THcDCHit* hit=fTrackSP[itrack].fHits[ihit];
Int_t plane=hit->GetPlaneNum()-1; Int_t plane=hit->GetPlaneNum()-1;
Double_t residual = hit->GetCoord() - fTrackSP[itrack].fCoords[plane]; Double_t residual = hit->GetCoord() - theDCTrack->GetCoord(plane);
fTrackSP[itrack].fResiduals[plane] = residual; theDCTrack->SetResidual(plane, residual);
chi2 += pow(residual/fSigma[plane],2); chi2 += pow(residual/fSigma[plane],2);
} }
if (fDebugDC) { if (fDebugDC) {
cout << "Residuals:" << endl; cout << "Residuals:" << endl;
for(Int_t ihit=0;ihit < fTrackSP[itrack].fNHits;ihit++) { for(Int_t ihit=0;ihit < theDCTrack->GetNHits();ihit++) {
THcDCHit* hit=fTrackSP[itrack].fHits[ihit]; THcDCHit* hit=theDCTrack->GetHit(ihit);
Int_t plane=hit->GetPlaneNum()-1; Int_t plane=hit->GetPlaneNum()-1;
cout << " " << plane << ": " << fTrackSP[itrack].fResiduals[ihit] << endl; cout << " " << plane << ": " << theDCTrack->GetResidual(plane) << endl;
} }
} }
fTrackSP[itrack].x_fp = dray[0]; theDCTrack->SetVector(dray[0], dray[1], 0.0, dray[2], dray[3]);
fTrackSP[itrack].y_fp = dray[1];
fTrackSP[itrack].z_fp = 0.0;
fTrackSP[itrack].xp_fp = dray[2];
fTrackSP[itrack].xp_fp = dray[3];
} }
fTrackSP[itrack].chi2_fp = chi2; theDCTrack->SetChisq(chi2);
} }
// Calculate residuals for each chamber if in single stub mode // Calculate residuals for each chamber if in single stub mode
...@@ -791,46 +794,42 @@ void THcDC::TrackFit() ...@@ -791,46 +794,42 @@ void THcDC::TrackFit()
// Specific for two chambers. Can/should it be generalized? // Specific for two chambers. Can/should it be generalized?
if(fSingleStub != 0) { if(fSingleStub != 0) {
if(ntracks_fp == 2) { if(fNDCTracks == 2) {
THcDCTrack *theDCTrack1 = static_cast<THcDCTrack*>( fDCTracks->At(0));
THcDCTrack *theDCTrack2 = static_cast<THcDCTrack*>( fDCTracks->At(1));
Int_t itrack=0; Int_t itrack=0;
Int_t ihit=0; Int_t ihit=0;
THcDCHit* hit=fTrackSP[itrack].fHits[ihit]; THcDCHit* hit=theDCTrack1->GetHit(ihit);
Int_t plane=hit->GetPlaneNum()-1; Int_t plane=hit->GetPlaneNum()-1;
Int_t chamber=fNChamber[plane]; Int_t chamber=fNChamber[plane];
if(chamber==1) { if(chamber==1) {
itrack=1; itrack=1;
hit=fTrackSP[itrack].fHits[ihit]; hit=theDCTrack2->GetHit(ihit);
plane=hit->GetPlaneNum()-1; plane=hit->GetPlaneNum()-1;
chamber=fNChamber[plane]; chamber=fNChamber[plane];
if(chamber==2) { if(chamber==2) {
Double_t ray1[4]; Double_t ray1[4];
Double_t ray2[4]; Double_t ray2[4];
ray1[0] = fTrackSP[0].x_fp; theDCTrack1->GetRay(ray1);
ray1[1] = fTrackSP[0].y_fp; theDCTrack2->GetRay(ray2);
ray1[2] = fTrackSP[0].xp_fp;
ray1[3] = fTrackSP[0].yp_fp;
ray2[0] = fTrackSP[1].x_fp;
ray2[1] = fTrackSP[1].y_fp;
ray2[2] = fTrackSP[1].xp_fp;
ray2[3] = fTrackSP[1].yp_fp;
itrack = 1; itrack = 1;
// Loop over hits in second chamber // Loop over hits in second chamber
for(Int_t ihit=0;ihit < fTrackSP[itrack].fNHits;ihit++) { for(Int_t ihit=0;ihit < theDCTrack2->GetNHits();ihit++) {
// Calculate residual in second chamber from first chamber track // Calculate residual in second chamber from first chamber track
THcDCHit* hit=fTrackSP[itrack].fHits[ihit]; THcDCHit* hit=theDCTrack2->GetHit(ihit);
Int_t plane=hit->GetPlaneNum()-1; Int_t plane=hit->GetPlaneNum()-1;
Double_t pos = DpsiFun(ray1,plane); Double_t pos = DpsiFun(ray1,plane);
fTrackSP[itrack-1].fDoubleResiduals[plane] = hit->GetCoord() - pos; theDCTrack1->SetDoubleResidual(plane,hit->GetCoord() - pos);
// hdc_dbl_res(pln) = hdc_double_residual(1,pln) for hists // hdc_dbl_res(pln) = hdc_double_residual(1,pln) for hists
} }
itrack=0; itrack=0;
// Loop over hits in first chamber // Loop over hits in first chamber
for(Int_t ihit=0;ihit < fTrackSP[itrack].fNHits;ihit++) { for(Int_t ihit=0;ihit < theDCTrack1->GetNHits();ihit++) {
// Calculate residual in first chamber from second chamber track // Calculate residual in first chamber from second chamber track
THcDCHit* hit=fTrackSP[itrack].fHits[ihit]; THcDCHit* hit=theDCTrack1->GetHit(ihit);
Int_t plane=hit->GetPlaneNum()-1; Int_t plane=hit->GetPlaneNum()-1;
Double_t pos = DpsiFun(ray1,plane); Double_t pos = DpsiFun(ray1,plane);
fTrackSP[itrack+1].fDoubleResiduals[plane] = hit->GetCoord() - pos; theDCTrack2->SetDoubleResidual(plane,hit->GetCoord() - pos);
// hdc_dbl_res(pln) = hdc_double_residual(1,pln) for hists // hdc_dbl_res(pln) = hdc_double_residual(1,pln) for hists
} }
} }
......
...@@ -75,6 +75,9 @@ public: ...@@ -75,6 +75,9 @@ public:
THcDC(); // for ROOT I/O THcDC(); // for ROOT I/O
protected: protected:
Int_t fDebugDC; Int_t fDebugDC;
Int_t fNDCTracks;
TClonesArray* fDCTracks; // Tracks found from stubs (THcDCTrack obj)
// Calibration // Calibration
// Per-event data // Per-event data
...@@ -129,19 +132,6 @@ protected: ...@@ -129,19 +132,6 @@ protected:
// Intermediate structure for building // Intermediate structure for building
static const char MAXTRACKS = 50; static const char MAXTRACKS = 50;
struct TrackSP {
Int_t nSP; /* Number of space points in this track */
Int_t spID[10]; /* List of space points in this track */
Int_t fNHits;
Int_t fNfree; /* Number of degrees of freedom */
std::vector<THcDCHit*> fHits; /* List of hits for this track */
std::vector<Double_t> fCoords; /* Coordinate on each plane */
std::vector<Double_t> fResiduals; /* Residual on each plane */
std::vector<Double_t> fDoubleResiduals; /* Residual on each plane for single stub mode */
Double_t x_fp, y_fp, z_fp;
Double_t xp_fp, yp_fp;
Double_t chi2_fp;
} fTrackSP[MAXTRACKS];
std::vector<THcDriftChamberPlane*> fPlanes; // List of plane objects std::vector<THcDriftChamberPlane*> fPlanes; // List of plane objects
std::vector<THcDriftChamber*> fChambers; // List of chamber objects std::vector<THcDriftChamber*> fChambers; // List of chamber objects
......
///////////////////////////////////////////////////////////////////////////////
// //
// THcDCTrack //
// //
// Class representing a track found from linking DC Space points //
///////////////////////////////////////////////////////////////////////////////
#include "THcDCHit.h"
#include "THcDCTrack.h"
THcDCTrack::THcDCTrack(Int_t nplanes) : fnSP(0), fNHits(0)
{
fHits.clear();
fCoords.resize(nplanes);
fResiduals.resize(nplanes);
fDoubleResiduals.resize(nplanes);
}
void THcDCTrack::AddHit(THcDCHit * hit)
{
// Add a hit to the track
fHits.push_back(hit);
fNHits++;
}
void THcDCTrack::AddSpacePoint( Int_t spid )
{
// Add to list of space points in this track
// Need a check for maximum spacepoints of 10
fspID[fnSP++] = spid;
}
void THcDCTrack::Clear( const Option_t* )
{
// Clear the space point and hit lists
fnSP = 0;
ClearHits();
// Need to set default values (0 or -100)
//fCoords.clear();
//fResiduals.clear();
//fDoubleResiduals.clear();
}
void THcDCTrack::ClearHits( )
{
fNHits = 0;
fHits.clear();
}
ClassImp(THcDCTrack)
///////////////////////////////////////////////////////////////////////////////
#ifndef ROOT_THcDCTrack
#define ROOT_THcDCTrack
///////////////////////////////////////////////////////////////////////////////
// //
// THcDCTrack //
// //
///////////////////////////////////////////////////////////////////////////////
//#include "THaCluster.h"
#include "TVector3.h"
#include "TObject.h"
//class THaVDCCluster;
class THcDCPlane;
class THaTrack;
class THcDCHit;
class THcDCTrack : public TObject {
public:
THcDCTrack(Int_t nplanes);
virtual ~THcDCTrack() {};
virtual void AddHit(THcDCHit * hit);
virtual void AddSpacePoint(Int_t spid);
//Get and Set Functions
Int_t* GetSpacePoints() {return fspID;}
Int_t GetNSpacePoints() const { return fnSP;}
Int_t GetSpacePointID(Int_t i) const {return fspID[i];}
THcDCHit* GetHit(Int_t i) const {return fHits[i];}
Int_t GetNHits() const {return fNHits;}
Int_t GetNFree() const {return fNfree;}
Double_t GetCoord(Int_t ip) const {return fCoords[ip];}
Double_t GetResidual(Int_t ip) const {return fResiduals[ip];}
void GetRay(Double_t *ray) const {ray[0]=fX_fp; ray[1]=fY_fp; ray[2]=fXp_fp; ray[3]=fYp_fp;}
void SetNFree(Int_t nfree) {fNfree = nfree;}
void SetCoord(Int_t ip, Double_t coord) {fCoords[ip] = coord;}
void SetResidual(Int_t ip, Double_t coord) {fResiduals[ip] = coord;}
void SetDoubleResidual(Int_t ip, Double_t coord) {fDoubleResiduals[ip] = coord;}
void SetVector(Double_t x, Double_t y, Double_t z,
Double_t xp, Double_t yp) {fX_fp = x; fY_fp=y; fZ_fp=z; fXp_fp=xp; fYp_fp=yp;}
void SetChisq(Double_t chi2) {fChi2_fp = chi2;}
virtual void ClearHits( );
// TObject functions redefined
virtual void Clear( Option_t* opt="" );
protected:
Int_t fnSP; /* Number of space points in this track */
/* Should we have a list of pointers to space points
instead of their indices? */
Int_t fspID[10]; /* List of space points in this track */
Int_t fNHits;
Int_t fNfree; /* Number of degrees of freedom */
std::vector<THcDCHit*> fHits; /* List of hits for this track */
std::vector<Double_t> fCoords; /* Coordinate on each plane */
std::vector<Double_t> fResiduals; /* Residual on each plane */
std::vector<Double_t> fDoubleResiduals; /* Residual on each plane for single stub mode */
Double_t fX_fp, fY_fp, fZ_fp;
Double_t fXp_fp, fYp_fp;
Double_t fChi2_fp;
private:
// Hide copy ctor and op=
THcDCTrack( const THcDCTrack& );
THcDCTrack& operator=( const THcDCTrack& );
ClassDef(THcDCTrack,0) // DCTrack class
};
#endif
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment