Hall-D Software  alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
JEventProcessor_MilleFieldOff.cc
Go to the documentation of this file.
1 // $Id$
2 //
3 // File: JEventProcessor_MilleFieldOff.cc
4 // Created: Tue Jan 17 19:32:32 Local time zone must be set--see zic manual page 2017
5 // Creator: mstaib (on Linux egbert 2.6.32-642.6.2.el6.x86_64 x86_64)
6 //
7 
11 #include "CDC/DCDCTrackHit.h"
12 #include "CDC/DCDCWire.h"
13 #include "TDirectory.h"
14 using namespace jana;
15 
16 
17 // Routine used to create our JEventProcessor
18 #include <JANA/JApplication.h>
19 #include <JANA/JFactory.h>
20 extern "C"{
21 void InitPlugin(JApplication *app){
22  InitJANAPlugin(app);
23  app->AddProcessor(new JEventProcessor_MilleFieldOff());
24 }
25 } // "C"
26 
27 
28 //------------------
29 // JEventProcessor_MilleFieldOff (Constructor)
30 //------------------
32 {
33 
34 }
35 
36 //------------------
37 // ~JEventProcessor_MilleFieldOff (Destructor)
38 //------------------
40 {
41 
42 }
43 
44 //------------------
45 // init
46 //------------------
48 {
49  // This is called once at program startup.
50  milleWriter = new Mille("nofield_mille_out.mil");
51 
52  gDirectory->mkdir("AlignmentConstants");
53  gDirectory->cd("AlignmentConstants");
54  // We need the constants used for this iteration
55  // Use a TProfile to avoid problems adding together multiple root files...
56  HistCurrentConstantsCDC = new TProfile("CDCAlignmentConstants", "Constants Used for CDC Alignment (In MILLEPEDE Order)", 16000 ,0.5, 16000.5);
57  HistCurrentConstantsFDC = new TProfile("FDCAlignmentConstants", "Constants Used for FDC Alignment (In MILLEPEDE Order)", 26000 ,0.5, 26000.5);
58 
59  gDirectory->cd("..");
60 
61  return NOERROR;
62 }
63 
64 //------------------
65 // brun
66 //------------------
67 jerror_t JEventProcessor_MilleFieldOff::brun(JEventLoop *eventLoop, int32_t runnumber)
68 {
69  // Get the current set of constants and sve them in the histogram
70  // This is called whenever the run number changes
71  // Check for magnetic field
72  DApplication* dapp=dynamic_cast<DApplication*>(eventLoop->GetJApplication());
73  bool dIsNoFieldFlag = (dynamic_cast<const DMagneticFieldMapNoField*>(dapp->GetBfield(runnumber)) != NULL);
74 
75  // This plugin is designed for field off data. If this is used for field on data, Abort...
76  if (!dIsNoFieldFlag){
77  jerr << " Plugin FDC_MilleFieldOff Must be run with zero magnetic field!!! Aborting " << endl;
78  jerr << " Use -PBFIELD_TYPE=NoField " << endl;
79  japp->Quit();
80  }
81 
82  // Store the current values of the alignment constants
83  JCalibration * jcalib = eventLoop->GetJCalibration();
84  vector<map<string,double> >vals;
85  if (jcalib->Get("CDC/global_alignment",vals)==false){
86  map<string,double> &row = vals[0];
87  // Get the offsets from the calibration database
88  HistCurrentConstantsCDC->Fill(1,row["dX"]);
89  HistCurrentConstantsCDC->Fill(2,row["dY"]);
90  HistCurrentConstantsCDC->Fill(3,row["dZ"]);
91  HistCurrentConstantsCDC->Fill(4,row["dPhiX"]);
92  HistCurrentConstantsCDC->Fill(5,row["dPhiY"]);
93  HistCurrentConstantsCDC->Fill(6,row["dPhiZ"]);
94  }
95 
96  if (jcalib->Get("CDC/wire_alignment",vals)==false){
97  for(unsigned int i=0; i<vals.size(); i++){
98  map<string,double> &row = vals[i];
99  // Get the offsets from the calibration database
100  HistCurrentConstantsCDC->Fill(1000+(i*4+1),row["dxu"]);
101  HistCurrentConstantsCDC->Fill(1000+(i*4+2),row["dyu"]);
102  HistCurrentConstantsCDC->Fill(1000+(i*4+3),row["dxd"]);
103  HistCurrentConstantsCDC->Fill(1000+(i*4+4),row["dyd"]);
104  }
105  }
106 
107  if (jcalib->Get("FDC/cathode_alignment",vals)==false){
108  for(unsigned int i=0; i<vals.size(); i++){
109  map<string,double> &row = vals[i];
110  // Get the offsets from the calibration database
111  HistCurrentConstantsFDC->Fill((i+1)*1000+101,row["dU"]);
112  HistCurrentConstantsFDC->Fill((i+1)*1000+102,row["dV"]);
113  HistCurrentConstantsFDC->Fill((i+1)*1000+103,row["dPhiU"]);
114  HistCurrentConstantsFDC->Fill((i+1)*1000+104,row["dPhiV"]);
115  }
116  }
117 
118  if (jcalib->Get("FDC/cell_offsets",vals)==false){
119  for(unsigned int i=0; i<vals.size(); i++){
120  map<string,double> &row = vals[i];
121  // Get the offsets from the calibration database
122  HistCurrentConstantsFDC->Fill((i+1)*1000+1,row["xshift"]);
123  HistCurrentConstantsFDC->Fill((i+1)*1000+100,row["yshift"]);
124  }
125  }
126 
127  if (jcalib->Get("FDC/cell_rotations",vals)==false){
128  for(unsigned int i=0; i<vals.size(); i++){
129  map<string,double> &row = vals[i];
130  HistCurrentConstantsFDC->Fill((i+1)*1000+2,row["dPhiX"]);
131  HistCurrentConstantsFDC->Fill((i+1)*1000+3,row["dPhiY"]);
132  HistCurrentConstantsFDC->Fill((i+1)*1000+4,row["dPhiZ"]);
133  }
134  }
135  if (jcalib->Get("FDC/wire_alignment",vals)==false){
136  for(unsigned int i=0; i<vals.size(); i++){
137  map<string,double> &row = vals[i];
138  // Get the offsets from the calibration database
139  HistCurrentConstantsFDC->Fill((i+1)*1000+5,row["dZ"]);
140  }
141  }
142 
143  if (jcalib->Get("FDC/package1/strip_gains_v2",vals)==false){
144  for(unsigned int i=0; i<vals.size(); i++){
145  map<string,double> &row = vals[i];
146  // Get the offsets from the calibration database
147  char name[32];
148  for (unsigned int j=1; j<= 216; j++)
149  {
150  sprintf(name,"strip%i", j);
151  if (i%2){ // V
152  HistCurrentConstantsFDC->Fill((i/2+1)*1000+600+j,row[name]);
153  }
154  else {// U
155  HistCurrentConstantsFDC->Fill((i/2+1)*1000+300+j,row[name]);
156  }
157  }
158  }
159  }
160 
161  if (jcalib->Get("FDC/package2/strip_gains_v2",vals)==false){
162  for(unsigned int i=0; i<vals.size(); i++){
163  map<string,double> &row = vals[i];
164  // Get the offsets from the calibration database
165  char name[32];
166  for (unsigned int j=1; j<= 216; j++)
167  {
168  sprintf(name,"strip%i", j);
169  if (i%2){ // V
170  HistCurrentConstantsFDC->Fill((i/2+7)*1000+600+j,row[name]);
171  }
172  else // U
173  HistCurrentConstantsFDC->Fill((i/2+7)*1000+300+j,row[name]);
174  }
175  }
176  }
177 
178  if (jcalib->Get("FDC/package3/strip_gains_v2",vals)==false){
179  for(unsigned int i=0; i<vals.size(); i++){
180  map<string,double> &row = vals[i];
181  // Get the offsets from the calibration database
182  char name[32];
183  for (unsigned int j=1; j<= 216; j++)
184  {
185  sprintf(name,"strip%i", j);
186  if (i%2){ // V
187  HistCurrentConstantsFDC->Fill((i/2+13)*1000+600+j,row[name]);
188  }
189  else // U
190  HistCurrentConstantsFDC->Fill((i/2+13)*1000+300+j,row[name]);
191  }
192  }
193  }
194 
195  if (jcalib->Get("FDC/package4/strip_gains_v2",vals)==false){
196  for(unsigned int i=0; i<vals.size(); i++){
197  map<string,double> &row = vals[i];
198  // Get the offsets from the calibration database
199  char name[32];
200  for (unsigned int j=1; j<= 216; j++)
201  {
202  sprintf(name,"strip%i", j);
203  if (i%2){ // V
204  HistCurrentConstantsFDC->Fill((i/2+19)*1000+600+j,row[name]);
205  }
206  else // U
207  HistCurrentConstantsFDC->Fill((i/2+19)*1000+300+j,row[name]);
208  }
209  }
210  }
211 
212  if (jcalib->Get("FDC/strip_pitches_v2",vals)==false){
213  for(unsigned int i=0; i<vals.size(); i++){
214  map<string,double> &row = vals[i];
215  HistCurrentConstantsFDC->Fill((i+1)*1000 + 200, row["U_SP_1"]);
216  HistCurrentConstantsFDC->Fill((i+1)*1000 + 201, row["U_G_1"]);
217  HistCurrentConstantsFDC->Fill((i+1)*1000 + 202, row["U_SP_2"]);
218  HistCurrentConstantsFDC->Fill((i+1)*1000 + 203, row["U_G_2"]);
219  HistCurrentConstantsFDC->Fill((i+1)*1000 + 204, row["U_SP_3"]);
220  HistCurrentConstantsFDC->Fill((i+1)*1000 + 205, row["V_SP_1"]);
221  HistCurrentConstantsFDC->Fill((i+1)*1000 + 206, row["V_G_1"]);
222  HistCurrentConstantsFDC->Fill((i+1)*1000 + 207, row["V_SP_2"]);
223  HistCurrentConstantsFDC->Fill((i+1)*1000 + 208, row["V_G_2"]);
224  HistCurrentConstantsFDC->Fill((i+1)*1000 + 209, row["V_SP_3"]);
225  }
226  }
227 
228  return NOERROR;
229 }
230 
231 //------------------
232 // evnt
233 //------------------
234 jerror_t JEventProcessor_MilleFieldOff::evnt(JEventLoop *loop, uint64_t eventnumber)
235 {
236  int straw_offset[29] = {0,0,42,84,138,192,258,324,404,484,577,670,776,882,1005,1128,1263,1398,1544,1690,1848,2006,2176,2346,2528,2710,2907,3104,3313};
237  // Loop over the tracks, get the tracking pulls, and fill some histograms. Easy peasy
238  vector<const DTrackCandidate *> trackVector;
239  loop->Get(trackVector,"StraightLine");
240 
241  for (size_t i = 0; i < trackVector.size(); i++){
242  const DTrackCandidate *track = trackVector[i];
243 
244  double trackingFOM = TMath::Prob(track->chisq, track->Ndof);
245 
246  // Some quality cuts for the tracks we will use
247  // Keep this minimal for now and investigate later
248  double trackingFOMCut = 0.001;
249  int trackingNDFCut = 18;
250 
251  if(trackingFOM < trackingFOMCut) continue;
252  if(track->Ndof < trackingNDFCut) continue;
253 
254  //double phi = bestHypothesis->momentum().Phi()*TMath::RadToDeg();
255  //double theta = bestHypothesis->momentum().Theta()*TMath::RadToDeg();
256 
257  if (!track->IsSmoothed) continue;
258 
259  vector<DTrackFitter::pull_t> pulls = track->pulls;
260  if (pulls.size() > 75) continue; // We don't want reallly long tracks
261 
262  // Determine TrackType
263  bool isCDCOnly=true; //bool isFDCOnly=true;
264  bool anyNaN = false;
265  for (size_t iPull = 0; iPull < pulls.size(); iPull++){
266  const DCDCTrackHit *cdc_hit = pulls[iPull].cdc_hit;
267  //const DFDCPseudo *fdc_hit = pulls[iPull].fdc_hit;
268  float resi = pulls[iPull].resi; // residual of measurement
269  float err = pulls[iPull].err;
270  float resic = pulls[iPull].resic; // residual for FDC cathode measuremtns
271  float errc = pulls[iPull].errc;
272  if (cdc_hit == NULL) isCDCOnly=false;
273  if (resi != resi || err != err || resic != resic || errc != errc || !isfinite(resi) || !isfinite(resic)) anyNaN=true;
274  // If any residual or error is truely huge, we dont really want this event. Cut at 10 cm.
275  if (fabs(resi)>10. || err > 10. || fabs(resic) > 10. || errc > 10.) anyNaN=true;
276  //if (fdc_hit == NULL) isFDCOnly=false;
277  }
278 
279  if (anyNaN) continue;
280 
281  // Need at least as many degrees of freedon as alignment parameters
282  if (!isCDCOnly) {
283  if(track->Ndof < 26) continue;
284  }
285 
286  japp->RootWriteLock(); // Just use the root lock as a temporary
287  for (size_t iPull = 0; iPull < pulls.size(); iPull++){
288  // Here is all of the information currently stored in the pulls from the fit
289  // From TRACKING/DTrackFitter.h
290  float resi = pulls[iPull].resi; // residual of measurement
291  float err = pulls[iPull].err; // estimated error of measurement
292  //double s = pulls[iPull].s;
293  //double tdrift = pulls[iPull].tdrift; // drift time of this measurement
294  //double d = pulls[iPull].d; // doca to wire
295  const DCDCTrackHit *cdc_hit = pulls[iPull].cdc_hit;
296  const DFDCPseudo *fdc_hit = pulls[iPull].fdc_hit;
297  //double docaphi = pulls[iPull].docaphi; // phi of doca in CDC straws
298  //double z = pulls[iPull].z;// z position at doca
299  //double tcorr = pulls[iPull].tcorr; // drift time with correction for B
300  float resic = pulls[iPull].resic; // residual for FDC cathode measuremtns
301  float errc = pulls[iPull].errc;
302 
303  vector<double> trackDerivatives = pulls[iPull].trackDerivatives;
304  if (trackDerivatives.size()==0) jerr << "Track derivative size == 0 ???" << endl;
305 
306  if (fdc_hit != NULL && fdc_hit->status == 6) {
307  // Add fdc hit
308  DFDCPseudo *thisHit = const_cast<DFDCPseudo *>(fdc_hit);
309 
310  vector<double> pseudoAlignmentDerivatives = thisHit->GetFDCPseudoAlignmentDerivatives();
311  vector<double> fdcStripGainDerivatives = thisHit->GetFDCStripGainDerivatives();
312  vector<double> fdcStripPitchDerivatives = thisHit->GetFDCStripPitchDerivatives();
313 
314  unsigned int layerOffset = 100000 + thisHit->wire->layer * 1000;
315 
316  //Now we need to fill the Mille structure once for our wire measurment and once for our cathode measurement
317  const int NLC = 4; // Number of local parameters
318  const int NGL_W = 5; // Number of global parameters for wire alignment
319  float localDerW[NLC];
320  float globalDerW[NGL_W];
321  int labelW[NGL_W];
322 
323  localDerW[0]=trackDerivatives[FDCTrackD::dDOCAW_dx]; localDerW[1]=trackDerivatives[FDCTrackD::dDOCAW_dy];
324  localDerW[2]=trackDerivatives[FDCTrackD::dDOCAW_dtx]; localDerW[3]=trackDerivatives[FDCTrackD::dDOCAW_dty];
325 
326  globalDerW[0] = trackDerivatives[FDCTrackD::dDOCAW_dDeltaX]; labelW[0] = layerOffset + 1;
327  globalDerW[1] = trackDerivatives[FDCTrackD::dDOCAW_dDeltaPhiX]; labelW[1] = layerOffset + 2;
328  globalDerW[2] = trackDerivatives[FDCTrackD::dDOCAW_dDeltaPhiY]; labelW[2] = layerOffset + 3;
329  globalDerW[3] = trackDerivatives[FDCTrackD::dDOCAW_dDeltaPhiZ]; labelW[3] = layerOffset + 4;
330  globalDerW[4] = trackDerivatives[FDCTrackD::dDOCAW_dDeltaZ]; labelW[4] = layerOffset + 5;
331 
332  milleWriter->mille(NLC, localDerW, NGL_W, globalDerW, labelW, resi, err);
333 
334  // Now for the cathode measurement, there are more alignment parameters.
335  const int NGL_C = 26; // Number of global parameters for wire alignment
336  float localDerC[NLC];
337  float globalDerC[NGL_C];
338  int labelC[NGL_C];
339 
340  localDerC[0]=trackDerivatives[FDCTrackD::dDOCAC_dx]; localDerC[1]=trackDerivatives[FDCTrackD::dDOCAC_dy];
341  localDerC[2]=trackDerivatives[FDCTrackD::dDOCAC_dtx]; localDerC[3]=trackDerivatives[FDCTrackD::dDOCAC_dty];
342 
343  globalDerC[0] = -1.0; labelC[0] = layerOffset + 100;
344  globalDerC[1] = trackDerivatives[FDCTrackD::dDOCAC_dDeltaX]; labelC[1] = layerOffset + 1;
345  globalDerC[2] = trackDerivatives[FDCTrackD::dDOCAC_dDeltaPhiX]; labelC[2] = layerOffset + 2;
346  globalDerC[3] = trackDerivatives[FDCTrackD::dDOCAC_dDeltaPhiY]; labelC[3] = layerOffset + 3;
347  globalDerC[4] = trackDerivatives[FDCTrackD::dDOCAC_dDeltaPhiZ]; labelC[4] = layerOffset + 4;
348  globalDerC[5] = trackDerivatives[FDCTrackD::dDOCAC_dDeltaZ]; labelC[5] = layerOffset + 5;
349 
350  // Cathode U and V offsets
351  globalDerC[6] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaU]; labelC[6] = layerOffset + 101;
352  globalDerC[7] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaV]; labelC[7] = layerOffset + 102;
353  globalDerC[8] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaPhiU]; labelC[8] = layerOffset + 103;
354  globalDerC[9] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaPhiV]; labelC[9] = layerOffset + 104;
355 
356  // Strip Pitch Calibration
357  globalDerC[10] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaU]*fdcStripPitchDerivatives[0]; labelC[10] = layerOffset + 200;
358  globalDerC[11] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaU]*fdcStripPitchDerivatives[1]; labelC[11] = layerOffset + 201;
359  globalDerC[12] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaU]*fdcStripPitchDerivatives[2]; labelC[12] = layerOffset + 202;
360  globalDerC[13] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaU]*fdcStripPitchDerivatives[3]; labelC[13] = layerOffset + 203;
361  globalDerC[14] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaU]*fdcStripPitchDerivatives[4]; labelC[14] = layerOffset + 204;
362  globalDerC[15] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaV]*fdcStripPitchDerivatives[5]; labelC[15] = layerOffset + 205;
363  globalDerC[16] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaV]*fdcStripPitchDerivatives[6]; labelC[16] = layerOffset + 206;
364  globalDerC[17] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaV]*fdcStripPitchDerivatives[7]; labelC[17] = layerOffset + 207;
365  globalDerC[18] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaV]*fdcStripPitchDerivatives[8]; labelC[18] = layerOffset + 208;
366  globalDerC[19] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaV]*fdcStripPitchDerivatives[9]; labelC[19] = layerOffset + 209;
367 
368  // Strip Gain Calibration
369  vector<const DFDCCathodeCluster *> cathodeClusters;
370  thisHit->Get(cathodeClusters);
371  unsigned int gainLabels[6]={};
372  for(unsigned int j = 0; j< cathodeClusters.size(); j++){
373  if(cathodeClusters[j]->plane == 3) { // U strips
374  unsigned int k=0;
375  for (; k < cathodeClusters[j]->members.size(); k++){
376  if (thisHit->cluster_u.index(0) == cathodeClusters[j]->members[k]->element) break;
377  }
378  gainLabels[0] = cathodeClusters[j]->members[k]->type != 3 ? thisHit->cluster_u.index(0) : thisHit->cluster_u.index(0)+108;
379  gainLabels[1] = cathodeClusters[j]->members[k+1]->type != 3 ? thisHit->cluster_u.index(1) : thisHit->cluster_u.index(1)+108;
380  gainLabels[2] = cathodeClusters[j]->members[k+2]->type != 3 ? thisHit->cluster_u.index(2) : thisHit->cluster_u.index(2)+108;
381  }
382  else if (cathodeClusters[j]->plane == 1) { // V strips
383  unsigned int k=0;
384  for (; k < cathodeClusters[j]->members.size(); k++){
385  if (thisHit->cluster_v.index(0) == cathodeClusters[j]->members[k]->element) break;
386  }
387  gainLabels[3] = cathodeClusters[j]->members[k]->type != 3 ? thisHit->cluster_v.index(0) : thisHit->cluster_v.index(0)+108;
388  gainLabels[4] = cathodeClusters[j]->members[k+1]->type != 3 ? thisHit->cluster_v.index(1) : thisHit->cluster_v.index(1)+108;
389  gainLabels[5] = cathodeClusters[j]->members[k+2]->type != 3 ? thisHit->cluster_v.index(2) : thisHit->cluster_v.index(2)+108;
390  }
391  }
392 
393  globalDerC[20] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaU]*fdcStripGainDerivatives[0]; labelC[20] = layerOffset + 300 + gainLabels[0];
394  globalDerC[21] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaU]*fdcStripGainDerivatives[1]; labelC[21] = layerOffset + 300 + gainLabels[1];
395  globalDerC[22] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaU]*fdcStripGainDerivatives[2]; labelC[22] = layerOffset + 300 + gainLabels[2];
396  globalDerC[23] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaV]*fdcStripGainDerivatives[3]; labelC[23] = layerOffset + 600 + gainLabels[3];
397  globalDerC[24] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaV]*fdcStripGainDerivatives[4]; labelC[24] = layerOffset + 600 + gainLabels[4];
398  globalDerC[25] = -pseudoAlignmentDerivatives[FDCPseudoD::dSddeltaV]*fdcStripGainDerivatives[5]; labelC[25] = layerOffset + 600 + gainLabels[5];
399 
400  milleWriter->mille(NLC, localDerC, NGL_C, globalDerC, labelC, resic, errc);
401  }
402 
403  if (cdc_hit != NULL){
404 
405  const DCDCWire *constWire = cdc_hit->wire;
406  DCDCWire *thisWire = const_cast<DCDCWire *>(constWire);
407 
408  vector<double> wireDerivatives = thisWire->GetCDCWireDerivatives();
409 
410  //Now we need to fill the Mille structure once for our wire measurment and once for our cathode measurement
411  const int NLC = 4; // Number of local parameters
412  const int NGL = 10; // Number of global parameters for wire alignment
413  float localDer[NLC];
414  float globalDer[NGL];
415  int label[NGL];
416 
417  localDer[0]=trackDerivatives[CDCTrackD::dDOCAdS0]; localDer[1]=trackDerivatives[CDCTrackD::dDOCAdS1];
418  localDer[2]=trackDerivatives[CDCTrackD::dDOCAdS2]; localDer[3]=trackDerivatives[CDCTrackD::dDOCAdS3];
419 
420  //if (isCDCOnly){ // Global shifts will not affect residuals
421  if (false) { // test
422  globalDer[0]=0.0; label[0]=1;
423  globalDer[1]=0.0; label[0]=2;
424  globalDer[2]=0.0; label[0]=3;
425  globalDer[3]=0.0; label[0]=4;
426  globalDer[4]=0.0; label[0]=5;
427  globalDer[5]=0.0; label[0]=6;
428  }
429  else{
430  globalDer[0]=trackDerivatives[CDCTrackD::dDOCAdOriginX]*wireDerivatives[CDCWireD::dOriginXddeltaX]
431  +trackDerivatives[CDCTrackD::dDOCAdOriginY]*wireDerivatives[CDCWireD::dOriginYddeltaX]
432  +trackDerivatives[CDCTrackD::dDOCAdOriginZ]*wireDerivatives[CDCWireD::dOriginZddeltaX]
433  +trackDerivatives[CDCTrackD::dDOCAdDirX]*wireDerivatives[CDCWireD::dDirXddeltaX]
434  +trackDerivatives[CDCTrackD::dDOCAdDirY]*wireDerivatives[CDCWireD::dDirYddeltaX]
435  +trackDerivatives[CDCTrackD::dDOCAdDirZ]*wireDerivatives[CDCWireD::dDirZddeltaX];
436  label[0]=1;
437 
438  globalDer[1]=trackDerivatives[CDCTrackD::dDOCAdOriginX]*wireDerivatives[CDCWireD::dOriginXddeltaY]
439  +trackDerivatives[CDCTrackD::dDOCAdOriginY]*wireDerivatives[CDCWireD::dOriginYddeltaY]
440  +trackDerivatives[CDCTrackD::dDOCAdOriginZ]*wireDerivatives[CDCWireD::dOriginZddeltaY]
441  +trackDerivatives[CDCTrackD::dDOCAdDirX]*wireDerivatives[CDCWireD::dDirXddeltaY]
442  +trackDerivatives[CDCTrackD::dDOCAdDirY]*wireDerivatives[CDCWireD::dDirYddeltaY]
443  +trackDerivatives[CDCTrackD::dDOCAdDirZ]*wireDerivatives[CDCWireD::dDirZddeltaY];
444  label[1]=2;
445 
446  globalDer[2]=trackDerivatives[CDCTrackD::dDOCAdOriginX]*wireDerivatives[CDCWireD::dOriginXddeltaZ]
447  +trackDerivatives[CDCTrackD::dDOCAdOriginY]*wireDerivatives[CDCWireD::dOriginYddeltaZ]
448  +trackDerivatives[CDCTrackD::dDOCAdOriginZ]*wireDerivatives[CDCWireD::dOriginZddeltaZ]
449  +trackDerivatives[CDCTrackD::dDOCAdDirX]*wireDerivatives[CDCWireD::dDirXddeltaZ]
450  +trackDerivatives[CDCTrackD::dDOCAdDirY]*wireDerivatives[CDCWireD::dDirYddeltaZ]
451  +trackDerivatives[CDCTrackD::dDOCAdDirZ]*wireDerivatives[CDCWireD::dDirZddeltaZ];
452  label[2]=3;
453  if (false){
454  jout << " Dumping deltaZ derivatives============ Wire " << thisWire->ring << " Straw " << thisWire->straw << endl;
455  jout << " Total = " << globalDer[2] << endl;
456  jout << "trackDerivatives[CDCTrackD::dDOCAdOriginX] " << trackDerivatives[CDCTrackD::dDOCAdOriginX] << " wireDerivatives[CDCWireD::dOriginXddeltaZ] " << wireDerivatives[CDCWireD::dOriginXddeltaZ] << endl;
457  jout << "trackDerivatives[CDCTrackD::dDOCAdOriginY] " << trackDerivatives[CDCTrackD::dDOCAdOriginY] << " wireDerivatives[CDCWireD::dOriginYddeltaZ] " << wireDerivatives[CDCWireD::dOriginYddeltaZ] << endl;
458  jout << "trackDerivatives[CDCTrackD::dDOCAdOriginZ] " << trackDerivatives[CDCTrackD::dDOCAdOriginZ] << " wireDerivatives[CDCWireD::dOriginZddeltaZ] " << wireDerivatives[CDCWireD::dOriginZddeltaZ] << endl;
459  jout << "trackDerivatives[CDCTrackD::dDOCAdDirX] " << trackDerivatives[CDCTrackD::dDOCAdDirX] << " wireDerivatives[CDCWireD::dDirXddeltaZ] " << wireDerivatives[CDCWireD::dDirXddeltaZ] << endl;
460  jout << "trackDerivatives[CDCTrackD::dDOCAdDirY] " << trackDerivatives[CDCTrackD::dDOCAdDirY] << " wireDerivatives[CDCWireD::dDirYddeltaZ] " << wireDerivatives[CDCWireD::dDirYddeltaZ] << endl;
461  jout << "trackDerivatives[CDCTrackD::dDOCAdDirZ] " << trackDerivatives[CDCTrackD::dDOCAdDirZ] << " wireDerivatives[CDCWireD::dDirZddeltaZ] " << wireDerivatives[CDCWireD::dDirZddeltaZ] << endl;
462  }
463 
464  globalDer[3]=trackDerivatives[CDCTrackD::dDOCAdOriginX]*wireDerivatives[CDCWireD::dOriginXddeltaPhiX]
465  +trackDerivatives[CDCTrackD::dDOCAdOriginY]*wireDerivatives[CDCWireD::dOriginYddeltaPhiX]
466  +trackDerivatives[CDCTrackD::dDOCAdOriginZ]*wireDerivatives[CDCWireD::dOriginZddeltaPhiX]
467  +trackDerivatives[CDCTrackD::dDOCAdDirX]*wireDerivatives[CDCWireD::dDirXddeltaPhiX]
468  +trackDerivatives[CDCTrackD::dDOCAdDirY]*wireDerivatives[CDCWireD::dDirYddeltaPhiX]
469  +trackDerivatives[CDCTrackD::dDOCAdDirZ]*wireDerivatives[CDCWireD::dDirZddeltaPhiX];
470  label[3]=4;
471 
472  globalDer[4]=trackDerivatives[CDCTrackD::dDOCAdOriginX]*wireDerivatives[CDCWireD::dOriginXddeltaPhiY]
473  +trackDerivatives[CDCTrackD::dDOCAdOriginY]*wireDerivatives[CDCWireD::dOriginYddeltaPhiY]
474  +trackDerivatives[CDCTrackD::dDOCAdOriginZ]*wireDerivatives[CDCWireD::dOriginZddeltaPhiY]
475  +trackDerivatives[CDCTrackD::dDOCAdDirX]*wireDerivatives[CDCWireD::dDirXddeltaPhiY]
476  +trackDerivatives[CDCTrackD::dDOCAdDirY]*wireDerivatives[CDCWireD::dDirYddeltaPhiY]
477  +trackDerivatives[CDCTrackD::dDOCAdDirZ]*wireDerivatives[CDCWireD::dDirZddeltaPhiY];
478  label[4]=5;
479 
480  globalDer[5]=trackDerivatives[CDCTrackD::dDOCAdOriginX]*wireDerivatives[CDCWireD::dOriginXddeltaPhiZ]
481  +trackDerivatives[CDCTrackD::dDOCAdOriginY]*wireDerivatives[CDCWireD::dOriginYddeltaPhiZ]
482  +trackDerivatives[CDCTrackD::dDOCAdOriginZ]*wireDerivatives[CDCWireD::dOriginZddeltaPhiZ]
483  +trackDerivatives[CDCTrackD::dDOCAdDirX]*wireDerivatives[CDCWireD::dDirXddeltaPhiZ]
484  +trackDerivatives[CDCTrackD::dDOCAdDirY]*wireDerivatives[CDCWireD::dDirYddeltaPhiZ]
485  +trackDerivatives[CDCTrackD::dDOCAdDirZ]*wireDerivatives[CDCWireD::dDirZddeltaPhiZ];
486  label[5]=6;
487  if (false){
488  jout << " Dumping deltaPhiZ derivatives============ Ring " << thisWire->ring << " Straw " << thisWire->straw << endl;
489  jout << " Total = " << globalDer[5] << endl;
490  jout << "trackDerivatives[CDCTrackD::dDOCAdOriginX] " << trackDerivatives[CDCTrackD::dDOCAdOriginX] << " wireDerivatives[CDCWireD::dOriginXddeltaPhiZ] " << wireDerivatives[CDCWireD::dOriginXddeltaPhiZ] << endl;
491  jout << "trackDerivatives[CDCTrackD::dDOCAdOriginY] " << trackDerivatives[CDCTrackD::dDOCAdOriginY] << " wireDerivatives[CDCWireD::dOriginYddeltaPhiZ] " << wireDerivatives[CDCWireD::dOriginYddeltaPhiZ] << endl;
492  jout << "trackDerivatives[CDCTrackD::dDOCAdOriginZ] " << trackDerivatives[CDCTrackD::dDOCAdOriginZ] << " wireDerivatives[CDCWireD::dOriginZddeltaPhiZ] " << wireDerivatives[CDCWireD::dOriginZddeltaPhiZ] << endl;
493  jout << "trackDerivatives[CDCTrackD::dDOCAdDirX] " << trackDerivatives[CDCTrackD::dDOCAdDirX] << " wireDerivatives[CDCWireD::dDirXddeltaPhiZ] " << wireDerivatives[CDCWireD::dDirXddeltaPhiZ] << endl;
494  jout << "trackDerivatives[CDCTrackD::dDOCAdDirY] " << trackDerivatives[CDCTrackD::dDOCAdDirY] << " wireDerivatives[CDCWireD::dDirYddeltaPhiZ] " << wireDerivatives[CDCWireD::dDirYddeltaPhiZ] << endl;
495  jout << "trackDerivatives[CDCTrackD::dDOCAdDirZ] " << trackDerivatives[CDCTrackD::dDOCAdDirZ] << " wireDerivatives[CDCWireD::dDirZddeltaPHiZ] " << wireDerivatives[CDCWireD::dDirZddeltaPhiZ] << endl;
496  }
497  }
498  globalDer[6]=trackDerivatives[CDCTrackD::dDOCAdOriginX]*wireDerivatives[CDCWireD::dOriginXddeltaXu]
499  +trackDerivatives[CDCTrackD::dDOCAdOriginY]*wireDerivatives[CDCWireD::dOriginYddeltaXu]
500  +trackDerivatives[CDCTrackD::dDOCAdOriginZ]*wireDerivatives[CDCWireD::dOriginZddeltaXu]
501  +trackDerivatives[CDCTrackD::dDOCAdDirX]*wireDerivatives[CDCWireD::dDirXddeltaXu]
502  +trackDerivatives[CDCTrackD::dDOCAdDirY]*wireDerivatives[CDCWireD::dDirYddeltaXu]
503  +trackDerivatives[CDCTrackD::dDOCAdDirZ]*wireDerivatives[CDCWireD::dDirZddeltaXu];
504  label[6]=1000 + (straw_offset[thisWire->ring]+(thisWire->straw-1))*4 + 1;
505 
506  if (false){
507  jout << " Dumping deltaXu derivatives============ Wire " << thisWire->ring << " Straw " << thisWire->straw << endl;
508  jout << " Total = " << globalDer[6] << endl;
509  jout << "trackDerivatives[CDCTrackD::dDOCAdOriginX] " << trackDerivatives[CDCTrackD::dDOCAdOriginX] << " wireDerivatives[CDCWireD::dOriginXddeltaXu]" << wireDerivatives[CDCWireD::dOriginXddeltaXu] << endl;
510  jout << "trackDerivatives[CDCTrackD::dDOCAdOriginY] " << trackDerivatives[CDCTrackD::dDOCAdOriginY] << " wireDerivatives[CDCWireD::dOriginYddeltaXu]" << wireDerivatives[CDCWireD::dOriginYddeltaXu] << endl;
511  jout << "trackDerivatives[CDCTrackD::dDOCAdOriginZ] " << trackDerivatives[CDCTrackD::dDOCAdOriginZ] << " wireDerivatives[CDCWireD::dOriginZddeltaXu]" << wireDerivatives[CDCWireD::dOriginZddeltaXu] << endl;
512  jout << "trackDerivatives[CDCTrackD::dDOCAdDirX] " << trackDerivatives[CDCTrackD::dDOCAdDirX] << " wireDerivatives[CDCWireD::dDirXddeltaXu]" << wireDerivatives[CDCWireD::dDirXddeltaXu] << endl;
513  jout << "trackDerivatives[CDCTrackD::dDOCAdDirY] " << trackDerivatives[CDCTrackD::dDOCAdDirY] << " wireDerivatives[CDCWireD::dDirYddeltaXu]" << wireDerivatives[CDCWireD::dDirYddeltaXu] << endl;
514  jout << "trackDerivatives[CDCTrackD::dDOCAdDirZ] " << trackDerivatives[CDCTrackD::dDOCAdDirZ] << " wireDerivatives[CDCWireD::dDirZddeltaXu]" << wireDerivatives[CDCWireD::dDirZddeltaXu] << endl;
515  }
516 
517  globalDer[7]=trackDerivatives[CDCTrackD::dDOCAdOriginX]*wireDerivatives[CDCWireD::dOriginXddeltaYu]
518  +trackDerivatives[CDCTrackD::dDOCAdOriginY]*wireDerivatives[CDCWireD::dOriginYddeltaYu]
519  +trackDerivatives[CDCTrackD::dDOCAdOriginZ]*wireDerivatives[CDCWireD::dOriginZddeltaYu]
520  +trackDerivatives[CDCTrackD::dDOCAdDirX]*wireDerivatives[CDCWireD::dDirXddeltaYu]
521  +trackDerivatives[CDCTrackD::dDOCAdDirY]*wireDerivatives[CDCWireD::dDirYddeltaYu]
522  +trackDerivatives[CDCTrackD::dDOCAdDirZ]*wireDerivatives[CDCWireD::dDirZddeltaYu];
523  label[7]=1000 + (straw_offset[thisWire->ring]+(thisWire->straw-1))*4 + 2;
524 
525  globalDer[8]=trackDerivatives[CDCTrackD::dDOCAdOriginX]*wireDerivatives[CDCWireD::dOriginXddeltaXd]
526  +trackDerivatives[CDCTrackD::dDOCAdOriginY]*wireDerivatives[CDCWireD::dOriginYddeltaXd]
527  +trackDerivatives[CDCTrackD::dDOCAdOriginZ]*wireDerivatives[CDCWireD::dOriginZddeltaXd]
528  +trackDerivatives[CDCTrackD::dDOCAdDirX]*wireDerivatives[CDCWireD::dDirXddeltaXd]
529  +trackDerivatives[CDCTrackD::dDOCAdDirY]*wireDerivatives[CDCWireD::dDirYddeltaXd]
530  +trackDerivatives[CDCTrackD::dDOCAdDirZ]*wireDerivatives[CDCWireD::dDirZddeltaXd];
531  label[8]=1000 + (straw_offset[thisWire->ring]+(thisWire->straw-1))*4 + 3;
532 
533  globalDer[9]=trackDerivatives[CDCTrackD::dDOCAdOriginX]*wireDerivatives[CDCWireD::dOriginXddeltaYd]
534  +trackDerivatives[CDCTrackD::dDOCAdOriginY]*wireDerivatives[CDCWireD::dOriginYddeltaYd]
535  +trackDerivatives[CDCTrackD::dDOCAdOriginZ]*wireDerivatives[CDCWireD::dOriginZddeltaYd]
536  +trackDerivatives[CDCTrackD::dDOCAdDirX]*wireDerivatives[CDCWireD::dDirXddeltaYd]
537  +trackDerivatives[CDCTrackD::dDOCAdDirY]*wireDerivatives[CDCWireD::dDirYddeltaYd]
538  +trackDerivatives[CDCTrackD::dDOCAdDirZ]*wireDerivatives[CDCWireD::dDirZddeltaYd];
539  label[9]=1000 + (straw_offset[thisWire->ring]+(thisWire->straw-1))*4 + 4;
540 
541  milleWriter->mille(NLC, localDer, NGL, globalDer, label, resi, err);
542  }
543 
544  }
545  milleWriter->end();
546 
547  japp->RootUnLock();
548 
549  }// End loop over tracks
550 
551  return NOERROR;
552  }
553 
554  //------------------
555  // erun
556  //------------------
558  {
559  // This is called whenever the run number changes, before it is
560  // changed to give you a chance to clean up before processing
561  // events from the next run number.
562  return NOERROR;
563  }
564 
565  //------------------
566  // fini
567  //------------------
569  {
570  // Called before program exit after event processing is finished.
571  delete milleWriter;
572  return NOERROR;
573  }
574 
Definition: track.h:16
DApplication * dapp
jerror_t evnt(jana::JEventLoop *eventLoop, uint64_t eventnumber)
Called every event.
const DCDCWire * wire
Definition: DCDCTrackHit.h:37
int status
status word for pseudopoint
Definition: DFDCPseudo.h:94
vector< double > GetCDCWireDerivatives()
Definition: DCDCWire.h:100
sprintf(text,"Post KinFit Cut")
DMagneticFieldMap * GetBfield(unsigned int run_number=1)
const DFDCWire * wire
DFDCWire for this wire.
Definition: DFDCPseudo.h:92
DMatrix3x1 index
Definition: DFDCPseudo.h:30
class DFDCPseudo: definition for a reconstructed point in the FDC
Definition: DFDCPseudo.h:74
vector< double > GetFDCPseudoAlignmentDerivatives()
Definition: DFDCPseudo.h:125
JApplication * japp
int layer
1-24
Definition: DFDCWire.h:16
jerror_t erun(void)
Called everytime run number changes, provided brun has been called.
centroid_t cluster_v
Cathode cluster centroids, Useful for gain/strip pitch calibration.
Definition: DFDCPseudo.h:88
InitPlugin_t InitPlugin
vector< double > GetFDCStripGainDerivatives()
Definition: DFDCPseudo.h:160
int straw
Definition: DCDCWire.h:81
vector< DTrackFitter::pull_t > pulls
&lt;A href=&quot;index.html#legend&quot;&gt; &lt;IMG src=&quot;CORE.png&quot; width=&quot;100&quot;&gt; &lt;/A&gt;
Class to write C binary file.
void cdc_hit(Int_t &, Int_t &, Int_t &, Int_t[], Int_t, Int_t, Int_t, Int_t, Int_t, Int_t)
Definition: fa125fns.h:55
float chisq
Chi-squared for the track (not chisq/dof!)
const vector< double > GetFDCStripPitchDerivatives()
Definition: DFDCPseudo.h:207
centroid_t cluster_u
Definition: DFDCPseudo.h:88
jerror_t fini(void)
Called after last event of last event source has been processed.
jerror_t brun(jana::JEventLoop *eventLoop, int32_t runnumber)
Called everytime a new run number is detected.
int ring
Definition: DCDCWire.h:80
jerror_t init(void)
Called once at program start.
int Ndof
Number of degrees of freedom in the fit.