diff --git a/examples/PARAM/general.param b/examples/PARAM/general.param index 898d4b3caeb8f21b64012dd0e0da900c920ca65d..c327a187a2e884bbfcd49a59f72417225575957f 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 c351efdbae4a030b810af1cefd78cb3530ec0651..7cfb0da2872f8e381dd3792f9ec7c4f24bd9412a 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 d8ad0110bd4a7fadd6763c861506e555fa292e49..0484b1ad39d1f2d1ee0171dfe0d2f5194f6a0ce0 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 3e5150e83521670b515c3d5e80d83ff31303bcac..79f88e7123f7a8d39d805719e7ce2b74ec831316 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 0485c1a308e0317f49f9fe17f426bc44757d766a..766545742f99b2742adf2818373980eb09d61570 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 3ae38a506f908be0eae7e6f9044d0762987ae933..b0091aeecff92a49d060c85b30f4c4238bf1ee15 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 2582249037c82987f176cca7eb13d04e9a5f5e33..960a5c373c8396a57c593c8657df78430d70bc2d 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 1c73a1e53bc82e0a276aba55f777002d7a2c79da..ea45f0678d544755809a999f8be920db55ac4a75 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 cae65633687f274352318e6f5f9ac3d5214d4a7a..816d27a650ea46865c1bd42f0655f5751de1ae1d 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 } };