Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
mvtxOM.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file mvtxOM.cc
1 // Updated by Xiaochun He on May 31, 2019 following Martin Purschke's
2 // suggestion correction
3 //
4 #include "mvtxOM.h"
5 
6 #include <TCanvas.h>
7 #include <TF1.h>
8 #include <TGaxis.h>
9 #include <TGraphAsymmErrors.h>
10 #include <TH1.h>
11 #include <TH2.h>
12 #include <TLatex.h>
13 #include <TLine.h>
14 #include <TStyle.h>
15 
16 #include <fstream>
17 #include <iostream>
18 #include <map>
19 
20 #define IDMVTXV1_MAXRUID 4
21 #define IDMVTXV1_MAXRUCHN 28
22 
23 int init_done = 0;
24 
25 using namespace std;
26 
27 const int NSTAVE = 4;
28 const bool chip_expected[4] = {true, true, true, true};
29 string stave_name[4] = {"E103", "C105", "C104", "A105"};
30 
32 
33 string outHitLocations = "/home/maps/meeg/felix/daq/felix_rcdaq/online_monitoring/hitLocations/locations.txt";
34 ofstream write_outHitLocations(outHitLocations.c_str());
35 
36 map<pair<int,int>,pair<int,int>> chipmap = {
37  {{1,1}, {0,0}},
38  {{1,2}, {0,1}},
39  {{1,3}, {0,2}},
40  {{1,4}, {0,3}},
41  {{1,5}, {0,4}},
42  {{1,6}, {0,5}},
43  {{1,7}, {0,6}},
44  {{1,8}, {0,7}},
45  {{1,9}, {3,8}},
46  {{1,10}, {3,7}},
47  {{1,11}, {3,6}},
48  {{1,12}, {3,5}},
49  {{1,13}, {3,4}},
50  {{1,14}, {3,3}},
51  {{1,15}, {3,2}},
52  {{1,16}, {3,1}},
53  {{1,17}, {3,0}},
54  {{1,18}, {2,8}},
55  {{1,19}, {2,7}},
56  {{1,20}, {2,6}},
57  {{1,21}, {2,5}},
58  {{1,22}, {2,4}},
59  {{1,23}, {2,3}},
60  {{1,24}, {2,2}},
61  {{1,25}, {2,1}},
62  {{1,26}, {2,0}},
63  {{1,27}, {0,8}},
64  {{2,1}, {1,2}},
65  {{2,2}, {1,1}},
66  {{2,3}, {1,0}},
67  {{2,4}, {1,3}},
68  {{2,5}, {1,4}},
69  {{2,6}, {1,5}},
70  {{2,7}, {1,6}},
71  {{2,8}, {1,7}},
72  {{2,9}, {1,8}}
73 }; //<ruid, ruchn> to <stave, chipID>
74 
76 int mvtx_verbose = 0;
77 int mvtx_refresh = -1;
78 int max_npixels = 512*1024; // cut to suppress plotting events with lots of hits
79 int mvtx_run = 0;
80 const bool flip_yaxis = false;
81 
82 //-- histograms filled in event loop
83 TH1F* hnevnt; // number of events
84 TH1F* hchip; // hits vs chip
85 TH1F* hwarn;
86 TH1F* herr;
87 TH1F* hnhit_chip[NSTAVE]; // nhit distribution for each chip
88 TH2F* h2d_chip[NSTAVE]; // 2D pixel hits for each chip
91 
92 //-- variables for online monitoring plotting
93 TCanvas* com; // canvas
94 TPad* phitrate;
95 TPad* pnhit;
96 TPad* p2d[NSTAVE];
97 TPad* pmean;
98 TPad* pdiffcol;
99 TPad* pdiffrow;
100 TH1F* haxis_chip; // axis for hits/event vs chip
101 TH1F* haxis_nhit; // axis for hits/event distribution for each chip
102 TH1F* haxis_2d; // axis for 2D pixel hits
103 TH1F* haxis_diff; // axis for diff pixel index (row or col)
104 TLatex* lnevents;
105 TLatex* ldiffcol[NSTAVE];
106 TLatex* ldiffrow[NSTAVE];
107 TLatex* lnhitmean[NSTAVE];
108 TGaxis* reversedaxis;
109 
110 //-- analysis histograms
116 TGraphAsymmErrors* gmpos_chip[NSTAVE];
120 
121 // The following line caused some problem of running ROOT6 in macro
122 //TF1* fg = new TF1("fg", "gaus", 0, 1024);
123 
124 TF1* fg;
125 
126 //-- some constants for different chips
127 int chipColor[] = {kBlue, kRed, kGreen+2, kMagenta+2};
128 int chipMarker[] = {kFullCircle, kFullSquare, kFullDiamond, kFullCross};
129 
130 // Show fit to beam center
131 TCanvas* cBeamCenter = nullptr;
132 const bool show_beam_fit = true;
133 
134 unsigned short decode_row(int hit)
135 {
136  return hit >> 16;
137 }
138 
139 unsigned short decode_col(int hit)
140 {
141  return hit & 0xffff;
142 }
143 
144 
145 //============================================================//
146 
147 void set_verbose(int v)
148 {
149  mvtx_verbose = v;
150 }
151 
152 //============================================================//
153 
154 void set_refresh(int r)
155 {
156  mvtx_refresh = r;
157 }
158 
159 //============================================================//
160 
161 int pinit()
162 {
163 
164  // Added by Xiaochun He following Martin's recommendation
165  fg = 0;
166 
167  if (init_done) return 1;
168  init_done = 1;
169 
170  // reset event counter
171  mvtx_evnts = 0;
172 
173  hnevnt = new TH1F("hnevent", "", 1, -0.5, 0.5);
174 
175  // hits vs chip index
176  hchip = new TH1F("hchip", ";chip index;N pixels", 4, -0.5, 3.5);
177  hchip->SetLineWidth(2);
178  hchip->SetLineColor(kBlue);
179 
180  // warnings vs chip index
181  hwarn = new TH1F("hwarn", ";chip index; warnings", 4, -0.5, 3.5);
182  hwarn->SetLineWidth(2);
183  hwarn->SetLineColor(kYellow+2);
184 
185  // errors vs chip index
186  herr = new TH1F("herr", ";chip index; errors", 4, -0.5, 3.5);
187  herr->SetLineWidth(2);
188  herr->SetLineColor(kRed);
189 
190  char name[500];
191  for (int i = 0; i < NSTAVE; i++)
192  {
193  // total hits/event distribution in each chip
194  sprintf(name, "hnhit_chip_%i", i);
195  hnhit_chip[i] = new TH1F(name, ";N pixels / event; events", 201, -0.5, 200.5);
196  hnhit_chip[i]->SetLineWidth(2);
197  hnhit_chip[i]->SetLineColor(chipColor[i]);
198 
199  // time distribution
200  sprintf(name, "hhittime_chip_%i", i);
201  hhittime_chip[i] = new TH1F(name, ";N pixels / 100 events; event number", 100, 0, 10000);
202  hhittime_chip[i]->SetLineWidth(2);
203  hhittime_chip[i]->SetLineColor(chipColor[i]);
204 
205  // 2D pixel hit dsitribution for each chip
206  sprintf(name, "h2d_chip_%i", i);
207  h2d_chip[i] = new TH2F(name, ";col;row;hits", 9216, -0.5, 9215.5, 512, -0.5, 511.5);
208 
209  // mean pixel index for each chip (aggregated over events)
210  sprintf(name, "gmpos_chip_%i", i);
211  gmpos_chip[i] = new TGraphAsymmErrors();
212  gmpos_chip[i]->SetName(name);
213  gmpos_chip[i]->SetMarkerStyle(chipMarker[i]);
214  gmpos_chip[i]->SetMarkerColor(chipColor[i]);
215  gmpos_chip[i]->SetLineColor(chipColor[i]);
216 
217  // difference in mean pixel index per event (from chip 0)
218  sprintf(name, "hdiffrow_chip_%i", i);
219  hdiffrow_chip[i] = new TH1F(name, ";row index diff", 1023, -511.5, 511.5);
220  hdiffrow_chip[i]->SetLineWidth(2);
221  hdiffrow_chip[i]->SetLineColor(chipColor[i]);
222 
223  sprintf(name, "hdiffcol_chip_%i", i);
224  hdiffcol_chip[i] = new TH1F(name, ";col index diff", 1023, -511.5, 511.5);
225  hdiffcol_chip[i]->SetLineWidth(2);
226  hdiffcol_chip[i]->SetLineColor(chipColor[i]);
227  }
228 
229 
230  // --------------------------------------------
231  // -- Setup canvas for Online Monitoring plots
232  // --------------------------------------------
233 
234  // -- setup axis
235  haxis_chip = new TH1F("haxis_chip",
236  ";chip index; N pixels / event",
237  NSTAVE, -0.5, NSTAVE-0.5);
238  haxis_chip->SetMinimum(0);
239  haxis_chip->SetMaximum(500);
240 
241  haxis_nhit = new TH1F("haxis_nhit",
242  ";N pixels / event; per event",
243  51, -0.5, 50.5);
244  haxis_nhit->SetMinimum(0);
245  haxis_nhit->SetMaximum(1.0);
246 
247  haxis_2d = new TH1F("haxis_2d",
248  ";row;col",
249  512, -0.5, 511.5);
250  haxis_2d->SetMinimum(-0.5);
251  haxis_2d->SetMaximum(511.5);
252  if (flip_yaxis)
253  {
254  haxis_2d->GetYaxis()->SetLabelOffset(999);
255  haxis_2d->GetYaxis()->SetTickLength(0);
256  }
257 
258  haxis_diff = new TH1F("haxis_diff",
259  ";mean idx diff",
260  1023, -511.5, 511.5);
261  haxis_diff->SetMinimum(0);
262  haxis_diff->SetMaximum(1);
263 
264  return 0;
265 
266 }
267 
268 //============================================================//
269 
271 {
272 
273  // Added by Xiaochu He following Martin's recommedation
274  if (!fg) fg = new TF1("fg", "gaus", 0, 1024);
275 
276  if ( e->getEvtType() == BEGRUNEVENT )
277  {
278  mvtx_run = e->getRunNumber();
279  mvtx_evnts = 0;
280  reset_histos();
281  OM();
282  return 0;
283  }
284  if ( e->getEvtType() != DATAEVENT )
285  return 0;
286 
287 
288  Packet *p = e->getPacket(2000);
289  if (p)
290  {
291 
292  bool evnt_err = false;
293 
294 
295  /*
296  int bad_ruchns = p->iValue(0, "BAD_RUCHNS");
297  if ( bad_ruchns > 0 )
298  {
299  if ( mvtx_verbose > 0 )
300  cout << "WARNING!! Event: " << mvtx_evnts << " Invalid RU channel IDs (really bad data)!"
301  << " BAD_RUCHNS:" << bad_ruchns << endl;
302  }
303 
304  int bad_chipids = p->iValue(0, "BAD_CHIPIDS");
305  if ( bad_chipids > 0 )
306  {
307  if ( mvtx_verbose > 0 )
308  cout << "WARNING!! Event: " << mvtx_evnts << " Invalid chip IDs (bad data)!"
309  << " BAD_CHIPIDS:" << bad_chipids << endl;
310  }
311 
312  int chipmax = p->iValue(0, "HIGHEST_CHIP") + 1;
313  if ( chipmax > NSTAVE )
314  {
315  if ( mvtx_verbose > 1 )
316  cout << "WARNING!! Event: " << mvtx_evnts << " More chips than expected!"
317  << " NSTAVE:" << NSTAVE << " HIGHEST_CHIP:" << chipmax << endl;
318 
319  chipmax = NSTAVE;
320  }
321  if ( mvtx_verbose > 2 )
322  {
323  cout << "Event:" << mvtx_evnts << " chipmax:" << chipmax << endl;
324  }
325  float mrow_chip0 = -1;
326  float mcol_chip0 = -1;
327 
328  int excess_data_bytes = p->iValue(0, "EXCESS_DATA_BYTES");
329  if ( excess_data_bytes>0 )
330  {
331  if ( mvtx_verbose > 1 )
332  cout << "WARNING!! Event: " << mvtx_evnts << " Data found past chip trailer"
333  << " EXCESS_DATA_BYTES: " << excess_data_bytes << endl;
334  }
335  bool evnt_err = false;
336  for ( int ichip = 0; ichip < NSTAVE; ichip++)
337  {
338  int header_found = p->iValue(ichip, "HEADER_FOUND");
339  int trailer_found = p->iValue(ichip, "TRAILER_FOUND");
340  int bunchcounter = p->iValue(ichip, "BUNCHCOUNTER");
341  int unexpected_bytes = p->iValue(ichip, "UNEXPECTED_BYTES");
342  int readout_flags = p->iValue(ichip, "READOUT_FLAGS");
343  //cout << "HEADER_FOUND: " << header_found << " TRAILER_FOUND: " << trailer_found << " BUNCHCOUNTER: " << bunchcounter << endl;
344 
345  bool has_warning = false;
346  bool has_error = false;
347  if (chip_expected[ichip])
348  {
349  if (header_found==0 || trailer_found==0)
350  {
351  if ( mvtx_verbose > 1 )
352  cout << "WARNING!! Event: " << mvtx_evnts << " Missing chip " << ichip << " HEADER_FOUND: " << header_found << " TRAILER_FOUND: " << trailer_found << " BUNCHCOUNTER: " << bunchcounter << endl;
353  has_warning = true;
354  }
355  if ( (header_found + trailer_found) == 1 )
356  {
357  if ( mvtx_verbose > 1 )
358  cout << "ERROR!! Event: " << mvtx_evnts << " header and trailer have different states for chip " << ichip << " HEADER_FOUND: " << header_found << " TRAILER_FOUND: " << trailer_found << " BUNCHCOUNTER: " << bunchcounter << endl;
359  has_error = true;
360  evnt_err = true;
361  }
362  }
363  else
364  {
365  if (header_found!=0 || trailer_found!=0)
366  {
367  if ( mvtx_verbose > 1 )
368  cout << "WARNING!! Event: " << mvtx_evnts << " Unexpected chip " << ichip << " HEADER_FOUND: " << header_found << " TRAILER_FOUND: " << trailer_found << " BUNCHCOUNTER: " << bunchcounter << endl;
369  has_warning = true;
370  }
371  }
372  if (unexpected_bytes!=0)
373  {
374  if ( mvtx_verbose > 0 )
375  cout << "WARNING!! Event: " << mvtx_evnts << " chip " << ichip << " UNEXPECTED_BYTES: " << unexpected_bytes << endl;
376  has_warning = true;
377  }
378  if (readout_flags > 0)
379  {
380  if ( mvtx_verbose > 1 )
381  cout << "WARNING!! Event: " << mvtx_evnts << " chip " << ichip << " READOUT_FLAGS: " << hex << readout_flags << dec << endl;
382  has_warning = true;
383  }
384 
385  if ( has_warning )
386  hwarn->Fill(ichip);
387  if ( has_error )
388  herr->Fill(ichip);
389  }
390  */
391 
392  int npixels[NSTAVE] = {0};
393  double mrow[NSTAVE] = {0};
394  double mcol[NSTAVE] = {0};
395  double mrow_refstave = -1;
396  double mcol_refstave = -1;
397 
398  if ( !evnt_err ) {
399  for (int ruid=0; ruid<IDMVTXV1_MAXRUID+1; ruid++)
400  {
401  if (p->iValue(ruid)!=-1)
402  {
403  for ( int ruchn = 0; ruchn < IDMVTXV1_MAXRUCHN+1; ruchn++)
404  {
405  if (p->iValue(ruid,ruchn)>0)
406  {
407  for (int i=0;i<p->iValue(ruid,ruchn);i++)
408  {
409  int hit = p->iValue(ruid,ruchn,i);
410  int irow = decode_row(hit);
411  int icol = decode_col(hit);
412  //cout << "(ruid " << ruid << ", ruchn " << ruchn << ") ";
413  //cout << "(row " << irow << ", col " << icol << ") ";
414  if (chipmap.count({ruid,ruchn}) != 1) {
415  cout << "invalid: (ruid " << ruid << ", ruchn " << ruchn << ") " << endl;
416  } else {
417  pair<int, int> chiplocation = chipmap[{ruid,ruchn}];
418  int istave = chiplocation.first;
419  int ichip = chiplocation.second;
420  //cout << "(stave " << istave << ", chip " << ichip << ") ";
421  npixels[istave]++;
422  mrow[istave]+=irow;
423  mcol[istave]+=icol+1024*ichip;
424  if (flip_yaxis)
425  h2d_chip[istave]->Fill(icol+1024*ichip,irow);
426  else
427  {
428  h2d_chip[istave]->Fill(icol+1024*ichip,511-irow);
429  if(ichip == 4){
430  write_outHitLocations<<e->getEvtSequence()<<", "<< istave<<", "<< icol<<", "<<irow<<endl;
431  //printf("There is a hit on stave %i at x = %i, y = %i for event %i\n", istave, icol, irow, e->getEvtSequence());
432  }
433  }
434  }
435  }
436  //cout << endl;
437  }
438  }
439  }
440  }
441  for (int istave=0;istave<NSTAVE;istave++) {
442  if (npixels[istave] > 0 && npixels[istave] < max_npixels)
443  {
444  mrow[istave] /= (float)npixels[istave];
445  mcol[istave] /= (float)npixels[istave];
446 
447  hchip->Fill(istave, npixels[istave]);
448  hhittime_chip[istave]->Fill(mvtx_evnts,npixels[istave]);
449 
450  if ( mvtx_verbose > 2 )
451  {
452  cout << " chip:" << istave << " npixels:" << npixels[istave]
453  << " mean:(" << mcol[istave] << ", " << mrow[istave] << ")" << endl;
454  }
455  }
456  }
457  if ( npixels[NSTAVE-1] > 0 && npixels[NSTAVE-1] < max_npixels )
458  {
459  mrow_refstave = mrow[NSTAVE-1];
460  mcol_refstave = mcol[NSTAVE-1];
461  for (int istave=0;istave<NSTAVE;istave++) {
462  if (npixels[istave] > 0 && npixels[istave] < max_npixels)
463  {
464  if ( mrow_refstave >= 0 && mcol_refstave >= 0 )
465  {
466  hdiffrow_chip[istave]->Fill(mrow[istave] - mrow_refstave);
467  hdiffcol_chip[istave]->Fill(mcol[istave] - mcol_refstave);
468  }
469 
470  }
471  }
472  }
473  }
474 
475  /*
476 
477  // fill regardless of the number of hits
478  hnhit_chip[ichip]->Fill(npixels);
479 
480  // only fill for nonzero hits
481  if (npixels > 0 && npixels < max_npixels)
482  {
483  mrow /= (float)npixels;
484  mcol /= (float)npixels;
485 
486  if ( ichip == NSTAVE - 1 )
487  {
488  mrow_chip0 = mrow;
489  mcol_chip0 = mcol;
490  }
491 
492  if ( mrow_chip0 >= 0 && mcol_chip0 >= 0 )
493  {
494  hdiffrow_chip[ichip]->Fill(mrow - mrow_chip0);
495  hdiffcol_chip[ichip]->Fill(mcol - mcol_chip0);
496  }
497 
498  hchip->Fill(ichip, npixels);
499  hhittime_chip[ichip]->Fill(mvtx_evnts,npixels);
500 
501  if ( mvtx_verbose > 2 )
502  {
503  cout << " chip:" << ichip << " npixels:" << npixels
504  << " mean:(" << mcol << ", " << mrow << ")" << endl;
505  }
506  }
507  } // ichip
508  */
509  hnevnt->Fill(0);
510  delete p;
511 }
512 
513 if ( mvtx_refresh > 0 && mvtx_evnts%mvtx_refresh == 0 ) OM();
514 
515 if ( mvtx_evnts%100000 == 0 )
516  cout << " processing event " << mvtx_evnts << endl;
517 
518  mvtx_evnts++;
519  return 0;
520  }
521 
522 //============================================================//
523 
524 int process_histos(float thresh)
525 {
526  int sum = 0;
527  int row,col;
528  for (int i = 0; i < NSTAVE; i++)
529  {
530  cout << endl;
531  cout << "-- stave " << i << endl;
532  int tot = 0;
533  for (int ix = 1; ix <= h2d_chip[i]->GetNbinsX(); ix++)
534  {
535  for (int iy = 1; iy <= h2d_chip[i]->GetNbinsY(); iy++)
536  {
537  if ( h2d_chip[i]->GetBinContent(ix, iy) > thresh )
538  {
539  tot++;
540  if (flip_yaxis)
541  {
542  row = 1023-(h2d_chip[i]->GetYaxis()->GetBinCenter(iy));
543  } else
544  {
545  row = h2d_chip[i]->GetYaxis()->GetBinCenter(iy);
546  }
547  col = h2d_chip[i]->GetXaxis()->GetBinCenter(ix);
548  cout << " row:" << row
549  << " col:" << col
550  << endl;
551  }
552  }
553  }
554  cout << "Total: " << tot << endl;
555  sum += tot;
556  }
557 
558  return sum;
559 }
560 
561 //============================================================//
562 
563 int mask_pixels(float thresh)
564 {
565  bool flip_yaxis = true; //Only uncomment if we need to flip the yAxis
566 
567  string file_name_tb1 = "/home/maps/git/RUv1_Test_sync2018-08/software/py/masklist_testbench1.txt";
568  string file_name_tb2 = "/home/maps/git/RUv1_Test_sync2018-08/software/py/masklist_testbench2.txt";
569  ofstream write_mask_file_tb1(file_name_tb1.c_str());
570  ofstream write_mask_file_tb2(file_name_tb2.c_str());
571 
572  write_mask_file_tb1<<"#Connector, ChipID, Col, Row"<<endl;
573  write_mask_file_tb2<<"#Connector, ChipID, Col, Row"<<endl;
574 
575  int sum = 0;
576  int chipid, row, col, tot, chip_tot, prev_tot;
577  Int_t stave_map[4] = {0, 4, 1, 2};
578  for (int i = 0; i < NSTAVE; i++)
579  {
580  tot = 0;
581  prev_tot = 0;
582  for (int ix = 1; ix <= h2d_chip[i]->GetNbinsX(); ix++)
583  {
584  for (int iy = 1; iy <= h2d_chip[i]->GetNbinsY(); iy++)
585  {
586  col = h2d_chip[i]->GetXaxis()->GetBinCenter(ix);
587  chipid = col/1024;
588  col = col%1024;
589 
590  if ( h2d_chip[i]->GetBinContent(ix, iy) > thresh )
591  {
592  tot++;
593  if (flip_yaxis)
594  {
595  row = 511-(h2d_chip[i]->GetYaxis()->GetBinCenter(iy));
596  }
597  else
598  {
599  row = h2d_chip[i]->GetYaxis()->GetBinCenter(iy);
600  }
601  if (i != 1) { write_mask_file_tb1 << stave_map[i] << "," << chipid << "," << col << "," << row << endl; }
602  else { write_mask_file_tb2 << stave_map[i] << "," << chipid << "," << col << "," << row << endl; }
603  }
604  }
605  chip_tot = tot - prev_tot;
606  if (col == 1023) {prev_tot = tot; printf("Total pixels masked on stave %i, chip %i: %i\n", i, chipid, chip_tot);}
607  }
608  sum += tot;
609  }
610  printf("Total pixels masked: %i\n", sum);
611  write_mask_file_tb1.close();
612  write_mask_file_tb2.close();
613  printf("New pixel mask has been created\n");
614 
615  return sum;
616 }
617 
618 //============================================================//
619 
620 int analysis()
621 {
622  //TH1D* hhitrate_chip;
623  //TH1D* hnhit_chip_norm[NSTAVE];
624  //TH2D* h2d_chip_norm[NSTAVE];
625  //TH1F* hdiffrow_chip_norm[NSTAVE];
626  //TH1F* hdiffcol_chip_norm[NSTAVE];
627  char name[500];
628 
629  //-- Check if we've already initialized the histograms
630  // if not, initialize
631  if ( !hhitrate_chip )
632  {
633  hhitrate_chip = (TH1F*) hchip->Clone("hhitrate_chip");
634  hwarnrate_chip = (TH1F*) hwarn->Clone("hwarnrate_chip");
635  herrrate_chip = (TH1F*) herr->Clone("herrrate_chip");
636 
637  for (int ichip = 0; ichip < NSTAVE; ichip++)
638  {
639  sprintf(name, "hnhit_chip_norm_%i", ichip);
640  hnhit_chip_norm[ichip] = (TH1F*) hnhit_chip[ichip]->Clone(name);
641 
642  sprintf(name, "h2d_chip_norm_%i", ichip);
643  h2d_chip_norm[ichip] = (TH2F*) h2d_chip[ichip]->Clone(name);
644 
645  sprintf(name, "hdiffrow_chip_norm_%i", ichip);
646  hdiffrow_chip_norm[ichip] = (TH1F*) hdiffrow_chip[ichip]->Clone(name);
647 
648  sprintf(name, "hdiffcol_chip_norm_%i", ichip);
649  hdiffcol_chip_norm[ichip] = (TH1F*) hdiffcol_chip[ichip]->Clone(name);
650  } // ichip
651  }
652 
653  //-- current number of events
654  double nevents = hnevnt->Integral();
655 
656  if ( nevents <= 0 )
657  return 1;
658 
659 
660  //-- hitrate vs chip index
661  hhitrate_chip->Reset();
662  hwarnrate_chip->Reset();
663  herrrate_chip->Reset();
664  for (int ib = 1; ib <= hchip->GetNbinsX(); ib++)
665  {
666  // hits
667  float bc = hchip->GetBinContent(ib);
668  hhitrate_chip->SetBinContent(ib, bc / nevents);
669  hhitrate_chip->SetBinError(ib, sqrt(bc) / nevents);
670 
671  // warnings
672  bc = hwarn->GetBinContent(ib);
673  hwarnrate_chip->SetBinContent(ib, bc / nevents);
674  hwarnrate_chip->SetBinError(ib, sqrt(bc) / nevents);
675 
676  // errors
677  bc = herr->GetBinContent(ib);
678  herrrate_chip->SetBinContent(ib, bc / nevents);
679  herrrate_chip->SetBinError(ib, sqrt(bc) / nevents);
680  } // ib
681 
682  if (show_beam_fit) { // create canvas for beamcenter plots
683  cBeamCenter = new TCanvas("cBeamCenter", "cBeamCenter", 1350, 900);
684  cBeamCenter->Divide(3,4);
685  }
686 
687  //-- per chip histograms
688  for (int ichip = 0; ichip < NSTAVE; ichip++)
689  {
690  //-- nhit distribution for each chip
691  hnhit_chip_norm[ichip]->Reset();
692  for (int ib = 1; ib <= hnhit_chip[ichip]->GetNbinsX(); ib++)
693  {
694  float bc = hnhit_chip[ichip]->GetBinContent(ib);
695  hnhit_chip_norm[ichip]->SetBinContent(ib, bc / nevents);
696  hnhit_chip_norm[ichip]->SetBinError(ib, sqrt(bc) / nevents);
697  } // ib
698 
699  //-- 2D distribution for each chip
700  h2d_chip_norm[ichip]->Reset();
701  for (int ix = 1; ix <= h2d_chip[ichip]->GetNbinsX(); ix++)
702  for (int iy = 1; iy <= h2d_chip[ichip]->GetNbinsY(); iy++)
703  {
704  float bc = h2d_chip[ichip]->GetBinContent(ix, iy);
705  h2d_chip_norm[ichip]->SetBinContent(ix, iy, bc / nevents);
706  }
707 
708  //-- mean pixel index for each chip
709  double m, s;
710  TH1D* h;
711 
712  if ( cBeamCenter )
713  {
714  cBeamCenter->cd(3*ichip+1);
715  gPad->SetLogz();
716  TH2D* h2d_beam = (TH2D*)h2d_chip_norm[ichip]->Clone("h2d_beam");
717  h2d_beam->GetXaxis()->SetTitle("Stave_Cols");
718  h2d_beam->GetYaxis()->SetTitle("Stave_Rows");
719  h2d_beam->Draw("colz");
720  TLatex lt;
721  lt.SetNDC();
722  lt.SetTextAlign(22);
723  lt.SetTextSize(1.5*lt.GetTextSize());
724  lt.SetTextColor(1);
725  lt.DrawText(0.5,0.96, stave_name[ichip].c_str());
726 
727  switch ( ichip )
728  {
729  case 0 :
730  dead_chip_forward[0]->Draw();
731  dead_chip_backward[0]->Draw();
732  break;
733 
734  case 2 :
735  dead_chip_forward[0]->Draw(); dead_chip_backward[0]->Draw();
736  dead_chip_forward[3]->Draw(); dead_chip_backward[3]->Draw();
737  dead_chip_forward[6]->Draw(); dead_chip_backward[6]->Draw();
738  dead_chip_forward[7]->Draw(); dead_chip_backward[7]->Draw();
739  dead_chip_forward[8]->Draw(); dead_chip_backward[8]->Draw();
740  break;
741 
742  case 3 :
743  dead_chip_forward[6]->Draw(); dead_chip_backward[6]->Draw();
744  break;
745 
746  default :
747  break;
748  }
749  }
750 
751  // cols projection
752  h = (TH1D*) h2d_chip_norm[ichip]->ProjectionX();
753  double m_col = h->GetMean();
754  double r_col = h->GetRMS();
755  fg->SetParameters(h->GetBinContent(h->GetMaximumBin()), m_col, r_col);
756  //YCM: initial range to full stave range
757  const int lastBinX = 9 * 1024 - 1;
758  const int lastBinY = 511;
759  fg->SetRange(0, lastBinX);
760  h->Fit(fg, "RQ0N");
761  m = fg->GetParameter(1);
762  s = fg->GetParameter(2);
763  //cout << " x1 mean:" << m << " sig:" << s << endl;
764  m = m - 2.5*s > lastBinX ? lastBinX : m;
765  fg->SetRange(m - 2.5*s, m + 2.5*s);
766  fg->SetParameters(fg->GetParameter(0), m, s);
767  h->Fit(fg, "RQ0N");
768  m = fg->GetParameter(1);
769  s = fg->GetParameter(2);
770  //cout << " x2 mean:" << m << " sig:" << s << endl;
771  if ( m > lastBinX )
772  {
773  s = s - (m - lastBinX);
774  m = lastBinX;
775  }
776  if ( m < 0 )
777  {
778  s = s - m;
779  m = 0;
780  }
781  //cout << " (m:" << m << " s:" << s << ")"<< endl;
782 
783  // save
784  m_col = m;
785  r_col = s;
786  if ( cBeamCenter )
787  {
788  cBeamCenter->cd( 3 * ichip + 2 );
789  h->GetXaxis()->SetTitle("Stave_Cols");
790  h->GetYaxis()->SetTitle("nClusters");
791  h->Draw();
792  TF1* f_col = (TF1*)fg->Clone("fx");
793  f_col->Draw("same");
794  TLatex lt;
795  lt.SetNDC();
796  lt.SetTextAlign(22);
797  lt.SetTextSize(1.5*lt.GetTextSize());
798  lt.SetTextColor(2);
799  lt.DrawLatex(.3, .7, Form("#mu_{COL} = %.2f #pm %.2f ", m_col, r_col));
800  lt.DrawLatex(.3, .6, Form("offset = %.2f", m_col));
801  }
802  else
803  delete h;
804 
805  // row projection
806  h = (TH1D*) h2d_chip_norm[ichip]->ProjectionY();
807  double m_row = h->GetMean();
808  double r_row = h->GetRMS();
809  fg->SetParameters(h->GetBinContent(h->GetMaximumBin()), m_row, r_row);
810  fg->SetRange(0, 1023);
811  h->Fit(fg, "RQ0N");
812  m = fg->GetParameter(1);
813  s = fg->GetParameter(2);
814  //cout << " y2 mean:" << m << " sig:" << s << endl;
815  fg->SetRange(m - 2.5*s, m + 2.5*s);
816  fg->SetParameters(fg->GetParameter(0), m, s);
817  h->Fit(fg, "RQ0N");
818  m = fg->GetParameter(1);
819  s = fg->GetParameter(2);
820  //cout << " y2 mean:" << m << " sig:" << s << endl;
821  if ( m > lastBinY )
822  {
823  s = s - (m - lastBinY);
824  m = lastBinY;
825  }
826  if ( m < 0 )
827  {
828  s = s - m;
829  m = 0;
830  }
831  //cout << " (m:" << m << " s:" << s << ")"<< endl;
832 
833  // save
834  m_row = m;
835  r_row = s;
836  if ( cBeamCenter )
837  {
838  cBeamCenter->cd( 3 * ichip + 3 );
839  h->GetXaxis()->SetTitle("Stave_Cols");
840  h->GetYaxis()->SetTitle("nClusters");
841  h->Draw();
842  TF1* f_row = (TF1*)fg->Clone("fy");
843  f_row->Draw("same");
844  TLatex lt;
845  lt.SetNDC();
846  lt.SetTextAlign(22);
847  lt.SetTextSize(1.5*lt.GetTextSize());
848  lt.SetTextColor(2);
849  lt.DrawLatex(.3, .7, Form("#mu_{ROW} = %.2f #pm %.2f ", m_row, r_row));
850  lt.DrawLatex(.3, .6, Form("offset = %.2f", lastBinY - m_row));
851  }
852  else
853  delete h;
854 
855  gmpos_chip[ichip]->SetPoint(0, m_col, m_row);
856  gmpos_chip[ichip]->SetPointError(0, r_col, r_col, r_row, r_row);
857 
858  //-- difference in mean pixel index per event for each chip
859  hdiffrow_chip_norm[ichip]->Reset();
860  hdiffcol_chip_norm[ichip]->Reset();
861  for (int ib = 1; ib <= hdiffrow_chip[ichip]->GetNbinsX(); ib++)
862  {
863  hdiffrow_chip_norm[ichip]->SetBinContent(ib, hdiffrow_chip[ichip]->GetBinContent(ib));
864  hdiffcol_chip_norm[ichip]->SetBinContent(ib, hdiffcol_chip[ichip]->GetBinContent(ib));
865  }
866 
867  float introw = hdiffrow_chip_norm[ichip]->Integral();
868  if ( introw > 0 )
869  hdiffrow_chip_norm[ichip]->Scale(1./introw);
870 
871  float intcol = hdiffcol_chip_norm[ichip]->Integral();
872  if ( intcol > 0 )
873  hdiffcol_chip_norm[ichip]->Scale(1./intcol);
874  } // ichip
875 
876  return 0;
877 }
878 
879 //============================================================//
880 
881 int OM()
882 {
883  //-- run analysis
884  analysis();
885 
886 
887 
888  gStyle->SetOptStat(0);
889  char name[500];
890  //-- setup canvas
891  static bool initialized = false;
892 
893  if ( !initialized )
894  {
895 
896  //com3 = new TCanvas("com","MVTX online monitoring", 1600, 800);
897  com = new TCanvas("com","MVTX online monitoring", 2000, 800);
898  com->SetMargin(0, 0, 0, 0);
899 
900  // hitrate vs chip index
901  //phitrate = new TPad("phitrate", "hit rate", 0.0, 0.45, 0.29, 0.9);
902  phitrate = new TPad("phitrate", "hit rate", 0.0, 0.60, 0.29, 0.89);
903  phitrate->SetMargin(0.15, 0.05, 0.15, 0.05);
904  phitrate->SetTicks(1, 1);
905  phitrate->Draw();
906 /*
907  // nhit distributions
908  com->cd();
909  pnhit = new TPad("pnhit", "n hit", 0.29, 0.60, 0.50, 0.90);
910  pnhit->SetMargin(0.12, 0.02, 0.15, 0.05);
911  pnhit->SetTicks(1, 1);
912  pnhit->Draw();
913 
914  // mean hit positions
915  com->cd();
916  pmean = new TPad("pmean", "mean", 0.00, 0.00, 0.29, 0.45);
917  pmean->SetMargin(0.15, 0.05, 0.15, 0.10);
918  pmean->SetTicks(1, 1);
919  pmean->Draw();
920 */
921  // diff col index
922  com->cd();
923  //pdiffcol = new TPad("pdiffcol", "pdiffcol", 0.29, 0.30, 0.50, 0.60);
924  pdiffcol = new TPad("pdiffcol", "pdiffcol", 0.00, 0.30, 0.29, 0.59);
925  //pdiffcol->SetMargin(0.12, 0.02, 0.12, 0.10);
926  pdiffcol->SetMargin(0.15, 0.05, 0.15, 0.05);
927  pdiffcol->SetTicks();
928  pdiffcol->Draw();
929 
930  // diff row index
931  com->cd();
932  //pdiffrow = new TPad("pdiffrow", "pdiffrow", 0.29, 0.00, 0.50, 0.30);
933  pdiffrow = new TPad("pdiffrow", "pdiffrow", 0.00, 0.00, 0.29, 0.29);
934  //pdiffrow->SetMargin(0.12, 0.02, 0.12, 0.10);
935  pdiffrow->SetMargin(0.15, 0.05, 0.15, 0.05);
936  pdiffrow->SetTicks();
937  pdiffrow->Draw();
938 
939  // 2D hit distributions
940  for (int i = 0; i < NSTAVE; i++)
941  {
942  com->cd();
943  sprintf(name, "p2d_%i", i);
944  //p2d[i] = new TPad(name, name,
945  // 0.50, 0.9-0.22*i,
946  // 1.00, 0.9-0.22*(i+1) );
947  p2d[i] = new TPad(name, name,
948  0.30, 0.9-0.22*i,
949  1.00, 0.9-0.22*(i+1) );
950  //p2d[i]->SetMargin(0.15, 0.15, 0.15, 0.15);
951  p2d[i]->SetMargin(0.01, 0.02, 0.25, 0.25);
952  p2d[i]->SetTicks(1, 1);
953  p2d[i]->Draw();
954  } // i
955 
956  //-- info
957  com->cd();
958  sprintf(name, "Number of Events: %.0f", hnevnt->Integral());
959  lnevents = new TLatex(0.5, 0.95, name);
960  lnevents->SetNDC();
961  lnevents->SetTextAlign(22);
962  lnevents->Draw();
963 
964  TLatex lt;
965  lt.SetNDC();
966  lt.SetTextAlign(22);
967 
968  //-- hitrate
969  phitrate->cd();
970  //gPad->SetLogy();
971  haxis_chip->GetYaxis()->SetRangeUser(0,5);
972  haxis_chip->Draw();
973  hhitrate_chip->Draw("hist e same");
974  hwarnrate_chip->Draw("hist e same");
975  herrrate_chip->Draw("hist e same");
976 
977  //-- nhit distributions
978  //pnhit->cd();
979  //gPad->SetLogy();
980  haxis_nhit->GetYaxis()->SetRangeUser(1e-4, 1);
981  //haxis_nhit->DrawCopy();
982  for (int i = 0; i < NSTAVE; i++)
983  {
984  hnhit_chip_norm[i]->Draw("hist e same");
985 
986  sprintf(name, "chip %i = %.2f", i, hnhit_chip_norm[i]->GetMean());
987  lnhitmean[i] = new TLatex(0.8, 0.85 - 0.05*i, name);
988  lnhitmean[i]->SetNDC();
989  lnhitmean[i]->SetTextColor(chipColor[i]);
990  lnhitmean[i]->Draw("same");
991  }
992 
993  //-- mean positions
994  //pmean->cd();
995  //haxis_2d->DrawCopy();
996  //gPad->Update();
997 /* if (flip_yaxis)
998  {
999  reversedaxis = new TGaxis(gPad->GetUxmin(),
1000  gPad->GetUymax(),
1001  gPad->GetUxmin()-0.001,
1002  gPad->GetUymin(),
1003  -0.5,
1004  1023.5,
1005  510, "+");
1006  reversedaxis->SetLabelOffset(-0.03);
1007  reversedaxis->Draw();
1008  }
1009  for (int i = 0; i < NSTAVE; i++)
1010  gmpos_chip[i]->Draw("P");
1011  lt.DrawLatex(0.5, 0.96, "Beam position & profile");
1012 */
1013  //-- diff col
1014  pdiffcol->cd();
1015  gPad->SetLogy();
1016  haxis_diff->GetYaxis()->SetRangeUser(1e-4, 1);
1017  //haxis_diff->GetXaxis()->SetRangeUser(-25, 25);
1018  haxis_diff->Draw();
1019  for (int i = 0; i < NSTAVE - 1; i++)
1020  {
1021  hdiffcol_chip_norm[i]->Draw("hist same");
1022 
1023  sprintf(name, "chip %i = %+.1f", i, hdiffcol_chip_norm[i]->GetMean());
1024  ldiffcol[i] = new TLatex(0.8, 0.90 - 0.05*i, name);
1025  ldiffcol[i]->SetNDC();
1026  ldiffcol[i]->SetTextColor(chipColor[i]);
1027  ldiffcol[i]->Draw("same");
1028 
1029  }
1030  lt.DrawLatex(0.5, 0.96, "Col");
1031 
1032  //-- diff row
1033  pdiffrow->cd();
1034  gPad->SetLogy();
1035  haxis_diff->GetYaxis()->SetRangeUser(1e-4, 1);
1036  //haxis_diff->GetXaxis()->SetRangeUser(-25, 25);
1037  haxis_diff->Draw();
1038  for (int i = 0; i < NSTAVE - 1; i++)
1039  {
1040  hdiffrow_chip_norm[i]->Draw("hist same");
1041 
1042  sprintf(name, "chip %i = %+.1f", i, hdiffrow_chip_norm[i]->GetMean());
1043  ldiffrow[i] = new TLatex(0.8, 0.90 - 0.05*i, name);
1044  ldiffrow[i]->SetNDC();
1045  ldiffrow[i]->SetTextColor(chipColor[i]);
1046  ldiffrow[i]->Draw("same");
1047 
1048  }
1049  lt.DrawLatex(0.5, 0.96, "Row");
1050 
1051  for (int i = 1; i < 10; ++i)
1052  {
1053  chip_edges.push_back(new TLine((1024*i)-1, 0, (1024*i)-1, 511));
1054  dead_chip_forward.push_back(new TLine((1024*(i-1))-1, 0, (1024*i)-1, 511));
1055  dead_chip_backward.push_back(new TLine((1024*(i-1))-1, 511, (1024*i)-1, 0));
1056  }
1057  //-- 2D distributions
1058  for (int i = 0; i < NSTAVE; i++)
1059  {
1060  p2d[i]->cd();
1061  gPad->SetLogz();
1062  //haxis_2d->DrawCopy();
1063  gPad->Update();
1064  if (flip_yaxis)
1065  {
1066  reversedaxis = new TGaxis(gPad->GetUxmin(),
1067  gPad->GetUymax(),
1068  gPad->GetUxmin()-0.001,
1069  gPad->GetUymin(),
1070  -0.5,
1071  511.5,
1072  510, "+");
1073  reversedaxis->SetLabelOffset(-0.03);
1074  //reversedaxis->Draw();
1075  }
1076  h2d_chip_norm[i]->Draw("colz same");
1077  h2d_chip_norm[i]->GetZaxis()->SetRangeUser(1e-6,1);
1078  for (int j = 0; j < 8; ++j){chip_edges[j]->Draw();}
1079  sprintf(name, "Stave %s", stave_name[i].c_str());
1080  lt.SetTextColor(chipColor[i]);
1081  lt.DrawLatex(0.5, 0.96, name);
1082  } // i
1083 
1084  p2d[0]->cd();
1085  dead_chip_forward[0]->Draw(); dead_chip_backward[0]->Draw();
1086  p2d[2]->cd();
1087  dead_chip_forward[0]->Draw(); dead_chip_backward[0]->Draw();
1088  dead_chip_forward[3]->Draw(); dead_chip_backward[3]->Draw();
1089  dead_chip_forward[6]->Draw(); dead_chip_backward[6]->Draw();
1090  dead_chip_forward[7]->Draw(); dead_chip_backward[7]->Draw();
1091  dead_chip_forward[8]->Draw(); dead_chip_backward[8]->Draw();
1092  p2d[3]->cd();
1093  dead_chip_forward[6]->Draw(); dead_chip_backward[6]->Draw();
1094 
1095  initialized = true;
1096  }
1097 
1098  //-- Update
1099  com->cd();
1100  sprintf(name, "Run %i, Number of Events: %.0f", mvtx_run, hnevnt->Integral());
1101  lnevents->SetText(0.5, 0.95, name);
1102  for (int i = 0; i < NSTAVE; i++)
1103  {
1104  sprintf(name, "stave %s = %.2f", stave_name[i].c_str(), hnhit_chip_norm[i]->GetMean());
1105  lnhitmean[i]->SetText(0.75, 0.88 - 0.05*i, name);
1106  }
1107  for (int i = 0; i < NSTAVE - 1; i++)
1108  {
1109  if ( hdiffcol_chip_norm[i]->Integral() > 0 )
1110  {
1111  double dcol = hdiffcol_chip_norm[i]->GetBinCenter(hdiffcol_chip_norm[i]->GetMaximumBin());
1112  sprintf(name, "stave %s = %+.0f", stave_name[i].c_str(), dcol);
1113  }
1114  else
1115  sprintf(name, "");
1116  ldiffcol[i]->SetText(0.75, 0.83 - 0.05*i, name);
1117 
1118  if ( hdiffrow_chip_norm[i]->Integral() > 0 )
1119  {
1120  double drow = hdiffrow_chip_norm[i]->GetBinCenter(hdiffrow_chip_norm[i]->GetMaximumBin());
1121  sprintf(name, "stave %s = %+.0f", stave_name[i].c_str(), drow);
1122  }
1123  else
1124  sprintf(name, "");
1125  ldiffrow[i]->SetText(0.75, 0.83 - 0.05*i, name);
1126  }
1127  com->Modified();
1128  phitrate->Modified();
1129  //pnhit->Modified();
1130  //pmean->Modified();
1131  pdiffcol->Modified();
1132  pdiffrow->Modified();
1133  for (int i = 0; i < NSTAVE; i++)
1134  p2d[i]->Modified();
1135  com->Update();
1136 
1137  return 0;
1138 }
1139 
1140 //============================================================//
1141 // Print plots on canvas to a png file
1142 
1144 {
1145  char outputFileName[128];
1146  sprintf(outputFileName, "Run%0.4i.png", mvtx_run);
1147  com->Print(outputFileName);
1148  if (show_beam_fit) {
1149  cBeamCenter->Print(Form("Beam%0.4d.png", mvtx_run));
1150  }
1151  return 0;
1152 }
1153 //============================================================//
1154 
1156 {
1157  analysis();
1158  cout << "==================================================" << endl;
1159  cout << "==================================================" << endl;
1160  cout << "== Number of events: " << hnevnt->Integral() << endl;
1161  for (int i = 0; i < NSTAVE; i++)
1162  {
1163  // only print chips with hits
1164  if ( !(hchip->GetBinContent(i + 1) > 0) )
1165  continue;
1166 
1167  cout << "==================================================" << endl;
1168  cout << "==== Chip " << i << endl;
1169  cout << "== events: " << hnhit_chip[i]->Integral() << endl;
1170  cout << "== Total hits: " << hchip->GetBinContent(i+1)
1171  << " (" << hhitrate_chip->GetBinContent(i+1) << ")"
1172  << endl;
1173  cout << "== Total Warnings: " << hwarn->GetBinContent(i+1)
1174  << " (" << hwarnrate_chip->GetBinContent(i+1) << ")"
1175  << endl;
1176  cout << "== Total Errors: " << herr->GetBinContent(i+1)
1177  << " (" << herrrate_chip->GetBinContent(i+1) << ")"
1178  << endl;
1179  cout << "== <hits/event>: " << hnhit_chip_norm[i]->GetMean() << endl;
1180  cout << "== RMS(hits/event): " << hnhit_chip_norm[i]->GetRMS() << endl;
1181 
1182  double mcol, mrow;
1183  gmpos_chip[i]->GetPoint(0, mcol, mrow);
1184  cout << "== <col>: " << mcol << endl;
1185  cout << "== <row>: " << mrow << endl;
1186 
1187  if ( i > 0 )
1188  {
1189  double dcol = hdiffcol_chip_norm[i]->GetMaximumBin();
1190  dcol = hdiffcol_chip_norm[i]->GetBinCenter(dcol);
1191  cout << "== <dcol>: " << dcol
1192  << " (" << hdiffcol_chip_norm[i]->GetMean() << ")"
1193  << endl;
1194 
1195  double drow = hdiffrow_chip_norm[i]->GetMaximumBin();
1196  drow = hdiffrow_chip_norm[i]->GetBinCenter(drow);
1197  cout << "== <drow>: " << drow
1198  << " (" << hdiffrow_chip_norm[i]->GetMean() << ")"
1199  << endl;
1200  }
1201  }
1202  cout << "==================================================" << endl;
1203  cout << "==================================================" << endl;
1204 }
1205 
1206 //============================================================//
1207 
1209 {
1210  hnevnt->Reset();
1211  hchip->Reset();
1212  hwarn->Reset();
1213  herr->Reset();
1214  for (int i = 0; i < NSTAVE; i++)
1215  {
1216  hnhit_chip[i]->Reset();
1217  h2d_chip[i]->Reset();
1218  hdiffrow_chip[i]->Reset();
1219  hdiffcol_chip[i]->Reset();
1220  } // i
1221 }
1222 
1223 //============================================================//
1224 
1226 {
1227  string fname(Form("beamcenter/beamcenter_%08d.txt", mvtx_run));
1228  ofstream fout(fname);
1229 
1230  for (int istave=0; istave<NSTAVE; ++istave)
1231  {
1232  double m_col, m_row;
1233  gmpos_chip[istave]->GetPoint(0, m_col, m_row);
1234  m_row = 511 - m_row; //chip row flipped wrt histo axis
1235 
1236  fout << istave << " 0 " << fixed << setprecision(2) << m_row << " " << 0 << " " << m_col << endl;
1237  cout << istave << " 0 " << fixed << setprecision(2) << m_row << " " << 0 << " " << m_col << endl;
1238  }
1239 
1240  fout.close();
1241  return;
1242 }