From 58c45eeeb61a5f6739213c7a3a17ab3f4aac5dc0 Mon Sep 17 00:00:00 2001
From: "Stephen A. Wood" <saw@jlab.org>
Date: Fri, 8 Mar 2013 16:47:08 -0500
Subject: [PATCH] Add computation of drift distance from drift time. It has
not been check that we are getting a good start time yet. Drift time and
distance added to tree Change DC plane names from 1, 2, 3, ... to 1x1, 1y1,
... This is so that the parameters holding the time to distance maps
can be found. (The parameters are e.g. hwc1x1fract) Need to find a way
not to have to hard code these plane names. Either use wire angles
(alpha) or some kind of parameter name mapping Changed output.def to match
new plane names
---
examples/PARAM/general.param | 2 +-
examples/output.def | 24 ++++-----
src/THcDCHit.cxx | 4 +-
src/THcDCHit.h | 6 ++-
src/THcDCLookupTTDConv.cxx | 95 +++++++++++-------------------------
src/THcDCLookupTTDConv.h | 26 +++-------
src/THcDCTimeToDistConv.h | 4 +-
src/THcDriftChamber.cxx | 41 ++++++++++++++--
src/THcDriftChamberPlane.cxx | 21 +++++++-
9 files changed, 115 insertions(+), 108 deletions(-)
diff --git a/examples/PARAM/general.param b/examples/PARAM/general.param
index 898d4b3..c327a18 100644
--- a/examples/PARAM/general.param
+++ b/examples/PARAM/general.param
@@ -25,4 +25,4 @@ raddeg=3.14159265/180
#include "PARAM/hhodo.param"
#include "PARAM/haero.param"
#include "PARAM/hdc.param"
-
+#include "PARAM/hdriftmap.param"
diff --git a/examples/output.def b/examples/output.def
index c351efd..7cfb0da 100644
--- a/examples/output.def
+++ b/examples/output.def
@@ -53,15 +53,15 @@ TH1F calposadc1 'HMS Cal ADC1' H.Cal.1pr.posadc1 150 -50 400
# Can we use variables for the constants. In CTP we used hdc_nwire(i)
#
-TH1F hdc1x1_wm 'HDC 1X1 Wiremap' H.dc.1.tdchits 113 0.5 113.5
-TH1F hdc1y1_wm 'HDC 1Y1 Wiremap' H.dc.2.tdchits 52 0.5 52.5
-TH1F hdc1u1_wm 'HDC 1U1 Wiremap' H.dc.3.tdchits 107 0.5 107.5
-TH1F hdc1v1_wm 'HDC 1V1 Wiremap' H.dc.4.tdchits 107 0.5 107.5
-TH1F hdc1y2_wm 'HDC 1Y2 Wiremap' H.dc.5.tdchits 52 0.5 52.5
-TH1F hdc1x2_wm 'HDC 1X2 Wiremap' H.dc.6.tdchits 113 0.5 113.5
-TH1F hdc2x1_wm 'HDC 2X1 Wiremap' H.dc.7.tdchits 113 0.5 113.5
-TH1F hdc2y1_wm 'HDC 2Y1 Wiremap' H.dc.8.tdchits 52 0.5 52.5
-TH1F hdc2u1_wm 'HDC 2U1 Wiremap' H.dc.9.tdchits 107 0.5 107.5
-TH1F hdc2v1_wm 'HDC 2V1 Wiremap' H.dc.10.tdchits 107 0.5 107.5
-TH1F hdc2y2_wm 'HDC 2Y2 Wiremap' H.dc.11.tdchits 52 0.5 52.5
-TH1F hdc2x2_wm 'HDC 2X2 Wiremap' H.dc.12.tdchits 113 0.5 113.5
+TH1F hdc1x1_wm 'HDC 1X1 Wiremap' H.dc.1x1.tdchits 113 0.5 113.5
+TH1F hdc1y1_wm 'HDC 1Y1 Wiremap' H.dc.1y1.tdchits 52 0.5 52.5
+TH1F hdc1u1_wm 'HDC 1U1 Wiremap' H.dc.1u1.tdchits 107 0.5 107.5
+TH1F hdc1v1_wm 'HDC 1V1 Wiremap' H.dc.1v1.tdchits 107 0.5 107.5
+TH1F hdc1y2_wm 'HDC 1Y2 Wiremap' H.dc.1y2.tdchits 52 0.5 52.5
+TH1F hdc1x2_wm 'HDC 1X2 Wiremap' H.dc.1x2.tdchits 113 0.5 113.5
+TH1F hdc2x1_wm 'HDC 2X1 Wiremap' H.dc.2x1.tdchits 113 0.5 113.5
+TH1F hdc2y1_wm 'HDC 2Y1 Wiremap' H.dc.2y1.tdchits 52 0.5 52.5
+TH1F hdc2u1_wm 'HDC 2U1 Wiremap' H.dc.2u1.tdchits 107 0.5 107.5
+TH1F hdc2v1_wm 'HDC 2V1 Wiremap' H.dc.2v1.tdchits 107 0.5 107.5
+TH1F hdc2y2_wm 'HDC 2Y2 Wiremap' H.dc.2y2.tdchits 52 0.5 52.5
+TH1F hdc2x2_wm 'HDC 2X2 Wiremap' H.dc.2x2.tdchits 113 0.5 113.5
diff --git a/src/THcDCHit.cxx b/src/THcDCHit.cxx
index d8ad011..0484b1a 100644
--- a/src/THcDCHit.cxx
+++ b/src/THcDCHit.cxx
@@ -14,7 +14,7 @@ ClassImp(THcDCHit)
const Double_t THcDCHit::kBig = 1.e38; // Arbitrary large value
//_____________________________________________________________________________
-Double_t THcDCHit::ConvertTimeToDist(Double_t slope)
+Double_t THcDCHit::ConvertTimeToDist()
{
// Converts TDC time to drift distance
// Takes the (estimated) slope of the track as an argument
@@ -24,7 +24,7 @@ Double_t THcDCHit::ConvertTimeToDist(Double_t slope)
if (ttdConv) {
// If a time to distance algorithm exists, use it to convert the TDC time
// to the drift distance
- fDist = ttdConv->ConvertTimeToDist(fTime, slope, &fdDist);
+ fDist = ttdConv->ConvertTimeToDist(fTime);
return fDist;
}
diff --git a/src/THcDCHit.h b/src/THcDCHit.h
index 3e5150e..79f88e7 100644
--- a/src/THcDCHit.h
+++ b/src/THcDCHit.h
@@ -15,10 +15,12 @@ class THcDCHit : public TObject {
public:
THcDCHit( THcDCWire* wire=NULL, Int_t rawtime=0, Double_t time=0.0 ) :
- fWire(wire), fRawTime(rawtime), fTime(time), fDist(0.0), ftrDist(kBig) {}
+ fWire(wire), fRawTime(rawtime), fTime(time), fDist(0.0), ftrDist(kBig) {
+ ConvertTimeToDist();
+ }
virtual ~THcDCHit() {}
- virtual Double_t ConvertTimeToDist(Double_t slope);
+ virtual Double_t ConvertTimeToDist();
Int_t Compare ( const TObject* obj ) const;
Bool_t IsSortable () const { return kTRUE; }
diff --git a/src/THcDCLookupTTDConv.cxx b/src/THcDCLookupTTDConv.cxx
index 0485c1a..7665457 100644
--- a/src/THcDCLookupTTDConv.cxx
+++ b/src/THcDCLookupTTDConv.cxx
@@ -1,6 +1,9 @@
///////////////////////////////////////////////////////////////////////////////
// //
-// THcDCLookupTTDConv //
+// THcDCLookupTTDConv //
+// //
+// Upon initialization needs to be given the lookup table for time //
+// to distance conversion. //
// //
///////////////////////////////////////////////////////////////////////////////
@@ -10,32 +13,21 @@ ClassImp(THcDCLookupTTDConv)
//______________________________________________________________________________
-THcDCLookupTTDConv::THcDCLookupTTDConv()
+THcDCLookupTTDConv::THcDCLookupTTDConv(Double_t T0, Double_t MaxDriftDistance,
+ Double_t BinSize, Int_t NumBins,
+ Double_t* Table) :
+fT0(T0), fMaxDriftDistance(MaxDriftDistance), fBinSize(BinSize),
+ fNumBins(NumBins)
{
//Normal constructor
-}
-//______________________________________________________________________________
-THcDCLookupTTDConv::THcDCLookupTTDConv( Double_t vel)
-{
- // Normal constructor
- fDriftVel = vel;
+ fTable = new Double_t[fNumBins];
+ for(Int_t i=0;i<fNumBins;i++) {
+ fTable[i] = Table[i];
+ }
- // TODO: This should be read from database!!
- fA1tdcCor[0] = 2.12e-3;
- fA1tdcCor[1] = 0.0;
- fA1tdcCor[2] = 0.0;
- fA1tdcCor[3] = 0.0;
- fA2tdcCor[0] = -4.20e-4;
- fA2tdcCor[1] = 1.3e-3;
- fA2tdcCor[2] = 1.06e-4;
- fA2tdcCor[3] = 0.0;
-
- fdtime = 4.e-9; // 4ns -> 200 microns
}
-
-
//______________________________________________________________________________
THcDCLookupTTDConv::~THcDCLookupTTDConv()
{
@@ -44,54 +36,25 @@ THcDCLookupTTDConv::~THcDCLookupTTDConv()
}
//______________________________________________________________________________
-Double_t THcDCLookupTTDConv::ConvertTimeToDist(Double_t time,
- Double_t tanTheta,
- Double_t *ddist)
+Double_t THcDCLookupTTDConv::ConvertTimeToDist(Double_t time)
{
- // Drift Velocity in m/s
- // time in s
- // Return m
-
-// printf("Converting Drift Time to Drift Distance!\n");
-
- Double_t a1 = 0.0, a2 = 0.0;
- // Find the values of a1 and a2 by evaluating the proper polynomials
- // a = A_3 * x^3 + A_2 * x^2 + A_1 * x + A_0
-
- tanTheta = 1.0 / tanTheta; // I assume this has to do w/ making the
- // polynomial have the proper variable...
-
- for (Int_t i = 3; i >= 1; i--) {
- a1 = tanTheta * (a1 + fA1tdcCor[i]);
- a2 = tanTheta * (a2 + fA2tdcCor[i]);
- }
- a1 += fA1tdcCor[0];
- a2 += fA2tdcCor[0];
-
-// printf("a1(%e) = %e\n", tanTheta, a1);
-// printf("a2(%e) = %e\n", tanTheta, a2);
-
- // ESPACE software includes corrections to the time for
- // 1. Cluster t0 (offset applied to entire cluster)
- // 2. Time of flight to scintillators
- Double_t dist = fDriftVel * time;
- Double_t unc = fDriftVel * fdtime; // watch uncertainty in the timing
- if (dist < 0) {
- // something screwy is going on
- } else if (dist < a1 ) {
- // dist = fDriftVel * time * (1 + 1 / (a1/a2 + 1));
- dist *= ( 1 + a2 / a1);
- unc *= ( 1 + a2 / a1);
- } else {
- dist += a2;
+ // Lookup in table
+ Int_t ib = (time-fT0)/fBinSize;
+ Double_t frac = 0;
+ if(ib >= 0 && ib+1 < fNumBins) {
+ Double_t tfrac = (time - (ib*fBinSize + fT0)) / fBinSize;
+ frac = fTable[ib]*(1-tfrac) + fTable[ib+1]*tfrac;
+ } else if (ib+1 >= fNumBins) {
+ frac = 1.0;
}
-
- if (ddist) *ddist = unc;
-// printf("D(%e) = %e\nUncorrected D = %e\n", time, dist, fDriftVel * time);
-
- return dist;
-}
+ Double_t drift_distance = fMaxDriftDistance * frac;
+
+ // Engine subtracts a hdc_card_delay from this. Seems
+ // to be zero in the PARAM files, bit is it always? Delay implies
+ // time? Whis is a time subtracted from a distance?
+ return(drift_distance);
+}
////////////////////////////////////////////////////////////////////////////////
diff --git a/src/THcDCLookupTTDConv.h b/src/THcDCLookupTTDConv.h
index 3ae38a5..b0091ae 100644
--- a/src/THcDCLookupTTDConv.h
+++ b/src/THcDCLookupTTDConv.h
@@ -13,31 +13,21 @@
class THcDCLookupTTDConv : public THcDCTimeToDistConv{
public:
- THcDCLookupTTDConv( );
- THcDCLookupTTDConv(Double_t vel);
+ THcDCLookupTTDConv(Double_t T0, Double_t MaxDriftDistance, Double_t BinSize,
+ Int_t NumBins, Double_t* Table);
virtual ~THcDCLookupTTDConv();
- virtual Double_t ConvertTimeToDist(Double_t time, Double_t tanTheta,
- Double_t *ddist=0);
+ virtual Double_t ConvertTimeToDist(Double_t time);
- // Get and Set Functions
- Double_t GetDriftVel() { return fDriftVel; }
-
- void SetDriftVel(Double_t v) {fDriftVel = v; }
-
protected:
- Double_t fDriftVel; // Drift velocity (m / s)
-
- // Coefficients for a polynomial yielding correction parameters
- // For now, hard code these values from db_eh845
- // Eventually, this need to be read directly from the database
- Double_t fA1tdcCor[4];
- Double_t fA2tdcCor[4];
-
- Double_t fdtime; // uncertainty in the measured time
+ Double_t fT0;
+ Double_t fBinSize;
+ Int_t fNumBins;
+ Double_t fMaxDriftDistance;
+ Double_t* fTable;
ClassDef(THcDCLookupTTDConv,0) // VDC Analytic TTD Conv class
};
diff --git a/src/THcDCTimeToDistConv.h b/src/THcDCTimeToDistConv.h
index 2582249..960a5c3 100644
--- a/src/THcDCTimeToDistConv.h
+++ b/src/THcDCTimeToDistConv.h
@@ -17,8 +17,8 @@ public:
THcDCTimeToDistConv() {}
virtual ~THcDCTimeToDistConv();
- virtual Double_t ConvertTimeToDist(Double_t time, Double_t tanTheta,
- Double_t *ddist=0) = 0;
+ virtual Double_t ConvertTimeToDist(Double_t time) = 0;
+
private:
THcDCTimeToDistConv( const THcDCTimeToDistConv& );
diff --git a/src/THcDriftChamber.cxx b/src/THcDriftChamber.cxx
index 1c73a1e..ea45f06 100644
--- a/src/THcDriftChamber.cxx
+++ b/src/THcDriftChamber.cxx
@@ -73,13 +73,48 @@ void THcDriftChamber::Setup(const char* name, const char* description)
cout << "Drift Chambers: " << fNPlanes << " planes in " << fNChambers << " chambers" << endl;
fPlaneNames = new char* [fNPlanes];
-
+ for(Int_t i=0;i<fNPlanes;i++) {fPlaneNames[i] = new char[4];}
+
+ // Big hack needed because the drift map parameters are in parameters with
+ // names like hwc1x1fract. Need to do some kind of parameter name mapping
+ // or generate the names from the wire angle (alpha)
+ if(prefix[0] == 'h') {
+ strcpy(fPlaneNames[0],"1x1");
+ strcpy(fPlaneNames[1],"1y1");
+ strcpy(fPlaneNames[2],"1u1");
+ strcpy(fPlaneNames[3],"1v1");
+ strcpy(fPlaneNames[4],"1y2");
+ strcpy(fPlaneNames[5],"1x2");
+ strcpy(fPlaneNames[6],"2x1");
+ strcpy(fPlaneNames[7],"2y1");
+ strcpy(fPlaneNames[8],"2u1");
+ strcpy(fPlaneNames[9],"2v1");
+ strcpy(fPlaneNames[10],"2y2");
+ strcpy(fPlaneNames[11],"2x2");
+ } else if (prefix[0] == 's') {
+ strcpy(fPlaneNames[0],"1u1");
+ strcpy(fPlaneNames[1],"1u2");
+ strcpy(fPlaneNames[2],"1x1");
+ strcpy(fPlaneNames[3],"1x2");
+ strcpy(fPlaneNames[4],"1v1");
+ strcpy(fPlaneNames[5],"1v2");
+ strcpy(fPlaneNames[6],"2u1");
+ strcpy(fPlaneNames[7],"2u2");
+ strcpy(fPlaneNames[8],"2x1");
+ strcpy(fPlaneNames[9],"2x2");
+ strcpy(fPlaneNames[10],"2v1");
+ strcpy(fPlaneNames[11],"2v2");
+ } else {
+ cout << "Unknown Spectrometer Prefix '" << prefix << "' Guessing names" << endl;
+ for(Int_t i=0;i<fNPlanes;i++) {
+ sprintf(fPlaneNames[i],"%d",i+1);
+ }
+ }
+
char *desc = new char[strlen(description)+100];
fPlanes = new THcDriftChamberPlane* [fNPlanes];
for(Int_t i=0;i<fNPlanes;i++) {
- fPlaneNames[i] = new char[5];
- sprintf(fPlaneNames[i],"%d",i+1);
strcpy(desc, description);
strcat(desc, " Plane ");
strcat(desc, fPlaneNames[i]);
diff --git a/src/THcDriftChamberPlane.cxx b/src/THcDriftChamberPlane.cxx
index cae6563..816d27a 100644
--- a/src/THcDriftChamberPlane.cxx
+++ b/src/THcDriftChamberPlane.cxx
@@ -80,11 +80,23 @@ Int_t THcDriftChamberPlane::ReadDatabase( const TDatime& date )
static const char* const here = "ReadDatabase()";
char prefix[2];
char parname[100];
+ Int_t NumDriftMapBins;
+ Double_t DriftMapFirstBin;
+ Double_t DriftMapBinSize;
+ Double_t DriftMap[1000];
prefix[0]=tolower(GetParent()->GetPrefix()[0]);
prefix[1]='\0';
+ DBRequest list[]={
+ {"driftbins", &NumDriftMapBins, kInt},
+ {"drift1stbin", &DriftMapFirstBin, kDouble},
+ {"driftbinsz", &DriftMapBinSize, kDouble},
+ {Form("wc%sfract",GetName()),DriftMap,kDouble,1000},
+ {0}
+ };
+ gHcParms->LoadParmValues((DBRequest*)&list,prefix);
- // Retrieve parameters we need
+ // Retrieve parameters we need from parent class
THcDriftChamber* fParent;
fParent = (THcDriftChamber*) GetParent();
@@ -103,7 +115,8 @@ Int_t THcDriftChamberPlane::ReadDatabase( const TDatime& date )
cout << fPlaneNum << " " << fNWires << endl;
- fTTDConv = new THcDCLookupTTDConv();// Need to pass the lookup table
+ fTTDConv = new THcDCLookupTTDConv(DriftMapFirstBin,fPitch/2,DriftMapBinSize,
+ NumDriftMapBins,DriftMap);
Int_t nWires = fParent->GetNWires(fPlaneNum);
// For HMS, wire numbers start with one, but arrays start with zero.
@@ -142,6 +155,10 @@ Int_t THcDriftChamberPlane::DefineVariables( EMode mode )
"fHits.THcDCHit.GetWireNum()"},
{"rawtdc", "Raw TDC Values",
"fHits.THcDCHit.GetRawTime()"},
+ {"time","Drift times",
+ "fHits.THcDCHit.GetTime()"},
+ {"dist","Drift distancess",
+ "fHits.THcDCHit.GetDist()"},
{ 0 }
};
--
GitLab