From b043a0c8344858849fccfde177ded09b829a5a79 Mon Sep 17 00:00:00 2001
From: Mark Jones <jones@jlab.org>
Date: Fri, 31 May 2019 12:52:47 -0400
Subject: [PATCH] Modify THcDC and THcHallCSpectrometer

1) THcDC
a) Add variable to tree called dc.InsideDipoleExit
which is a Bool_t flag that is kTRUE if golden track is inside the
dipole exit. Uses THcHallCSpectrometer::InsideDipoleExitWindow
in THcDC::SetFocalPlaneBestTrack to set variable.

2) THcHallCSpectrometer
a) Add method Bool_t InsideDipoleExitWindow which checks whether
the track is within the dipole exit window. Checks
prefix used for database to set flag
fUseSHMSDipoleExitWindow or fUseHMSDipoleExitWindow
to be true which is used in InsideDipoleExitWindow
to select either SHMSDipoleExitWindow or
HMSDipoleExitWindow method to be used as the test.
b) Hardcoded z distance from spectrometer focal plane to
dipole exit window.
c) Methods Bool_t SHMSDipoleExitWindow and
Bool_t HMSDipoleExitWindow return true/false if
the track was inside/outside the dipole exit window.
d) Added dipole exit window test to the BestTrackUsingPrune
method. Need to set parameter hprune_DipoleExit=1 or
pprune_DipoleExit=1 to use the test in the method.
By default it is turned off.
---
 src/THcDC.cxx                |  5 +++
 src/THcDC.h                  |  1 +
 src/THcHallCSpectrometer.cxx | 71 ++++++++++++++++++++++++++++++++++++
 src/THcHallCSpectrometer.h   |  8 +++-
 4 files changed, 84 insertions(+), 1 deletion(-)

diff --git a/src/THcDC.cxx b/src/THcDC.cxx
index 85d027c..37efbb3 100644
--- a/src/THcDC.cxx
+++ b/src/THcDC.cxx
@@ -27,6 +27,7 @@ the number of parameters per plane.
 #include "TMath.h"
 #include "TVectorD.h"
 #include "THaApparatus.h"
+#include "THcHallCSpectrometer.h"
 
 #include <cstring>
 #include <cstdio>
@@ -413,6 +414,7 @@ Int_t THcDC::DefineVariables( EMode mode )
     { "chisq", "chisq/dof (golden track) ", "fChisq_best"},
     { "sp1_id", " (golden track) ", "fSp1_ID_best"},
     { "sp2_id", " (golden track) ", "fSp2_ID_best"},
+    { "InsideDipoleExit", " ","fInSideDipoleExit_best"},
     { "gtrack_nsp", " Number of space points in golden track ", "fNsp_best"},
     { "residual", "Residuals", "fResiduals"},
     { "residualExclPlane", "Residuals", "fResidualsExclPlane"},
@@ -507,6 +509,7 @@ void THcDC::ClearEvent()
   fYp_fp_best=-10000.;
   fChisq_best=kBig;
   fNsp_best=0;
+  fInSideDipoleExit_best = kTRUE;
   for(UInt_t i=0;i<fNChambers;i++) {
     fChambers[i]->Clear();
   }
@@ -648,6 +651,8 @@ void THcDC::SetFocalPlaneBestTrack(Int_t golden_track_index)
       fY_fp_best=tr1->GetY();
       fXp_fp_best=tr1->GetXP();
       fYp_fp_best=tr1->GetYP();
+      THcHallCSpectrometer *app = dynamic_cast<THcHallCSpectrometer*>(GetApparatus());
+      fInSideDipoleExit_best = app->InsideDipoleExitWindow(fX_fp_best, fXp_fp_best ,fY_fp_best,fYp_fp_best);
       fSp1_ID_best=tr1->GetSp1_ID();
       fSp2_ID_best=tr1->GetSp2_ID();
       fChisq_best=tr1->GetChisq();
