From b1b202d266694c3e029f7ef6088afa0a601675a6 Mon Sep 17 00:00:00 2001
From: "Stephen A. Wood" <zviwood@gmail.com>
Date: Sun, 20 Jan 2013 21:20:46 -0500
Subject: [PATCH] Drift chamber code work. Modeled on podd and h_trans_dc.f
Start adding to the drift chamber code. Setup a class structure
similar to what podd uses for the VDCs. Raw hit class renamed
to THcRawDCHit. Signal hits are how THcDCHit, modeled after
podd hits. Add THcDCWire, and classes for time to distance conversion.
Time to distnace conversion doesn't do anything yet.
Find some improved DC parameter files so that min and max TDC limits
are good for the run we are using.
THcDriftChamberPlane asks ThcDriftChamber for many parameters.
---
Makefile | 4 +-
examples/PARAM/general.param | 1 +
examples/PARAM/hdc.param | 32 +++++------
examples/dchitmaps.C | 2 +-
src/HallC_LinkDef.h | 4 ++
src/THcDCHit.cxx | 104 ++++++++++++++++-------------------
src/THcDCHit.h | 73 +++++++++++++++---------
src/THcDCLookupTTDConv.cxx | 97 ++++++++++++++++++++++++++++++++
src/THcDCLookupTTDConv.h | 48 ++++++++++++++++
src/THcDCTimeToDistConv.cxx | 21 +++++++
src/THcDCTimeToDistConv.h | 32 +++++++++++
src/THcDCWire.cxx | 14 +++++
src/THcDCWire.h | 51 +++++++++++++++++
src/THcDriftChamber.cxx | 53 ++++++++++++++++--
src/THcDriftChamber.h | 49 ++++++++++++++++-
src/THcDriftChamberPlane.cxx | 85 ++++++++++++++++++++++++----
src/THcDriftChamberPlane.h | 34 ++++++++++--
src/THcRawDCHit.cxx | 78 ++++++++++++++++++++++++++
src/THcRawDCHit.h | 37 +++++++++++++
19 files changed, 692 insertions(+), 127 deletions(-)
create mode 100644 src/THcDCLookupTTDConv.cxx
create mode 100644 src/THcDCLookupTTDConv.h
create mode 100644 src/THcDCTimeToDistConv.cxx
create mode 100644 src/THcDCTimeToDistConv.h
create mode 100644 src/THcDCWire.cxx
create mode 100644 src/THcDCWire.h
create mode 100644 src/THcRawDCHit.cxx
create mode 100644 src/THcRawDCHit.h
diff --git a/Makefile b/Makefile
index 226591b..fadfb29 100644
--- a/Makefile
+++ b/Makefile
@@ -16,7 +16,9 @@ SRC = src/THcInterface.cxx src/THcParmList.cxx src/THcAnalyzer.cxx \
src/THcHodoscope.cxx src/THcScintillatorPlane.cxx \
src/THcHodoscopeHit.cxx \
src/THcDriftChamber.cxx src/THcDriftChamberPlane.cxx \
- src/THcDCHit.cxx \
+ src/THcRawDCHit.cxx src/THcDCHit.cxx \
+ src/THcDCWire.cxx \
+ src/THcDCLookupTTDConv.cxx src/THcDCTimeToDistConv.cxx \
src/THcShower.cxx src/THcShowerPlane.cxx \
src/THcShowerHit.cxx \
src/THcAerogel.cxx src/THcAerogelHit.cxx
diff --git a/examples/PARAM/general.param b/examples/PARAM/general.param
index 6cb0731..898d4b3 100644
--- a/examples/PARAM/general.param
+++ b/examples/PARAM/general.param
@@ -24,4 +24,5 @@ raddeg=3.14159265/180
#include "PARAM/scal.pos"
#include "PARAM/hhodo.param"
#include "PARAM/haero.param"
+#include "PARAM/hdc.param"
diff --git a/examples/PARAM/hdc.param b/examples/PARAM/hdc.param
index 926a370..e1bcd4a 100644
--- a/examples/PARAM/hdc.param
+++ b/examples/PARAM/hdc.param
@@ -16,26 +16,26 @@
0.020
0.020
; hms dc tdc minimum tdc value array allowed for a good hit
- hdc_tdc_min_win = 1500, 1500, 1500, 1500, 1500, 1500
- 1500, 1500, 1500, 1500, 1500, 1500
+ hdc_tdc_min_win = 2900, 2900, 2900, 2900, 2900, 2900
+ 2900, 2900, 2900, 2900, 2900, 2900
; hms dc tdc maximum tdc value array allowed for a good hit
- hdc_tdc_max_win = 3200, 3200, 3200, 3200, 3200, 3200
- 3200, 3200, 3200, 3200, 3200, 3200
+ hdc_tdc_max_win = 3400, 3400, 3400, 3400, 3400, 3400
+ 3400, 3400, 3400, 3400, 3400, 3400
; hms drift chamber tdc's time per channel
hdc_tdc_time_per_channel = 0.50
; hms zero time for drift chambers !DECREASING this number moves the hdtime plots to LOWER time.
- hdc_plane_time_zero = (1670+14+18-220-10)
- (1670+12+20-220-10)
- (1670+13+18-220-10)
- (1670+13+20-220-10)
- (1670+12+20-220-10)
- (1670+14+16-220-10)
- (1670+15+16-220-10)
- (1670+11.5+16-220-10)
- (1670+13+12-220-10)
- (1670+13+12-220-10)
- (1670+10.5+16-220-10)
- (1670+13+18-220-10)
+ hdc_plane_time_zero = (1670+14+18-10)
+ (1670+12+20-10)
+ (1670+13+18-10)
+ (1670+13+20-10)
+ (1670+12+20-10)
+ (1670+14+16-10)
+ (1670+15+16-10)
+ (1670+11.5+16-10)
+ (1670+13+12-10)
+ (1670+13+12-10)
+ (1670+10.5+16-10)
+ (1670+13+18-10)
; Dave Abbott's wire velocity correction
hdc_wire_velocity = 12.0
diff --git a/examples/dchitmaps.C b/examples/dchitmaps.C
index e291872..ce968cd 100644
--- a/examples/dchitmaps.C
+++ b/examples/dchitmaps.C
@@ -1,7 +1,7 @@
{
TFile* f = new TFile("hodtest.root");
- TCanvas *c1 = new TCanvas("c1", "Scintillator Hit Maps", 800, 800);
+ TCanvas *c1 = new TCanvas("c1", "Drift Chamber Hit Maps", 800, 800);
c1->Divide(2, 6);
TH1F* h[12];
diff --git a/src/HallC_LinkDef.h b/src/HallC_LinkDef.h
index 95765c1..608c004 100644
--- a/src/HallC_LinkDef.h
+++ b/src/HallC_LinkDef.h
@@ -20,7 +20,11 @@
#pragma link C++ class THcHodoscopeHit+;
#pragma link C++ class THcDriftChamber+;
#pragma link C++ class THcDriftChamberPlane+;
+#pragma link C++ class THcRawDCHit+;
#pragma link C++ class THcDCHit+;
+#pragma link C++ class THcDCWire+;
+#pragma link C++ class THcDCLookupTTDConv+;
+#pragma link C++ class THcDCTimeToDistConv+;
#pragma link C++ class THcShower+;
#pragma link C++ class THcShowerPlane+;
#pragma link C++ class THcShowerHit+;
diff --git a/src/THcDCHit.cxx b/src/THcDCHit.cxx
index 261245a..d8ad011 100644
--- a/src/THcDCHit.cxx
+++ b/src/THcDCHit.cxx
@@ -1,78 +1,66 @@
///////////////////////////////////////////////////////////////////////////////
// //
-// THcDCHit //
+// THcDCHit //
// //
-// Class representing for drift chamber wire (or other device with //
-// a single multihit TDC channel per detector element //
+// Class representing a single hit for the VDC //
// //
///////////////////////////////////////////////////////////////////////////////
#include "THcDCHit.h"
+#include "THcDCTimeToDistConv.h"
-using namespace std;
-
+ClassImp(THcDCHit)
-void THcDCHit::SetData(Int_t signal, Int_t data) {
- fTDC[fNHits++] = data;
-}
+const Double_t THcDCHit::kBig = 1.e38; // Arbitrary large value
-// Return just the first hit
-Int_t THcDCHit::GetData(Int_t signal) {
- if(fNHits>0) {
- return(fTDC[0]);
- } else {
- return(-1);
+//_____________________________________________________________________________
+Double_t THcDCHit::ConvertTimeToDist(Double_t slope)
+{
+ // Converts TDC time to drift distance
+ // Takes the (estimated) slope of the track as an argument
+
+ THcDCTimeToDistConv* ttdConv = (fWire) ? fWire->GetTTDConv() : NULL;
+
+ 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);
+ return fDist;
}
-}
+
+ Error("ConvertTimeToDist()", "No Time to dist algorithm available");
+ return 0.0;
-// Return a requested hit
-Int_t THcDCHit::GetData(Int_t signal, Int_t ihit) {
- if(ihit >=0 && ihit< fNHits) {
- return(fTDC[ihit]);
- } else {
- return(-1);
- }
}
-
-Int_t THcDCHit::Compare(const TObject* obj) const
-{
- // Compare to sort by plane and counter
- // Should we be able to move this into THcRawHit
-
- const THcDCHit* hit = dynamic_cast<const THcDCHit*>(obj);
-
- if(!hit) return -1;
- Int_t p1 = fPlane;
- Int_t p2 = hit->fPlane;
- if(p1 < p2) return -1;
- else if(p1 > p2) return 1;
- else {
- Int_t c1 = fCounter;
- Int_t c2 = hit->fCounter;
- if(c1 < c2) return -1;
- else if (c1 == c2) return 0;
- else return 1;
- }
-}
//_____________________________________________________________________________
-THcDCHit& THcDCHit::operator=( const THcDCHit& rhs )
+Int_t THcDCHit::Compare( const TObject* obj ) const
{
- // Assignment operator.
+ // Used to sort hits
+ // A hit is "less than" another hit if it occurred on a lower wire number.
+ // Also, for hits on the same wire, the first hit on the wire (the one with
+ // the smallest time) is "less than" one with a higher time. If the hits
+ // are sorted according to this scheme, they will be in order of increasing
+ // wire number and, for each wire, will be in the order in which they hit
+ // the wire
- THcRawHit::operator=(rhs);
- if ( this != &rhs ) {
- fPlane = rhs.fPlane;
- fCounter = rhs.fCounter;
- fNHits = rhs.fNHits;
- for(Int_t ihit=0;ihit<fNHits;ihit++) {
- fTDC[ihit] = rhs.fTDC[ihit];
- }
+ if( !obj || IsA() != obj->IsA() || !fWire )
+ return -1;
+
+ const THcDCHit* hit = static_cast<const THcDCHit*>( obj );
+
+ Int_t myWireNum = fWire->GetNum();
+ Int_t hitWireNum = hit->GetWire()->GetNum();
+ // Compare wire numbers
+ if (myWireNum < hitWireNum) return -1;
+ if (myWireNum > hitWireNum) return 1;
+ if (myWireNum == hitWireNum) {
+ // If wire numbers are the same, compare times
+ Double_t hitTime = hit->GetTime();
+ if (fTime < hitTime) return -1;
+ if (fTime > hitTime) return 1;
}
- return *this;
+ return 0;
}
-
-//////////////////////////////////////////////////////////////////////////
-ClassImp(THcDCHit)
-
+////////////////////////////////////////////////////////////////////////////////
diff --git a/src/THcDCHit.h b/src/THcDCHit.h
index 5523c79..3e5150e 100644
--- a/src/THcDCHit.h
+++ b/src/THcDCHit.h
@@ -1,37 +1,60 @@
#ifndef ROOT_THcDCHit
#define ROOT_THcDCHit
-#include "THcRawHit.h"
+///////////////////////////////////////////////////////////////////////////////
+// //
+// THcDCHit //
+// //
+///////////////////////////////////////////////////////////////////////////////
-#define MAXHITS 16
+#include "TObject.h"
+#include "THcDCWire.h"
+#include <cstdio>
-class THcDCHit : public THcRawHit {
+class THcDCHit : public TObject {
- public:
-
- THcDCHit(Int_t plane=0, Int_t counter=0) : THcRawHit(plane, counter),
- fNHits(0) {
- }
- THcDCHit& operator=( const THcDCHit& );
+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) {}
virtual ~THcDCHit() {}
- virtual void Clear( Option_t* opt="" ) { fNHits=0; }
-
- void SetData(Int_t signal, Int_t data);
- Int_t GetData(Int_t signal);
- Int_t GetData(Int_t signal, Int_t ihit);
-
- virtual Bool_t IsSortable () const {return kTRUE; }
- virtual Int_t Compare(const TObject* obj) const;
-
- Int_t fNHits;
- Int_t fTDC[MAXHITS];
-
- protected:
-
+ virtual Double_t ConvertTimeToDist(Double_t slope);
+ Int_t Compare ( const TObject* obj ) const;
+ Bool_t IsSortable () const { return kTRUE; }
+
+ // Get and Set Functions
+ THcDCWire* GetWire() const { return fWire; }
+ Int_t GetWireNum() const { return fWire->GetNum(); }
+ Int_t GetRawTime() const { return fRawTime; }
+ Double_t GetTime() const { return fTime; }
+ Double_t GetDist() const { return fDist; }
+ Double_t GetPos() const { return fWire->GetPos(); } //Position of hit wire
+ Double_t GetdDist() const { return fdDist; }
+
+ void SetWire(THcDCWire * wire) { fWire = wire; }
+ void SetRawTime(Int_t time) { fRawTime = time; }
+ void SetTime(Double_t time) { fTime = time; }
+ void SetDist(Double_t dist) { fDist = dist; }
+ void SetdDist(Double_t ddist) { fdDist = ddist; }
+ void SetFitDist(Double_t dist) { ftrDist = dist; }
+
+protected:
+ static const Double_t kBig; //!
+
+ THcDCWire* fWire; // Wire on which the hit occurred
+ Int_t fRawTime; // TDC value (channels)
+ Double_t fTime; // Time corrected for time offset of wire (s)
+ Double_t fDist; // (Perpendicular) Drift Distance
+ Double_t fdDist; // uncertainty in fDist (for chi2 calc)
+ Double_t ftrDist; // (Perpendicular) distance from the track
+
private:
+ THcDCHit( const THcDCHit& );
+ THcDCHit& operator=( const THcDCHit& );
+
+ ClassDef(THcDCHit,2) // VDCHit class
+};
- ClassDef(THcDCHit, 0); // DC hit class
-};
+////////////////////////////////////////////////////////////////////////////////
#endif
diff --git a/src/THcDCLookupTTDConv.cxx b/src/THcDCLookupTTDConv.cxx
new file mode 100644
index 0000000..0485c1a
--- /dev/null
+++ b/src/THcDCLookupTTDConv.cxx
@@ -0,0 +1,97 @@
+///////////////////////////////////////////////////////////////////////////////
+// //
+// THcDCLookupTTDConv //
+// //
+///////////////////////////////////////////////////////////////////////////////
+
+#include "THcDCLookupTTDConv.h"
+
+ClassImp(THcDCLookupTTDConv)
+
+
+//______________________________________________________________________________
+THcDCLookupTTDConv::THcDCLookupTTDConv()
+{
+ //Normal constructor
+}
+
+//______________________________________________________________________________
+THcDCLookupTTDConv::THcDCLookupTTDConv( Double_t vel)
+{
+ // Normal constructor
+ fDriftVel = vel;
+
+ // 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()
+{
+ // Destructor. Remove variables from global list.
+
+}
+
+//______________________________________________________________________________
+Double_t THcDCLookupTTDConv::ConvertTimeToDist(Double_t time,
+ Double_t tanTheta,
+ Double_t *ddist)
+{
+ // 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;
+ }
+
+ if (ddist) *ddist = unc;
+// printf("D(%e) = %e\nUncorrected D = %e\n", time, dist, fDriftVel * time);
+
+ return dist;
+
+}
+
+
+////////////////////////////////////////////////////////////////////////////////
diff --git a/src/THcDCLookupTTDConv.h b/src/THcDCLookupTTDConv.h
new file mode 100644
index 0000000..3ae38a5
--- /dev/null
+++ b/src/THcDCLookupTTDConv.h
@@ -0,0 +1,48 @@
+#ifndef ROOT_THcDCLookupTTDConv
+#define ROOT_THcDCLookupTTDConv
+
+///////////////////////////////////////////////////////////////////////////////
+// //
+// THcDCLookupTTDConv //
+// //
+// Uses a drift velocity (um/ns) to convert time (ns) into distance (cm) //
+// //
+///////////////////////////////////////////////////////////////////////////////
+#include "THcDCTimeToDistConv.h"
+
+class THcDCLookupTTDConv : public THcDCTimeToDistConv{
+
+public:
+ THcDCLookupTTDConv( );
+ THcDCLookupTTDConv(Double_t vel);
+
+ virtual ~THcDCLookupTTDConv();
+
+ virtual Double_t ConvertTimeToDist(Double_t time, Double_t tanTheta,
+ Double_t *ddist=0);
+
+
+ // 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
+
+ ClassDef(THcDCLookupTTDConv,0) // VDC Analytic TTD Conv class
+};
+
+
+////////////////////////////////////////////////////////////////////////////////
+
+#endif
diff --git a/src/THcDCTimeToDistConv.cxx b/src/THcDCTimeToDistConv.cxx
new file mode 100644
index 0000000..1596589
--- /dev/null
+++ b/src/THcDCTimeToDistConv.cxx
@@ -0,0 +1,21 @@
+///////////////////////////////////////////////////////////////////////////////
+// //
+// THcDCTimeToDistConv //
+// //
+///////////////////////////////////////////////////////////////////////////////
+
+#include "THcDCTimeToDistConv.h"
+
+
+ClassImp(THcDCTimeToDistConv)
+
+
+//______________________________________________________________________________
+THcDCTimeToDistConv::~THcDCTimeToDistConv()
+{
+ // Destructor.
+
+}
+
+
+////////////////////////////////////////////////////////////////////////////////
diff --git a/src/THcDCTimeToDistConv.h b/src/THcDCTimeToDistConv.h
new file mode 100644
index 0000000..2582249
--- /dev/null
+++ b/src/THcDCTimeToDistConv.h
@@ -0,0 +1,32 @@
+#ifndef ROOT_THcDCTimeToDistConv
+#define ROOT_THcDCTimeToDistConv
+
+///////////////////////////////////////////////////////////////////////////////
+// //
+// THcDCTimeToDistConv //
+// //
+// Base class for algorithms for converting TDC time into perpendicular //
+// drift distance //
+///////////////////////////////////////////////////////////////////////////////
+
+#include "TObject.h"
+
+class THcDCTimeToDistConv : public TObject {
+
+public:
+ THcDCTimeToDistConv() {}
+ virtual ~THcDCTimeToDistConv();
+
+ virtual Double_t ConvertTimeToDist(Double_t time, Double_t tanTheta,
+ Double_t *ddist=0) = 0;
+private:
+
+ THcDCTimeToDistConv( const THcDCTimeToDistConv& );
+ THcDCTimeToDistConv& operator=( const THcDCTimeToDistConv& );
+
+ ClassDef(THcDCTimeToDistConv,0) // DCTimeToDistConv class
+};
+
+////////////////////////////////////////////////////////////////////////////////
+
+#endif
diff --git a/src/THcDCWire.cxx b/src/THcDCWire.cxx
new file mode 100644
index 0000000..3b6520b
--- /dev/null
+++ b/src/THcDCWire.cxx
@@ -0,0 +1,14 @@
+///////////////////////////////////////////////////////////////////////////////
+// //
+// THcDCWire //
+// //
+// Class to represent a drift chamber wire //
+// //
+///////////////////////////////////////////////////////////////////////////////
+
+#include "THcDCWire.h"
+
+ClassImp(THcDCWire)
+
+
+///////////////////////////////////////////////////////////////////////////////
diff --git a/src/THcDCWire.h b/src/THcDCWire.h
new file mode 100644
index 0000000..7f2ab5b
--- /dev/null
+++ b/src/THcDCWire.h
@@ -0,0 +1,51 @@
+#ifndef ROOT_THcDCWire
+#define ROOT_THcDCWire
+
+///////////////////////////////////////////////////////////////////////////////
+// //
+// THcDCWire //
+// //
+///////////////////////////////////////////////////////////////////////////////
+#include "TObject.h"
+
+class THcDCTimeToDistConv;
+
+class THcDCWire : public TObject {
+
+public:
+
+ THcDCWire( Int_t num=0, Double_t pos=0.0, Double_t offset=0.0,
+ THcDCTimeToDistConv* ttd=NULL ) :
+ fNum(num), fFlag(0), fPos(pos), fTOffset(offset), fTTDConv(ttd) {}
+ virtual ~THcDCWire() {}
+
+ // Get and Set Functions
+ Int_t GetNum() const { return fNum; }
+ Int_t GetFlag() const { return fFlag; }
+ Double_t GetPos() const { return fPos; }
+ Double_t GetTOffset() const { return fTOffset; }
+ THcDCTimeToDistConv * GetTTDConv() { return fTTDConv; }
+
+ void SetNum (Int_t num) {fNum = num;}
+ void SetFlag (Int_t flag) {fFlag = flag;}
+ void SetPos (Double_t pos) { fPos = pos; }
+ void SetTOffset (Double_t tOffset){ fTOffset = tOffset; }
+ void SetTTDConv (THcDCTimeToDistConv * ttdConv){ fTTDConv = ttdConv;}
+
+protected:
+ Int_t fNum; //Wire Number
+ Int_t fFlag; //Flag for errors (e.g. Bad wire)
+ Double_t fPos; //Position within the plane
+ Double_t fTOffset; //Timing Offset
+ THcDCTimeToDistConv* fTTDConv; //!Time to Distance Converter
+
+private:
+ THcDCWire( const THcDCWire& );
+ THcDCWire& operator=( const THcDCWire& );
+
+ ClassDef(THcDCWire,0) // VDCWire class
+};
+
+////////////////////////////////////////////////////////////////////////////////
+
+#endif
diff --git a/src/THcDriftChamber.cxx b/src/THcDriftChamber.cxx
index 0b2ad19..1c73a1e 100644
--- a/src/THcDriftChamber.cxx
+++ b/src/THcDriftChamber.cxx
@@ -62,12 +62,14 @@ void THcDriftChamber::Setup(const char* name, const char* description)
DBRequest list[]={
{"dc_num_planes",&fNPlanes, kInt},
- {"dc_num_chambers",&fNChambers, kDouble},
+ {"dc_num_chambers",&fNChambers, kInt},
+ {"dc_tdc_time_per_channel",&fNSperChan, kDouble},
+ {"dc_wire_velocity",&fWireVelocity,kDouble},
{0}
};
gHcParms->LoadParmValues((DBRequest*)&list,prefix);
-
+
cout << "Drift Chambers: " << fNPlanes << " planes in " << fNChambers << " chambers" << endl;
fPlaneNames = new char* [fNPlanes];
@@ -105,7 +107,7 @@ THaAnalysisObject::EStatus THcDriftChamber::Init( const TDatime& date )
// Should probably put this in ReadDatabase as we will know the
// maximum number of hits after setting up the detector map
- THcHitList::InitHitList(fDetMap, "THcDCHit", 1000);
+ THcHitList::InitHitList(fDetMap, "THcRawDCHit", 1000);
EStatus status;
// This triggers call of ReadDatabase and DefineVariables
@@ -169,9 +171,48 @@ Int_t THcDriftChamber::ReadDatabase( const TDatime& date )
prefix[1]='\0';
+ fXCenter = new Double_t [fNChambers];
+ fYCenter = new Double_t [fNChambers];
+
+ fTdcWinMin = new Int_t [fNPlanes];
+ fTdcWinMax = new Int_t [fNPlanes];
+ fCentralTime = new Int_t [fNPlanes];
fNWires = new Int_t [fNPlanes];
+ fNChamber = new Int_t [fNPlanes]; // Which chamber is this plane
+ fWireOrder = new Int_t [fNPlanes]; // Wire readout order
+ fDriftTimeSign = new Int_t [fNPlanes];
+
+ fZPos = new Double_t [fNPlanes];
+ fAlphaAngle = new Double_t [fNPlanes];
+ fBetaAngle = new Double_t [fNPlanes];
+ fGammaAngle = new Double_t [fNPlanes];
+ fPitch = new Double_t [fNPlanes];
+ fCentralWire = new Double_t [fNPlanes];
+ fPlaneTimeZero = new Double_t [fNPlanes];
+
+
DBRequest list[]={
- {"dc_nrwire",fNWires, kInt, fNPlanes},
+ {"dc_tdc_time_per_channel",&fNSperChan, kDouble},
+ {"dc_wire_velocity",&fWireVelocity,kDouble},
+
+ {"dc_xcenter", fXCenter, kDouble, fNChambers},
+ {"dc_ycenter", fYCenter, kDouble, fNChambers},
+
+ {"dc_tdc_min_win", fTdcWinMin, kInt, fNPlanes},
+ {"dc_tdc_max_win", fTdcWinMax, kInt, fNPlanes},
+ {"dc_central_time", fCentralTime, kInt, fNPlanes},
+ {"dc_nrwire", fNWires, kInt, fNPlanes},
+ {"dc_chamber_planes", fNChamber, kInt, fNPlanes},
+ {"dc_wire_counting", fWireOrder, kInt, fNPlanes},
+ {"dc_drifttime_sign", fDriftTimeSign, kInt, fNPlanes},
+
+ {"dc_zpos", fZPos, kDouble, fNPlanes},
+ {"dc_alpha_angle", fAlphaAngle, kDouble, fNPlanes},
+ {"dc_beta_angle", fBetaAngle, kDouble, fNPlanes},
+ {"dc_gamma_angle", fGammaAngle, kDouble, fNPlanes},
+ {"dc_pitch", fPitch, kDouble, fNPlanes},
+ {"dc_central_wire", fCentralWire, kDouble, fNPlanes},
+ {"dc_plane_time_zero", fPlaneTimeZero, kDouble, fNPlanes},
{0}
};
gHcParms->LoadParmValues((DBRequest*)&list,prefix);
@@ -306,9 +347,9 @@ Int_t THcDriftChamber::Decode( const THaEvData& evdata )
}
#if 0
- // fRawHitList is TClones array of THcDCHit objects
+ // fRawHitList is TClones array of THcRawDCHit objects
for(Int_t ihit = 0; ihit < fNRawHits ; ihit++) {
- THcDCHit* hit = (THcDCHit *) fRawHitList->At(ihit);
+ THcRawDCHit* hit = (THcRawDCHit *) fRawHitList->At(ihit);
// cout << ihit << " : " << hit->fPlane << ":" << hit->fCounter << " : "
// << endl;
for(Int_t imhit = 0; imhit < hit->fNHits; imhit++) {
diff --git a/src/THcDriftChamber.h b/src/THcDriftChamber.h
index f96363d..a973c46 100644
--- a/src/THcDriftChamber.h
+++ b/src/THcDriftChamber.h
@@ -9,8 +9,9 @@
#include "THaTrackingDetector.h"
#include "THcHitList.h"
-#include "THcDCHit.h"
+#include "THcRawDCHit.h"
#include "THcDriftChamberPlane.h"
+#include "TMath.h"
//class THaScCalib;
class TClonesArray;
@@ -33,7 +34,25 @@ public:
Int_t GetNTracks() const { return fTrackProj->GetLast()+1; }
const TClonesArray* GetTrackHits() const { return fTrackProj; }
-
+
+ Int_t GetNWires(Int_t plane) const { return fNWires[plane-1];}
+ Int_t GetNChamber(Int_t plane) const { return fNChamber[plane-1];}
+ Int_t GetWireOrder(Int_t plane) const { return fWireOrder[plane-1];}
+ Int_t GetPitch(Int_t plane) const { return fPitch[plane-1];}
+ Int_t GetCentralWire(Int_t plane) const { return fCentralWire[plane-1];}
+ Int_t GetTdcWinMin(Int_t plane) const { return fTdcWinMin[plane-1];}
+ Int_t GetTdcWinMax(Int_t plane) const { return fTdcWinMax[plane-1];}
+
+ Double_t GetPlaneTimeZero(Int_t plane) const { return fPlaneTimeZero[plane-1];}
+
+ Double_t GetNSperChan() const { return fNSperChan;}
+
+ Double_t GetCenter(Int_t plane) const {
+ Int_t chamber = GetNChamber(plane)-1;
+ return
+ fXCenter[chamber]*sin(fAlphaAngle[plane-1]) +
+ fYCenter[chamber]*cos(fAlphaAngle[plane-1]);
+ }
// friend class THaScCalib;
THcDriftChamber(); // for ROOT I/O
@@ -43,12 +62,36 @@ protected:
// Per-event data
-
// Potential Hall C parameters. Mostly here for demonstration
Int_t fNPlanes;
char** fPlaneNames;
Int_t fNChambers;
+
+ Double_t fNSperChan; /* TDC bin size */
+ Double_t fWireVelocity;
+
+ // Each of these will be dimensioned with the number of chambers
+ Double_t* fXCenter;
+ Double_t* fYCenter;
+
+ // Each of these will be dimensioned with the number of planes
+ // A THcDriftChamberPlane class object will need to access the value for
+ // its plane number. Should we have a Get method for each or
+ Int_t* fTdcWinMin;
+ Int_t* fTdcWinMax;
+ Int_t* fCentralTime;
Int_t* fNWires; // Number of wires per plane
+ Int_t* fNChamber;
+ Int_t* fWireOrder;
+ Int_t* fDriftTimeSign;
+
+ Double_t* fZPos;
+ Double_t* fAlphaAngle;
+ Double_t* fBetaAngle;
+ Double_t* fGammaAngle;
+ Double_t* fPitch;
+ Double_t* fCentralWire;
+ Double_t* fPlaneTimeZero;
THcDriftChamberPlane** fPlanes; // List of plane objects
diff --git a/src/THcDriftChamberPlane.cxx b/src/THcDriftChamberPlane.cxx
index a4316ce..5408b4f 100644
--- a/src/THcDriftChamberPlane.cxx
+++ b/src/THcDriftChamberPlane.cxx
@@ -7,7 +7,9 @@
//////////////////////////////////////////////////////////////////////////
#include "THcDriftChamberPlane.h"
-#include "TClonesArray.h"
+#include "THcDCWire.h"
+#include "THcDCHit.h"
+#include "THcDCLookupTTDConv.h"
#include "THcSignalHit.h"
#include "THcGlobals.h"
#include "THcParmList.h"
@@ -32,7 +34,9 @@ THcDriftChamberPlane::THcDriftChamberPlane( const char* name,
: THaSubDetector(name,description,parent)
{
// Normal constructor with name and description
- fTDCHits = new TClonesArray("THcSignalHit",100);
+ fHits = new TClonesArray("THcDCHit",100);
+ fWires = new TClonesArray("THcDCWire", 100);
+
fPlaneNum = planenum;
}
@@ -40,7 +44,9 @@ THcDriftChamberPlane::THcDriftChamberPlane( const char* name,
THcDriftChamberPlane::~THcDriftChamberPlane()
{
// Destructor
- delete fTDCHits;
+ delete fWires;
+ delete fHits;
+ delete fTTDConv;
}
THaAnalysisObject::EStatus THcDriftChamberPlane::Init( const TDatime& date )
@@ -78,6 +84,36 @@ Int_t THcDriftChamberPlane::ReadDatabase( const TDatime& date )
prefix[1]='\0';
// Retrieve parameters we need
+ THcDriftChamber* fParent;
+
+ fParent = (THcDriftChamber*) GetParent();
+ // These are single variables here, but arrays in THcDriftChamber.
+ fNChamber = fParent->GetNChamber(fPlaneNum);
+ fNWires = fParent->GetNWires(fPlaneNum);
+ fWireOrder = fParent->GetWireOrder(fPlaneNum);
+ fPitch = fParent->GetPitch(fPlaneNum);
+ fCentralWire = fParent->GetCentralWire(fPlaneNum);
+ fTdcWinMin = fParent->GetTdcWinMin(fPlaneNum);
+ fTdcWinMax = fParent->GetTdcWinMax(fPlaneNum);
+ fPlaneTimeZero = fParent->GetPlaneTimeZero(fPlaneNum);
+ fCenter = fParent->GetCenter(fPlaneNum);
+
+ fNSperChan = fParent->GetNSperChan();
+
+ cout << fPlaneNum << " " << fNWires << endl;
+
+ fTTDConv = new THcDCLookupTTDConv();// Need to pass the lookup table
+
+ Int_t nWires = fParent->GetNWires(fPlaneNum);
+ // For HMS, wire numbers start with one, but arrays start with zero.
+ // So wire number is index+1
+ for (int i=0; i<nWires; i++) {
+ Double_t pos = fPitch*( (fWireOrder==0?(i+1):fNWires-i)
+ - fCentralWire) - fCenter;
+ THcDCWire* wire = new((*fWires)[i])
+ THcDCWire( i+1, pos , 0.0, fTTDConv);
+ //if( something < 0 ) wire->SetFlag(1);
+ }
return kOK;
}
@@ -94,7 +130,9 @@ Int_t THcDriftChamberPlane::DefineVariables( EMode mode )
// Register variables in global list
RVarDef vars[] = {
{"tdchits", "List of TDC hits",
- "fTDCHits.THcSignalHit.GetPaddleNumber()"},
+ "fHits.THcDCHit.GetWireNum()"},
+ {"rawtdc", "Raw TDC Values",
+ "fHits.THcDCHit.GetRawTime()"},
{ 0 }
};
@@ -106,7 +144,7 @@ void THcDriftChamberPlane::Clear( Option_t* )
{
//cout << " Calling THcDriftChamberPlane::Clear " << GetName() << endl;
// Clears the hit lists
- fTDCHits->Clear();
+ fHits->Clear();
}
//_____________________________________________________________________________
@@ -137,22 +175,45 @@ Int_t THcDriftChamberPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
// Assumes that the hit list is sorted by plane, so we stop when the
// plane doesn't agree and return the index for the next hit.
- Int_t nTDCHits=0;
- fTDCHits->Clear();
+ //Int_t nTDCHits=0;
+ fHits->Clear();
Int_t nrawhits = rawhits->GetLast()+1;
// cout << "THcDriftChamberPlane::ProcessHits " << fPlaneNum << " " << nexthit << "/" << nrawhits << endl;
Int_t ihit = nexthit;
+ Int_t nextHit = 0;
while(ihit < nrawhits) {
- THcDCHit* hit = (THcDCHit *) rawhits->At(ihit);
+ THcRawDCHit* hit = (THcRawDCHit *) rawhits->At(ihit);
if(hit->fPlane > fPlaneNum) {
break;
}
- // Just put in the first hit for now
- if(hit->fNHits > 0) { // Should always be the case
- THcSignalHit *sighit = (THcSignalHit*) fTDCHits->ConstructedAt(nTDCHits++);
- sighit->Set(hit->fCounter, hit->fTDC[0]);
+ Int_t wireNum = hit->fCounter;
+ THcDCWire* wire = GetWire(wireNum);
+ Int_t wire_last = -1;
+ for(Int_t mhit=0; mhit<hit->fNHits; mhit++) {
+ /* Sort into early, late and ontime */
+ Int_t rawtdc = hit->fTDC[mhit];
+ if(rawtdc < fTdcWinMin) {
+ // Increment early counter (Actually late because TDC is backward)
+ } else if (rawtdc > fTdcWinMax) {
+ // Increment late count
+ } else {
+ // A good hit
+ if(wire_last == wireNum) {
+ // Increment extra hit counter
+ // Are we choosing the correct hit in the case of multiple hits?
+ // Are we choose the same hit that ENGINE chooses?
+ // cout << "Extra hit " << fPlaneNum << " " << wireNum << " " << rawtdc << endl;
+ } else {
+ Double_t time = // -hstart_time (comes from h_trans_scin
+ - rawtdc*fNSperChan + fPlaneTimeZero;
+ // How do we get this start time from the hodoscope to here
+ // (or at least have it ready by coarse process)
+ new( (*fHits)[nextHit++] ) THcDCHit(wire, rawtdc, time);
+ }
+ wire_last = wireNum;
+ }
}
ihit++;
}
diff --git a/src/THcDriftChamberPlane.h b/src/THcDriftChamberPlane.h
index 29f4bb7..152796d 100644
--- a/src/THcDriftChamberPlane.h
+++ b/src/THcDriftChamberPlane.h
@@ -14,9 +14,14 @@
#include "THaSubDetector.h"
#include "TClonesArray.h"
+#include <cassert>
class THaEvData;
-class THaSignalHit;
+class THcDCWire;
+class THcDCHit;
+class THcDCTimeToDistConv;
+
+/*class THaSignalHit;*/
class THcDriftChamberPlane : public THaSubDetector {
@@ -36,19 +41,38 @@ class THcDriftChamberPlane : public THaSubDetector {
virtual Int_t ProcessHits(TClonesArray* rawhits, Int_t nexthit);
- Double_t fSpacing;
-
- TClonesArray* fParentHitList;
+ // Get and Set functions
+ Int_t GetNWires() const { return fWires->GetLast()+1; }
+ THcDCWire* GetWire(Int_t i) const
+ { assert( i>=1 && i<=GetNWires() );
+ return (THcDCWire*)fWires->UncheckedAt(i-1); }
protected:
- TClonesArray* fTDCHits;
+ TClonesArray* fParentHitList;
+
+ TClonesArray* fHits;
+ TClonesArray* fWires;
Int_t fPlaneNum;
+ Int_t fNChamber;
+ Int_t fNWires;
+ Int_t fWireOrder;
+ Int_t fTdcWinMin;
+ Int_t fTdcWinMax;
+ Double_t fPitch;
+ Double_t fCentralWire;
+ Double_t fPlaneTimeZero;
+
+ Double_t fCenter;
+
+ Double_t fNSperChan; /* TDC bin size */
virtual Int_t ReadDatabase( const TDatime& date );
virtual Int_t DefineVariables( EMode mode = kDefine );
+ THcDCTimeToDistConv* fTTDConv; // Time-to-distance converter for this plane's wires
+
ClassDef(THcDriftChamberPlane,0)
};
#endif
diff --git a/src/THcRawDCHit.cxx b/src/THcRawDCHit.cxx
new file mode 100644
index 0000000..ac317e1
--- /dev/null
+++ b/src/THcRawDCHit.cxx
@@ -0,0 +1,78 @@
+///////////////////////////////////////////////////////////////////////////////
+// //
+// THcRawDCHit //
+// //
+// Class representing for drift chamber wire (or other device with //
+// a single multihit TDC channel per detector element //
+// //
+///////////////////////////////////////////////////////////////////////////////
+
+#include "THcRawDCHit.h"
+
+using namespace std;
+
+
+void THcRawDCHit::SetData(Int_t signal, Int_t data) {
+ fTDC[fNHits++] = data;
+}
+
+// Return just the first hit
+Int_t THcRawDCHit::GetData(Int_t signal) {
+ if(fNHits>0) {
+ return(fTDC[0]);
+ } else {
+ return(-1);
+ }
+}
+
+// Return a requested hit
+Int_t THcRawDCHit::GetData(Int_t signal, Int_t ihit) {
+ if(ihit >=0 && ihit< fNHits) {
+ return(fTDC[ihit]);
+ } else {
+ return(-1);
+ }
+}
+
+
+Int_t THcRawDCHit::Compare(const TObject* obj) const
+{
+ // Compare to sort by plane and counter
+ // Should we be able to move this into THcRawHit
+
+ const THcRawDCHit* hit = dynamic_cast<const THcRawDCHit*>(obj);
+
+ if(!hit) return -1;
+ Int_t p1 = fPlane;
+ Int_t p2 = hit->fPlane;
+ if(p1 < p2) return -1;
+ else if(p1 > p2) return 1;
+ else {
+ Int_t c1 = fCounter;
+ Int_t c2 = hit->fCounter;
+ if(c1 < c2) return -1;
+ else if (c1 == c2) return 0;
+ else return 1;
+ }
+}
+//_____________________________________________________________________________
+THcRawDCHit& THcRawDCHit::operator=( const THcRawDCHit& rhs )
+{
+ // Assignment operator.
+
+ THcRawHit::operator=(rhs);
+ if ( this != &rhs ) {
+ fPlane = rhs.fPlane;
+ fCounter = rhs.fCounter;
+ fNHits = rhs.fNHits;
+ for(Int_t ihit=0;ihit<fNHits;ihit++) {
+ fTDC[ihit] = rhs.fTDC[ihit];
+ }
+ }
+ return *this;
+}
+
+
+//////////////////////////////////////////////////////////////////////////
+ClassImp(THcRawDCHit)
+
diff --git a/src/THcRawDCHit.h b/src/THcRawDCHit.h
new file mode 100644
index 0000000..7aca52b
--- /dev/null
+++ b/src/THcRawDCHit.h
@@ -0,0 +1,37 @@
+#ifndef ROOT_THcRawDCHit
+#define ROOT_THcRawDCHit
+
+#include "THcRawHit.h"
+
+#define MAXHITS 16
+
+class THcRawDCHit : public THcRawHit {
+
+ public:
+
+ THcRawDCHit(Int_t plane=0, Int_t counter=0) : THcRawHit(plane, counter),
+ fNHits(0) {
+ }
+ THcRawDCHit& operator=( const THcRawDCHit& );
+ virtual ~THcRawDCHit() {}
+
+ virtual void Clear( Option_t* opt="" ) { fNHits=0; }
+
+ void SetData(Int_t signal, Int_t data);
+ Int_t GetData(Int_t signal);
+ Int_t GetData(Int_t signal, Int_t ihit);
+
+ virtual Bool_t IsSortable () const {return kTRUE; }
+ virtual Int_t Compare(const TObject* obj) const;
+
+ Int_t fNHits;
+ Int_t fTDC[MAXHITS];
+
+ protected:
+
+ private:
+
+ ClassDef(THcRawDCHit, 0); // DC hit class
+};
+
+#endif
--
GitLab