Hall-D Software  alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
JEventProcessor_FDC_InternalAlignment.cc
Go to the documentation of this file.
1 // $Id$
2 //
3 // File: JEventProcessor_FDC_InternalAlignment.cc
4 // Created: Sun Nov 27 16:10:26 EST 2016
5 // Creator: mstaib (on Linux ifarm1102 2.6.32-431.el6.x86_64 x86_64)
6 //
7 
9 using namespace jana;
10 #include "HistogramTools.h"
11 #include "FDC/DFDCPseudo.h"
12 
13 // Routine used to create our JEventProcessor
14 #include <JANA/JApplication.h>
15 #include <JANA/JFactory.h>
16 extern "C"{
17 void InitPlugin(JApplication *app){
18  InitJANAPlugin(app);
19  app->AddProcessor(new JEventProcessor_FDC_InternalAlignment());
20 }
21 } // "C"
22 
23 
24 //------------------
25 // JEventProcessor_FDC_InternalAlignment (Constructor)
26 //------------------
28 {
29 
30 }
31 
32 //------------------
33 // ~JEventProcessor_FDC_InternalAlignment (Destructor)
34 //------------------
36 {
37 
38 }
39 
40 //------------------
41 // init
42 //------------------
44 {
45  // This is called once at program startup.
46 
47  gDirectory->mkdir("FDC_InternalAlignment");
48  gDirectory->cd("FDC_InternalAlignment");
49  for (int i=1; i<=24;i++){
50 
51  char thisName[256];
52  sprintf(thisName, "Plane %.2i Wire Position Vs XY", i);
53 
54  Hist3D[i-1]=new TH3I(thisName, thisName, 100, -50., 50., 100, -50., 50., 100, -0.5, 0.5);
55  }
56  // We need the constants used for this iteration
57  // Use a TProfile to avoid problems adding together multiple root files...
58  HistCurrentConstants = new TProfile("CathodeAlignmentConstants", "Constants Used for Cathode Alignment (In CCDB Order)", 450,0.5,450.5);
59 
60  gDirectory->cd("..");
61 
62  return NOERROR;
63 }
64 
65 //------------------
66 // brun
67 //------------------
68 jerror_t JEventProcessor_FDC_InternalAlignment::brun(JEventLoop *eventLoop, int32_t runnumber)
69 {
70  // This is called whenever the run number changes
71  // Get curent cathode alignment constants
72 
73  JCalibration * jcalib = eventLoop->GetJCalibration();
74  vector<map<string,double> >vals;
75  if (jcalib->Get("FDC/cathode_alignment",vals)==false){
76  for(unsigned int i=0; i<vals.size(); i++){
77  map<string,double> &row = vals[i];
78  // Get the offsets from the calibration database
79  HistCurrentConstants->Fill(i*4+1,row["dPhiU"]);
80  HistCurrentConstants->Fill(i*4+2,row["dU"]);
81  HistCurrentConstants->Fill(i*4+3,row["dPhiV"]);
82  HistCurrentConstants->Fill(i*4+4,row["dV"]);
83  }
84  }
85 
86  if (jcalib->Get("FDC/strip_pitches_v2",vals)==false){
87  for(unsigned int i=0; i<vals.size(); i++){
88  map<string,double> &row = vals[i];
89  // Get the offsets from the calibration database
90  HistCurrentConstants->Fill(i*10+101,row["U_SP_1"]);
91  HistCurrentConstants->Fill(i*10+102,row["U_G_1"]);
92  HistCurrentConstants->Fill(i*10+103,row["U_SP_2"]);
93  HistCurrentConstants->Fill(i*10+104,row["U_G_2"]);
94  HistCurrentConstants->Fill(i*10+105,row["U_SP_3"]);
95  HistCurrentConstants->Fill(i*10+106,row["V_SP_1"]);
96  HistCurrentConstants->Fill(i*10+107,row["V_G_1"]);
97  HistCurrentConstants->Fill(i*10+108,row["V_SP_2"]);
98  HistCurrentConstants->Fill(i*10+109,row["V_G_2"]);
99  HistCurrentConstants->Fill(i*10+110,row["V_SP_3"]);
100  }
101  }
102 
103  if (jcalib->Get("FDC/cell_offsets",vals)==false){
104  for(unsigned int i=0; i<vals.size(); i++){
105  map<string,double> &row = vals[i];
106  // Get the offsets from the calibration database
107  HistCurrentConstants->Fill(i*2+401,row["xshift"]);
108  HistCurrentConstants->Fill(i*2+402,row["yshift"]);
109  }
110  }
111 
112  return NOERROR;
113 }
114 
115 //------------------
116 // evnt
117 //------------------
118 jerror_t JEventProcessor_FDC_InternalAlignment::evnt(JEventLoop *loop, uint64_t eventnumber)
119 {
120  vector <const DFDCPseudo*> fdcPseudoVector;
121  loop->Get(fdcPseudoVector);
122 
123  for (unsigned int i = 0; i<fdcPseudoVector.size(); i++){
124  const DFDCPseudo* thisPseudo = fdcPseudoVector[i];
125  if (thisPseudo->status != 6) continue;
126 
127  japp->RootWriteLock();
128  Hist3D[thisPseudo->wire->layer - 1]->Fill(thisPseudo->w, thisPseudo->s, thisPseudo->w_c - thisPseudo->w);
129  japp->RootUnLock();
130 
131  // Plot the wire times
132  char thisName[256];
133  sprintf(thisName, "Plane %.2i Wire t Vs Wire Number", thisPseudo->wire->layer);
134  Fill2DHistogram("FDC_InternalAlignment","Wire t0",thisName,
135  thisPseudo->wire->wire, thisPseudo->time,
136  thisName,
137  96, 0.5, 96.5, 250, -50.0, 200.0);
138 
139  sprintf(thisName, "Plane %.2i Wire Position Vs XY", thisPseudo->wire->layer);
140  Fill2DProfile("FDC_InternalAlignment", "Profile2D", thisName,
141  thisPseudo->w, thisPseudo->s, thisPseudo->w_c - thisPseudo->w,
142  thisName,
143  100, -50.,50., 100,-50., 50.);
144 
145  // Calculate the projection of u,w onto v
146  double sinphiu = sin(thisPseudo->phi_u), sinphiv = sin(thisPseudo->phi_v);
147  double sinphiumphiv = sin(thisPseudo->phi_u-thisPseudo->phi_v);
148  double deltaX = HistCurrentConstants->GetBinContent((thisPseudo->wire->layer-1)*2+401);
149  double deltaU = HistCurrentConstants->GetBinContent((thisPseudo->wire->layer-1)*4+2);
150  double deltaV = HistCurrentConstants->GetBinContent((thisPseudo->wire->layer-1)*4+4);
151  double upred = ((thisPseudo->w-deltaX)*sinphiumphiv + thisPseudo->v*sinphiu)/sinphiv;
152  double vpred = -((thisPseudo->w-deltaX)*sinphiumphiv-thisPseudo->u*sinphiv)/sinphiu;
153 
154  sprintf(thisName, "Plane %.2i u_res vs u", thisPseudo->wire->layer);
155  Fill2DHistogram("FDC_InternalAlignment","Cathode Projections", thisName,
156  thisPseudo->u - deltaU, thisPseudo->u - upred,
157  "u_res Vs. u;u_{local};u-upred",
158  192,-47.5,47.5,200, -0.2,0.2);
159 
160  sprintf(thisName, "Plane %.2i v_res vs v", thisPseudo->wire->layer);
161  Fill2DHistogram("FDC_InternalAlignment","Cathode Projections", thisName,
162  thisPseudo->v - deltaV, thisPseudo->v - vpred,
163  "v_res Vs. v;v_{local};v-vpred",
164  192,-47.5,47.5,200, -0.2,0.2);
165 
166  // For the gains we need to break up the two halves of the cathode planes
167  if((thisPseudo->w-deltaX)<0.){
168  sprintf(thisName, "Plane %.2i u_res vs u", thisPseudo->wire->layer);
169  Fill2DHistogram("FDC_InternalAlignment","Cathode Projections Negative", thisName,
170  thisPseudo->u - deltaU, thisPseudo->u - upred,
171  "u_res Vs. u;u_{local};u-upred",
172  192,-47.5,47.5,200, -0.2,0.2);
173 
174  sprintf(thisName, "Plane %.2i v_res vs v", thisPseudo->wire->layer);
175  Fill2DHistogram("FDC_InternalAlignment","Cathode Projections Negative", thisName,
176  thisPseudo->v - deltaV, thisPseudo->v - vpred,
177  "v_res Vs. v;v_{local};v-vpred",
178  192,-47.5,47.5,200, -0.2,0.2);
179  }
180  else{
181  sprintf(thisName, "Plane %.2i u_res vs u", thisPseudo->wire->layer);
182  Fill2DHistogram("FDC_InternalAlignment","Cathode Projections Positive", thisName,
183  thisPseudo->u - deltaU, thisPseudo->u - upred,
184  "u_res Vs. u;u_{local};u-upred",
185  192,-47.5,47.5,200, -0.2,0.2);
186 
187  sprintf(thisName, "Plane %.2i v_res vs v", thisPseudo->wire->layer);
188  Fill2DHistogram("FDC_InternalAlignment","Cathode Projections Positive", thisName,
189  thisPseudo->v - deltaV, thisPseudo->v - vpred,
190  "v_res Vs. v;v_{local};v-vpred",
191  192,-47.5,47.5,200, -0.2,0.2);
192  }
193 
194  //jout << "==upred " << upred << " u " << thisPseudo->u << endl;
195  //jout << "vpred " << vpred << " v " << thisPseudo->v << endl;
196 
197  }
198 
199 
200  // Test For Lubomir...
201  /*
202  double constants[192] = { -0.0043154, -0.1691660, 0.0007645 , -0.0002993 , 5.0059700 , 5.0065600 , -0.2965900 , 0.5452100 ,
203  0.3594690, 0.2248910 , 0.0004884 , 0.0001668 , 5.0040600 , 5.0063800 , 0.3391100 , 0.5998900 ,
204  -0.2362050, -0.1891460, 0.0005648 , -0.0005000 , 5.0055400 , 5.0069200 , 0.5505000 , 0.1661200 ,
205  -0.1562530, -0.1499240, 0.0000257 , 0.0000318 , 5.0071500 , 5.0048300 , 0.1676100 , -0.2970800 ,
206  -0.2040510, -0.5034870, -0.0002976, -0.0009806 , 5.0056300 , 5.0060900 , -0.4721900 , -0.252970 ,
207  -0.3032560, -0.2500840, 0.0002768 , -0.0003581 , 5.0128400 , 5.0098100 , -0.6313500 , 0.4893600 ,
208  0.1563260, -0.0514443, 0.0001716 , 0.0002051 , 5.0078500 , 5.0084900 , -0.3773300 , 0.2163800 ,
209  0.7422930, 0.3868080 , -0.0004678 , -0.0006898 , 5.0062200 , 5.0058600 , 0.2586500 , 0.0381900 ,
210  -0.1862160, -0.3921980, -0.0003152, 0.0005399 , 5.0058100 , 5.0070700 , 0.4594900 , 0.3248800 ,
211  -0.8467560, 0.0483033 , -0.0003287 , 0.0005756 , 5.0100700 , 5.0068000 , 0.0434000 , -0.2726000 ,
212  -0.1340560, -0.0933582, 0.0012310 , 0.0000364 , 5.0065800 , 5.0076400 , -0.3061300 , -0.277280 ,
213  0.0450700, -0.0604768, -0.0010056, 0.0001712 , 5.0093500 , 5.0060800 , -0.2725000 , 0.1725300 ,
214  -0.2025840, -0.1141940, 0.0007098 , -0.0001474 , 5.0055900 , 5.0052500 , -0.0582500 , 0.3434700 ,
215  -0.0267204, -0.0445084, 0.0002086 , 0.0011202 , 5.0062000 , 5.0040900 , 0.2537100 , 0.2302600 ,
216  0.1929130, -0.2376780, 0.0002160 , 0.0005313 , 5.0057000 , 5.0057500 , 0.2489000 , 0.0218399 ,
217  -0.2891600, 0.1794350 , 0.0011607 , 0.0007243 , 5.0054800 , 5.0059400 , 0.1409400 , -0.0889300 ,
218  0.3780580, -0.0605780, 0.0003887 , 0.0006382 , 5.0094100 , 5.0071100 , -0.1349100 , 0.0094002 ,
219  -0.1119370, -0.0205372, 0.0003893 , 0.0000034 , 5.0099600 , 5.0073600 , 0.0673701 , 0.0904398 ,
220  -0.1003550, 0.1026090 , -0.0003057 , 0.0000636 , 5.0070700 , 5.0071000 , -0.1667400 , -0.3741500 ,
221  0.4315830, 0.2629620 , 0.0005938 , 0.0007167 , 5.0070100 , 5.0063000 , 0.1928000 , -0.1254600 ,
222  0.0109732, 0.0740713 , 0.0006309 , 0.0000190 , 5.0058600 , 5.0063100 , -0.2017900 , 0.1079100 ,
223  -0.7136310, -0.1934810, 0.0006850 , 0.0000467 , 5.0024400 , 5.0048500 , -0.0082700 , 0.4051300 ,
224  0.2190650, -0.0636425, 0.0006076 , 0.0011785 , 5.0050100 , 5.0049200 , 0.0951700 , 0.1409100 ,
225  -0.0848179, 0.1256750 , -0.0000582 , -0.0000461 , 5.0078600 , 5.0059700 , 0.2675400 , -0.1691500};
226 
227  int npseudo=fdcPseudoVector.size();
228  for (unsigned int i = 0; i<npseudo; i++){
229  const DFDCPseudo* thisPseudo = fdcPseudoVector[i];
230 
231  int lpseudo = thisPseudo->wire->layer;
232  double zpseudo = thisPseudo->wire->origin.Z();
233  double upseudo = (thisPseudo->u)*10.;
234  double vpseudo = (thisPseudo->v)*10.;
235  double spseudo = (thisPseudo->s)*10.;
236  double wpseudo = (thisPseudo->w)*10.;
237 
238  int layer=lpseudo;
239 
240 // get constants for this layer
241 
242 float pxwl=wpseudo;
243 float pycl=-spseudo;
244 double u0 = constants[(layer-1)*8+0];
245 double d0 = constants[(layer-1)*8+1];
246 double udel = constants[(layer-1)*8+2];
247 double ddel = constants[(layer-1)*8+3];
248 double strip_pitch_u = constants[(layer-1)*8+4];
249 double strip_pitch_d = constants[(layer-1)*8+5];
250 double xshift = constants[(layer-1)*8+6];
251 double yshift = constants[(layer-1)*8+7];
252 
253 float v=vpseudo*strip_pitch_u/5.005-u0;
254 float u=upseudo*strip_pitch_d/5.005-d0;
255 float myangle=15.*M_PI/180.;
256 float strip_angle=75.*M_PI/180.;
257 float phi_u=strip_angle+ddel;
258 float phi_v=M_PI-strip_angle+udel;
259 float cosPhiU=cos(phi_u);
260 float cosPhiV=cos(phi_v);
261 float sinPhiU=sin(phi_u);
262 float sinPhiV=sin(phi_v);
263 float pxcl=(u*sinPhiV-v*sinPhiU)/(cosPhiV*sinPhiU-cosPhiU*sinPhiV);//-xshift;
264 pycl=(u*cosPhiV-v*cosPhiU)/(cosPhiU*sinPhiV-cosPhiV*sinPhiU);//+yshift;
265 
266 // The wire position projected from the wire
267 japp->RootWriteLock();
268 Hist3D[layer-1]->Fill(pxwl/10.,pycl/10.,(pxwl-pxcl)/10.);
269 japp->RootUnLock();
270 
271 char thisName[256];
272 sprintf(thisName, "Plane %.2i Wire Position Vs XY", layer);
273 Fill2DProfile("FDC_InternalAlignment", "Profile2D", thisName,
274  pxwl/10., pycl/10., (pxwl-pxcl)/10.,
275  thisName,
276  100, -50.,50., 100,-50., 50.);
277 }
278 */
279 return NOERROR;
280 }
281 
282 //------------------
283 // erun
284 //------------------
286 {
287  // This is called whenever the run number changes, before it is
288  // changed to give you a chance to clean up before processing
289  // events from the next run number.
290  return NOERROR;
291 }
292 
293 //------------------
294 // fini
295 //------------------
297 {
298  // Called before program exit after event processing is finished.
299  return NOERROR;
300 }
301 
int status
status word for pseudopoint
Definition: DFDCPseudo.h:94
sprintf(text,"Post KinFit Cut")
jerror_t evnt(jana::JEventLoop *eventLoop, uint64_t eventnumber)
Called every event.
const DFDCWire * wire
DFDCWire for this wire.
Definition: DFDCPseudo.h:92
class DFDCPseudo: definition for a reconstructed point in the FDC
Definition: DFDCPseudo.h:74
JApplication * japp
int layer
1-24
Definition: DFDCWire.h:16
double phi_v
rotation angles for cathode planes
Definition: DFDCPseudo.h:87
double w_c
Definition: DFDCPseudo.h:90
double u
Definition: DFDCPseudo.h:85
jerror_t brun(jana::JEventLoop *eventLoop, int32_t runnumber)
Called everytime a new run number is detected.
InitPlugin_t InitPlugin
jerror_t fini(void)
Called after last event of last event source has been processed.
double s
&lt; wire position computed from cathode data, assuming the avalanche occurs at the wire ...
Definition: DFDCPseudo.h:91
void Fill2DHistogram(const char *plugin, const char *directoryName, const char *name, const double valueX, const double valueY, const char *title, int nBinsX, double xmin, double xmax, int nBinsY, double ymin, double ymax, bool print=false)
double w
Definition: DFDCPseudo.h:89
double time
time corresponding to this pseudopoint.
Definition: DFDCPseudo.h:93
double sin(double)
jerror_t erun(void)
Called everytime run number changes, provided brun has been called.
int wire
1-N
Definition: DFDCWire.h:17
double phi_u
Definition: DFDCPseudo.h:87
void Fill2DProfile(const char *plugin, const char *directoryName, const char *name, const double valueX, const double valueY, const double valueZ, const char *title, int nBinsX, double xmin, double xmax, int nBinsY, double ymin, double ymax, bool print=false)
double v
centroid positions in the two cathode views
Definition: DFDCPseudo.h:85
jerror_t init(void)
Called once at program start.