diff --git a/src/THcDC.h b/src/THcDC.h
index 8fbd3a6..f0100cd 100644
--- a/src/THcDC.h
+++ b/src/THcDC.h
@@ -176,6 +176,7 @@ protected:
   Double_t fChisq_best;
   Int_t fSp1_ID_best;
   Int_t fSp2_ID_best;
+  Bool_t fInSideDipoleExit_best;
  // For accumulating statitics for efficiencies
   Int_t fTotEvents;
   Int_t* fNChamHits;
diff --git a/src/THcHallCSpectrometer.cxx b/src/THcHallCSpectrometer.cxx
index 254ddc7..0324434 100644
--- a/src/THcHallCSpectrometer.cxx
+++ b/src/THcHallCSpectrometer.cxx
@@ -275,10 +275,12 @@ Int_t THcHallCSpectrometer::ReadDatabase( const TDatime& date )
     {"prune_chibeta",         &fPruneChiBeta,          kDouble,         0,  1},
     {"prune_npmt",            &fPruneNPMT,           kDouble,         0,  1},
     {"prune_fptime",          &fPruneFpTime,             kDouble,         0,  1},
+    {"prune_DipoleExit",          &fPruneDipoleExit,             kDouble,         0,  1},
     {0}
   };
 
   // Default values
+  fPruneDipoleExit=0;
   fSelUsingScin = 0;
   fSelUsingPrune = 0;
   fPruneXp = .2;
@@ -295,6 +297,10 @@ Int_t THcHallCSpectrometer::ReadDatabase( const TDatime& date )
   fMispointing_x=999.;
   fMispointing_y=999.;
   gHcParms->LoadParmValues((DBRequest*)&list,prefix);
+  fUseHMSDipoleExitWindow=kFALSE;
+  fUseSHMSDipoleExitWindow=kFALSE;
+  if (prefix[0]=='h') fUseHMSDipoleExitWindow=kTRUE;
+  if (prefix[0]=='p') fUseSHMSDipoleExitWindow=kTRUE;
  
   //  mispointing in transport system y is horizontal and +x is vertical down
   if (fMispointing_y == 999.) {
@@ -862,6 +868,32 @@ Int_t THcHallCSpectrometer::BestTrackUsingPrune()
     PruneSelect++;
     if (nGood==1 && fPruneSelect ==0 && fNtracks>1) fPruneSelect=PruneSelect;
 
+   // !     Prune on dipole exit
+    nGood = 0;
+    for (Int_t ptrack = 0; ptrack < fNtracks; ptrack++ ){
+      Double_t xfp=testTracks[ptrack]->GetX();
+      Double_t yfp=testTracks[ptrack]->GetY();
+      Double_t xpfp=testTracks[ptrack]->GetTheta();
+       Double_t ypfp=testTracks[ptrack]->GetPhi();
+       if ( fPruneDipoleExit==1 && InsideDipoleExitWindow(xfp,xpfp,yfp,ypfp) && ( keep[ptrack] ) ){
+	nGood ++;
+      }
+    }
+    if (nGood > 0 ) {
+      for (Int_t ptrack = 0; ptrack < fNtracks; ptrack++ ){
+      Double_t xfp=testTracks[ptrack]->GetX();
+      Double_t yfp=testTracks[ptrack]->GetY();
+      Double_t xpfp=testTracks[ptrack]->GetTheta();
+       Double_t ypfp=testTracks[ptrack]->GetPhi();
+	if (!InsideDipoleExitWindow(xfp,xpfp,yfp,ypfp)  ){
+	  keep[ptrack] = kFALSE;
+	  reject[ptrack] += 30;
+	}
+      }
+    }
+    PruneSelect++;
+    if (nGood==1 && fPruneSelect ==0 && fNtracks>1) fPruneSelect=PruneSelect;
+
     // !     Prune on beta
     nGood = 0;
     for (Int_t ptrack = 0; ptrack < fNtracks; ptrack++ ){
@@ -1074,6 +1106,45 @@ Bool_t THcHallCSpectrometer::IsMyEvent(Int_t evtype) const
 
   return kFALSE;
 }
+//
+Bool_t THcHallCSpectrometer::InsideDipoleExitWindow(Double_t x_fp, Double_t xp_fp, Double_t y_fp, Double_t yp_fp) {
+  Bool_t inside=kTRUE;
+  Double_t DipoleExitWindowZpos=0.;
+  if (fUseSHMSDipoleExitWindow)  DipoleExitWindowZpos=-307.;
+  if (fUseHMSDipoleExitWindow)  DipoleExitWindowZpos=-147.48;
+  Double_t xdip = x_fp + xp_fp*DipoleExitWindowZpos;
+  Double_t ydip = y_fp + yp_fp*DipoleExitWindowZpos;
+  if (fUseSHMSDipoleExitWindow) inside = SHMSDipoleExitWindow(xdip,ydip);
+  if (fUseHMSDipoleExitWindow) inside = HMSDipoleExitWindow(xdip,ydip);
+  return inside;
+}
+//
+Bool_t THcHallCSpectrometer::SHMSDipoleExitWindow(Double_t xdip,Double_t ydip ) {
+    Bool_t insideSHMS=kTRUE;
+    Double_t crad=23.81; // radius of semicircle
+    Double_t voffset= crad-24.035;
+    Double_t hwid=11.549/2.;
+    if ( TMath::Abs(ydip) < hwid) {
+      if (TMath::Abs(xdip) > (crad+voffset)) insideSHMS=kFALSE;
+    } else {
+      if ( ydip >=hwid) {
+      if ( ((xdip-voffset)*(xdip-voffset)+(ydip-hwid)*(ydip-hwid)) > crad*crad) insideSHMS=kFALSE;
+      }
+      if ( ydip <=-hwid) {
+      if ( ((xdip-voffset)*(xdip-voffset)+(ydip+hwid)*(ydip+hwid)) > crad*crad) insideSHMS=kFALSE;
+      }
+    }
+      return insideSHMS;
+ }
+//
+Bool_t  THcHallCSpectrometer::HMSDipoleExitWindow(Double_t xdip,Double_t ydip) {
+  Bool_t insideHMS=kTRUE;
+    Double_t xpipe_offset = 2.8;
+    Double_t ypipe_offset = 0.0;
+    Double_t pipe_rad=46.507;
+    if ( ((xdip-xpipe_offset)*(xdip-xpipe_offset)+(ydip-ypipe_offset)*(ydip-ypipe_offset)) > pipe_rad*pipe_rad) insideHMS=kFALSE; 
+    return insideHMS ;
+}
 
 //_____________________________________________________________________________
 ClassImp(THcHallCSpectrometer)
diff --git a/src/THcHallCSpectrometer.h b/src/THcHallCSpectrometer.h
index beedc84..ee20876 100644
--- a/src/THcHallCSpectrometer.h
+++ b/src/THcHallCSpectrometer.h
@@ -68,13 +68,19 @@ public:
   virtual Int_t GetNumTypes() { return eventtypes.size(); };
   virtual Bool_t IsPresent() {return fPresent;};
 
+  Bool_t InsideDipoleExitWindow(Double_t x_fp, Double_t xp_fp, Double_t y_fp, Double_t yp_fp);
 
 protected:
   void InitializeReconstruction();
 
-  //  Bool_t*      fKeep;
+  Bool_t SHMSDipoleExitWindow(Double_t x_dip, Double_t y_dip);
+  Bool_t HMSDipoleExitWindow(Double_t x_dip, Double_t y_dip);
+  Bool_t fUseSHMSDipoleExitWindow;
+  Bool_t fUseHMSDipoleExitWindow;
+ //  Bool_t*      fKeep;
   //  Int_t*       fReject;
 
+  Double_t     fPruneDipoleExit;
   Double_t     fPartMass;
   Double_t     fPruneXp;
   Double_t     fPruneYp;
-- 
GitLab