Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:20:43

0001 // Updated by Xiaochun He on May 31, 2019 following Martin Purschke's
0002 // suggestion correction
0003 //
0004 #include "mvtxOM.h"
0005 
0006 #include <TCanvas.h>
0007 #include <TF1.h>
0008 #include <TGaxis.h>
0009 #include <TGraphAsymmErrors.h>
0010 #include <TH1.h>
0011 #include <TH2.h>
0012 #include <TLatex.h>
0013 #include <TLine.h>
0014 #include <TStyle.h>
0015 
0016 #include <fstream>
0017 #include <iostream>
0018 #include <map>
0019 
0020 #define IDMVTXV1_MAXRUID       4
0021 #define IDMVTXV1_MAXRUCHN      28
0022 
0023 int init_done = 0;
0024 
0025 using namespace std;
0026 
0027 const int NSTAVE = 4;
0028 const bool chip_expected[4] = {true, true, true, true};
0029 string stave_name[4] = {"E103", "C105", "C104", "A105"};
0030 
0031 vector<TLine*> chip_edges, dead_chip_forward, dead_chip_backward;
0032 
0033 string outHitLocations = "/home/maps/meeg/felix/daq/felix_rcdaq/online_monitoring/hitLocations/locations.txt";
0034 ofstream write_outHitLocations(outHitLocations.c_str());
0035 
0036 map<pair<int,int>,pair<int,int>> chipmap = {
0037   {{1,1}, {0,0}},
0038   {{1,2}, {0,1}},
0039   {{1,3}, {0,2}},
0040   {{1,4}, {0,3}},
0041   {{1,5}, {0,4}},
0042   {{1,6}, {0,5}},
0043   {{1,7}, {0,6}},
0044   {{1,8}, {0,7}},
0045   {{1,9}, {3,8}},
0046   {{1,10}, {3,7}},
0047   {{1,11}, {3,6}},
0048   {{1,12}, {3,5}},
0049   {{1,13}, {3,4}},
0050   {{1,14}, {3,3}},
0051   {{1,15}, {3,2}},
0052   {{1,16}, {3,1}},
0053   {{1,17}, {3,0}},
0054   {{1,18}, {2,8}},
0055   {{1,19}, {2,7}},
0056   {{1,20}, {2,6}},
0057   {{1,21}, {2,5}},
0058   {{1,22}, {2,4}},
0059   {{1,23}, {2,3}},
0060   {{1,24}, {2,2}},
0061   {{1,25}, {2,1}},
0062   {{1,26}, {2,0}},
0063   {{1,27}, {0,8}},
0064   {{2,1}, {1,2}},
0065   {{2,2}, {1,1}},
0066   {{2,3}, {1,0}},
0067   {{2,4}, {1,3}},
0068   {{2,5}, {1,4}},
0069   {{2,6}, {1,5}},
0070   {{2,7}, {1,6}},
0071   {{2,8}, {1,7}},
0072   {{2,9}, {1,8}}
0073 }; //<ruid, ruchn> to <stave, chipID>
0074 
0075 int mvtx_evnts;
0076 int mvtx_verbose = 0;
0077 int mvtx_refresh = -1;
0078 int max_npixels = 512*1024; // cut to suppress plotting events with lots of hits
0079 int mvtx_run = 0;
0080 const bool flip_yaxis = false;
0081 
0082 //-- histograms filled in event loop
0083 TH1F* hnevnt; // number of events
0084 TH1F* hchip;  // hits vs chip
0085 TH1F* hwarn;
0086 TH1F* herr;
0087 TH1F* hnhit_chip[NSTAVE]; // nhit distribution for each chip
0088 TH2F* h2d_chip[NSTAVE]; // 2D pixel hits for each chip
0089 TH1F* hdiffrow_chip[NSTAVE];
0090 TH1F* hdiffcol_chip[NSTAVE];
0091 
0092 //-- variables for online monitoring plotting
0093 TCanvas* com; // canvas
0094 TPad* phitrate;
0095 TPad* pnhit;
0096 TPad* p2d[NSTAVE];
0097 TPad* pmean;
0098 TPad* pdiffcol;
0099 TPad* pdiffrow;
0100 TH1F* haxis_chip; // axis for hits/event vs chip
0101 TH1F* haxis_nhit; // axis for hits/event distribution for each chip
0102 TH1F* haxis_2d;   // axis for 2D pixel hits
0103 TH1F* haxis_diff; // axis for diff pixel index (row or col)
0104 TLatex* lnevents;
0105 TLatex* ldiffcol[NSTAVE];
0106 TLatex* ldiffrow[NSTAVE];
0107 TLatex* lnhitmean[NSTAVE];
0108 TGaxis* reversedaxis;
0109 
0110 //-- analysis histograms
0111 TH1F* hhitrate_chip;
0112 TH1F* hwarnrate_chip;
0113 TH1F* herrrate_chip;
0114 TH1F* hnhit_chip_norm[NSTAVE];
0115 TH2F* h2d_chip_norm[NSTAVE];
0116 TGraphAsymmErrors* gmpos_chip[NSTAVE];
0117 TH1F* hdiffrow_chip_norm[NSTAVE];
0118 TH1F* hdiffcol_chip_norm[NSTAVE];
0119 TH1F* hhittime_chip[NSTAVE];
0120 
0121 // The following line caused some problem of running ROOT6 in macro
0122 //TF1* fg = new TF1("fg", "gaus", 0, 1024);
0123 
0124 TF1* fg;
0125 
0126 //-- some constants for different chips
0127 int chipColor[] = {kBlue, kRed, kGreen+2, kMagenta+2};
0128 int chipMarker[] = {kFullCircle, kFullSquare, kFullDiamond, kFullCross};
0129 
0130 // Show fit to beam center
0131 TCanvas* cBeamCenter = nullptr;
0132 const bool show_beam_fit = true;
0133 
0134 unsigned short decode_row(int hit)
0135 {
0136     return hit >> 16;
0137 }
0138 
0139 unsigned short decode_col(int hit)
0140 {
0141     return hit & 0xffff;
0142 }
0143 
0144 
0145 //============================================================//
0146 
0147 void set_verbose(int v)
0148 {
0149     mvtx_verbose = v;
0150 }
0151 
0152 //============================================================//
0153 
0154 void set_refresh(int r)
0155 {
0156     mvtx_refresh = r;
0157 }
0158 
0159 //============================================================//
0160 
0161 int pinit()
0162 {
0163 
0164   // Added by Xiaochun He following Martin's recommendation
0165   fg = 0;
0166 
0167   if (init_done) return 1;
0168   init_done = 1;
0169 
0170   // reset event counter
0171   mvtx_evnts = 0;
0172 
0173   hnevnt = new TH1F("hnevent", "", 1, -0.5, 0.5);
0174 
0175   // hits vs chip index
0176   hchip = new TH1F("hchip", ";chip index;N pixels", 4, -0.5, 3.5);
0177   hchip->SetLineWidth(2);
0178   hchip->SetLineColor(kBlue);
0179 
0180   // warnings vs chip index
0181   hwarn = new TH1F("hwarn", ";chip index; warnings", 4, -0.5, 3.5);
0182   hwarn->SetLineWidth(2);
0183   hwarn->SetLineColor(kYellow+2);
0184 
0185   // errors vs chip index
0186   herr = new TH1F("herr", ";chip index; errors", 4, -0.5, 3.5);
0187   herr->SetLineWidth(2);
0188   herr->SetLineColor(kRed);
0189 
0190   char name[500];
0191   for (int i = 0; i < NSTAVE; i++)
0192     {
0193       // total hits/event distribution in each chip
0194       sprintf(name, "hnhit_chip_%i", i);
0195       hnhit_chip[i] = new TH1F(name, ";N pixels / event; events", 201, -0.5, 200.5);
0196       hnhit_chip[i]->SetLineWidth(2);
0197       hnhit_chip[i]->SetLineColor(chipColor[i]);
0198 
0199       // time distribution
0200       sprintf(name, "hhittime_chip_%i", i);
0201       hhittime_chip[i] = new TH1F(name, ";N pixels / 100 events; event number", 100, 0, 10000);
0202       hhittime_chip[i]->SetLineWidth(2);
0203       hhittime_chip[i]->SetLineColor(chipColor[i]);
0204 
0205       // 2D pixel hit dsitribution for each chip
0206       sprintf(name, "h2d_chip_%i", i);
0207       h2d_chip[i] = new TH2F(name, ";col;row;hits", 9216, -0.5, 9215.5, 512, -0.5, 511.5);
0208 
0209       // mean pixel index for each chip (aggregated over events)
0210       sprintf(name, "gmpos_chip_%i", i);
0211       gmpos_chip[i] = new TGraphAsymmErrors();
0212       gmpos_chip[i]->SetName(name);
0213       gmpos_chip[i]->SetMarkerStyle(chipMarker[i]);
0214       gmpos_chip[i]->SetMarkerColor(chipColor[i]);
0215       gmpos_chip[i]->SetLineColor(chipColor[i]);
0216 
0217       // difference in mean pixel index per event (from chip 0)
0218       sprintf(name, "hdiffrow_chip_%i", i);
0219       hdiffrow_chip[i] = new TH1F(name, ";row index diff", 1023, -511.5, 511.5);
0220       hdiffrow_chip[i]->SetLineWidth(2);
0221       hdiffrow_chip[i]->SetLineColor(chipColor[i]);
0222 
0223       sprintf(name, "hdiffcol_chip_%i", i);
0224       hdiffcol_chip[i] = new TH1F(name, ";col index diff", 1023, -511.5, 511.5);
0225       hdiffcol_chip[i]->SetLineWidth(2);
0226       hdiffcol_chip[i]->SetLineColor(chipColor[i]);
0227     }
0228 
0229 
0230   // --------------------------------------------
0231   // -- Setup canvas for Online Monitoring plots
0232   // --------------------------------------------
0233 
0234   // -- setup axis
0235   haxis_chip = new TH1F("haxis_chip",
0236             ";chip index; N pixels / event",
0237             NSTAVE, -0.5, NSTAVE-0.5);
0238   haxis_chip->SetMinimum(0);
0239   haxis_chip->SetMaximum(500);
0240 
0241   haxis_nhit = new TH1F("haxis_nhit",
0242             ";N pixels / event; per event",
0243             51, -0.5, 50.5);
0244   haxis_nhit->SetMinimum(0);
0245   haxis_nhit->SetMaximum(1.0);
0246 
0247   haxis_2d = new TH1F("haxis_2d",
0248               ";row;col",
0249               512, -0.5, 511.5);
0250   haxis_2d->SetMinimum(-0.5);
0251   haxis_2d->SetMaximum(511.5);
0252   if (flip_yaxis)
0253     {
0254       haxis_2d->GetYaxis()->SetLabelOffset(999);
0255       haxis_2d->GetYaxis()->SetTickLength(0);
0256     }
0257 
0258   haxis_diff = new TH1F("haxis_diff",
0259             ";mean idx diff",
0260             1023, -511.5, 511.5);
0261   haxis_diff->SetMinimum(0);
0262   haxis_diff->SetMaximum(1);
0263 
0264   return 0;
0265 
0266 }
0267 
0268 //============================================================//
0269 
0270 int process_event (Event * e)
0271 {
0272 
0273   // Added by Xiaochu He following Martin's recommedation
0274   if (!fg) fg = new TF1("fg", "gaus", 0, 1024);
0275 
0276   if ( e->getEvtType() == BEGRUNEVENT )
0277     {
0278       mvtx_run = e->getRunNumber();
0279       mvtx_evnts = 0;
0280       reset_histos();
0281       OM();
0282       return 0;
0283     }
0284   if ( e->getEvtType() != DATAEVENT )
0285     return 0;
0286 
0287 
0288   Packet *p = e->getPacket(2000);
0289   if (p)
0290     {
0291 
0292       bool evnt_err = false;
0293 
0294 
0295       /*
0296     int bad_ruchns = p->iValue(0, "BAD_RUCHNS");
0297     if ( bad_ruchns > 0 )
0298     {
0299     if ( mvtx_verbose > 0 )
0300     cout << "WARNING!! Event: " << mvtx_evnts << " Invalid RU channel IDs (really bad data)!"
0301     << " BAD_RUCHNS:" << bad_ruchns << endl;
0302     }
0303 
0304     int bad_chipids = p->iValue(0, "BAD_CHIPIDS");
0305     if ( bad_chipids > 0 )
0306     {
0307     if ( mvtx_verbose > 0 )
0308     cout << "WARNING!! Event: " << mvtx_evnts << " Invalid chip IDs (bad data)!"
0309     << " BAD_CHIPIDS:" << bad_chipids << endl;
0310     }
0311 
0312     int chipmax = p->iValue(0, "HIGHEST_CHIP") + 1;
0313     if ( chipmax > NSTAVE )
0314     {
0315     if ( mvtx_verbose > 1 )
0316     cout << "WARNING!! Event: " << mvtx_evnts << " More chips than expected!"
0317     << " NSTAVE:" << NSTAVE << " HIGHEST_CHIP:" << chipmax << endl;
0318 
0319     chipmax = NSTAVE;
0320     }
0321     if ( mvtx_verbose > 2 )
0322     {
0323     cout << "Event:" << mvtx_evnts << " chipmax:" << chipmax << endl;
0324     }
0325     float mrow_chip0 = -1;
0326     float mcol_chip0 = -1;
0327 
0328     int excess_data_bytes = p->iValue(0, "EXCESS_DATA_BYTES");
0329     if ( excess_data_bytes>0 )
0330     {
0331     if ( mvtx_verbose > 1 )
0332     cout << "WARNING!! Event: " << mvtx_evnts << " Data found past chip trailer"
0333     << " EXCESS_DATA_BYTES: " << excess_data_bytes << endl;
0334     }
0335     bool evnt_err = false;
0336     for ( int ichip = 0; ichip < NSTAVE; ichip++)
0337     {
0338     int header_found = p->iValue(ichip, "HEADER_FOUND");
0339     int trailer_found = p->iValue(ichip, "TRAILER_FOUND");
0340     int bunchcounter = p->iValue(ichip, "BUNCHCOUNTER");
0341     int unexpected_bytes = p->iValue(ichip, "UNEXPECTED_BYTES");
0342     int readout_flags = p->iValue(ichip, "READOUT_FLAGS");
0343         //cout << "HEADER_FOUND: " << header_found << " TRAILER_FOUND: " << trailer_found << " BUNCHCOUNTER: " << bunchcounter << endl;
0344 
0345         bool has_warning = false;
0346         bool has_error = false;
0347         if (chip_expected[ichip])
0348         {
0349         if (header_found==0 || trailer_found==0)
0350         {
0351         if ( mvtx_verbose > 1 )
0352         cout << "WARNING!! Event: " << mvtx_evnts << " Missing chip " << ichip << " HEADER_FOUND: " << header_found << " TRAILER_FOUND: " << trailer_found << " BUNCHCOUNTER: " << bunchcounter << endl;
0353         has_warning = true;
0354         }
0355         if ( (header_found + trailer_found) == 1 )
0356         {
0357         if ( mvtx_verbose > 1 )
0358         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;
0359         has_error = true;
0360         evnt_err = true;
0361         }
0362         }
0363         else
0364         {
0365         if (header_found!=0 || trailer_found!=0)
0366         {
0367             if ( mvtx_verbose > 1 )
0368                 cout << "WARNING!! Event: " << mvtx_evnts << " Unexpected chip " << ichip << " HEADER_FOUND: " << header_found << " TRAILER_FOUND: " << trailer_found << " BUNCHCOUNTER: " << bunchcounter << endl;
0369             has_warning = true;
0370         }
0371     }
0372     if (unexpected_bytes!=0)
0373     {
0374         if ( mvtx_verbose > 0 )
0375             cout << "WARNING!! Event: " << mvtx_evnts << " chip " << ichip << " UNEXPECTED_BYTES: " << unexpected_bytes << endl;
0376         has_warning = true;
0377     }
0378     if (readout_flags > 0)
0379     {
0380         if ( mvtx_verbose > 1 )
0381             cout << "WARNING!! Event: " << mvtx_evnts << " chip " << ichip << " READOUT_FLAGS: " << hex << readout_flags << dec << endl;
0382         has_warning = true;
0383     }
0384 
0385     if ( has_warning )
0386         hwarn->Fill(ichip);
0387     if ( has_error )
0388         herr->Fill(ichip);
0389     }
0390     */
0391 
0392     int npixels[NSTAVE] = {0};
0393     double mrow[NSTAVE] = {0};
0394     double mcol[NSTAVE] = {0};
0395     double mrow_refstave = -1;
0396     double mcol_refstave = -1;
0397 
0398     if ( !evnt_err ) {
0399         for (int ruid=0; ruid<IDMVTXV1_MAXRUID+1; ruid++)
0400         {
0401             if (p->iValue(ruid)!=-1)
0402             {
0403                 for ( int ruchn = 0; ruchn < IDMVTXV1_MAXRUCHN+1; ruchn++)
0404                 {
0405                     if (p->iValue(ruid,ruchn)>0)
0406                     {
0407                         for (int i=0;i<p->iValue(ruid,ruchn);i++)
0408                         {
0409                             int hit = p->iValue(ruid,ruchn,i);
0410                             int irow = decode_row(hit);
0411                             int icol = decode_col(hit);
0412                             //cout << "(ruid " << ruid << ", ruchn " << ruchn << ") ";
0413                             //cout << "(row " << irow << ", col " << icol << ") ";
0414                             if (chipmap.count({ruid,ruchn}) != 1) {
0415                                 cout << "invalid: (ruid " << ruid << ", ruchn " << ruchn << ") " << endl;
0416                             } else {
0417                                 pair<int, int> chiplocation = chipmap[{ruid,ruchn}];
0418                                 int istave = chiplocation.first;
0419                                 int ichip = chiplocation.second;
0420                                 //cout << "(stave " << istave << ", chip " << ichip << ") ";
0421                                 npixels[istave]++;
0422                                 mrow[istave]+=irow;
0423                                 mcol[istave]+=icol+1024*ichip;
0424                                 if (flip_yaxis)
0425                                     h2d_chip[istave]->Fill(icol+1024*ichip,irow);
0426                                 else
0427                                 {
0428                                     h2d_chip[istave]->Fill(icol+1024*ichip,511-irow);
0429                                     if(ichip == 4){
0430                                         write_outHitLocations<<e->getEvtSequence()<<", "<< istave<<", "<< icol<<", "<<irow<<endl;
0431                                         //printf("There is a hit on stave %i at x = %i, y = %i for event %i\n", istave, icol, irow, e->getEvtSequence());
0432                                         }
0433                                  }
0434                             }
0435                         }
0436                         //cout << endl;
0437                     }
0438                 }
0439             }
0440         }
0441         for (int istave=0;istave<NSTAVE;istave++) {
0442             if (npixels[istave] > 0  && npixels[istave] < max_npixels)
0443             {
0444                 mrow[istave] /= (float)npixels[istave];
0445                 mcol[istave] /= (float)npixels[istave];
0446 
0447                 hchip->Fill(istave, npixels[istave]);
0448                 hhittime_chip[istave]->Fill(mvtx_evnts,npixels[istave]);
0449 
0450                 if ( mvtx_verbose > 2 )
0451                 {
0452                     cout << "    chip:" << istave << " npixels:" << npixels[istave]
0453                         << " mean:(" << mcol[istave] << ", " << mrow[istave] << ")" << endl;
0454                 }
0455             }
0456         }
0457         if ( npixels[NSTAVE-1] > 0 && npixels[NSTAVE-1] < max_npixels )
0458         {
0459             mrow_refstave = mrow[NSTAVE-1];
0460             mcol_refstave = mcol[NSTAVE-1];
0461             for (int istave=0;istave<NSTAVE;istave++) {
0462                 if (npixels[istave] > 0  && npixels[istave] < max_npixels)
0463                 {
0464                     if ( mrow_refstave >= 0 && mcol_refstave >= 0 )
0465                     {
0466                         hdiffrow_chip[istave]->Fill(mrow[istave] - mrow_refstave);
0467                         hdiffcol_chip[istave]->Fill(mcol[istave] - mcol_refstave);
0468                     }
0469 
0470                 }
0471             }
0472         }
0473     }
0474 
0475     /*
0476 
0477     // fill regardless of the number of hits
0478     hnhit_chip[ichip]->Fill(npixels);
0479 
0480     // only fill for nonzero hits
0481     if (npixels > 0  && npixels < max_npixels)
0482     {
0483     mrow /= (float)npixels;
0484     mcol /= (float)npixels;
0485 
0486     if ( ichip == NSTAVE - 1 )
0487     {
0488     mrow_chip0 = mrow;
0489     mcol_chip0 = mcol;
0490     }
0491 
0492     if ( mrow_chip0 >= 0 && mcol_chip0 >= 0 )
0493     {
0494     hdiffrow_chip[ichip]->Fill(mrow - mrow_chip0);
0495     hdiffcol_chip[ichip]->Fill(mcol - mcol_chip0);
0496     }
0497 
0498     hchip->Fill(ichip, npixels);
0499     hhittime_chip[ichip]->Fill(mvtx_evnts,npixels);
0500 
0501     if ( mvtx_verbose > 2 )
0502     {
0503     cout << "    chip:" << ichip << " npixels:" << npixels
0504     << " mean:(" << mcol << ", " << mrow << ")" << endl;
0505     }
0506     }
0507     } // ichip
0508     */
0509     hnevnt->Fill(0);
0510     delete p;
0511 }
0512 
0513 if ( mvtx_refresh > 0 &&  mvtx_evnts%mvtx_refresh == 0 ) OM();
0514 
0515 if ( mvtx_evnts%100000 == 0 )
0516     cout << " processing event " << mvtx_evnts << endl;
0517 
0518     mvtx_evnts++;
0519     return 0;
0520     }
0521 
0522 //============================================================//
0523 
0524 int process_histos(float thresh)
0525 {
0526     int sum = 0;
0527     int row,col;
0528     for (int i = 0; i < NSTAVE; i++)
0529     {
0530         cout << endl;
0531         cout << "-- stave " << i << endl;
0532         int tot = 0;
0533         for (int ix = 1; ix <= h2d_chip[i]->GetNbinsX(); ix++)
0534         {
0535             for (int iy = 1; iy <= h2d_chip[i]->GetNbinsY(); iy++)
0536             {
0537                 if ( h2d_chip[i]->GetBinContent(ix, iy) > thresh )
0538                 {
0539                     tot++;
0540                     if (flip_yaxis)
0541                     {
0542                         row = 1023-(h2d_chip[i]->GetYaxis()->GetBinCenter(iy));
0543                     } else
0544                     {
0545                         row = h2d_chip[i]->GetYaxis()->GetBinCenter(iy);
0546                     }
0547                     col = h2d_chip[i]->GetXaxis()->GetBinCenter(ix);
0548                     cout << " row:" << row
0549                         << " col:" << col
0550                         << endl;
0551                 }
0552             }
0553         }
0554         cout << "Total: " << tot << endl;
0555         sum += tot;
0556     }
0557 
0558     return sum;
0559 }
0560 
0561 //============================================================//
0562 
0563 int mask_pixels(float thresh)
0564 {
0565     bool flip_yaxis = true; //Only uncomment if we need to flip the yAxis
0566 
0567     string file_name_tb1 = "/home/maps/git/RUv1_Test_sync2018-08/software/py/masklist_testbench1.txt";
0568     string file_name_tb2 = "/home/maps/git/RUv1_Test_sync2018-08/software/py/masklist_testbench2.txt";
0569     ofstream write_mask_file_tb1(file_name_tb1.c_str());
0570     ofstream write_mask_file_tb2(file_name_tb2.c_str());
0571 
0572     write_mask_file_tb1<<"#Connector, ChipID, Col, Row"<<endl;
0573     write_mask_file_tb2<<"#Connector, ChipID, Col, Row"<<endl;
0574 
0575     int sum = 0;
0576     int chipid, row, col, tot, chip_tot, prev_tot;
0577     Int_t stave_map[4] = {0, 4, 1, 2};
0578     for (int i = 0; i < NSTAVE; i++)
0579     {
0580         tot = 0;
0581         prev_tot = 0;
0582         for (int ix = 1; ix <= h2d_chip[i]->GetNbinsX(); ix++)
0583         {
0584             for (int iy = 1; iy <= h2d_chip[i]->GetNbinsY(); iy++)
0585             {
0586                 col = h2d_chip[i]->GetXaxis()->GetBinCenter(ix);
0587                 chipid = col/1024;
0588                 col = col%1024;
0589 
0590                 if ( h2d_chip[i]->GetBinContent(ix, iy) > thresh )
0591                 {
0592                     tot++;
0593                     if (flip_yaxis)
0594                     {
0595                         row = 511-(h2d_chip[i]->GetYaxis()->GetBinCenter(iy));
0596                     }
0597                     else
0598                     {
0599                         row = h2d_chip[i]->GetYaxis()->GetBinCenter(iy);
0600                     }
0601                     if (i != 1) { write_mask_file_tb1 << stave_map[i] << "," << chipid  << "," << col << "," << row << endl; }
0602                     else { write_mask_file_tb2 << stave_map[i] << "," << chipid  << "," << col << "," << row << endl; }
0603                 }
0604             }
0605             chip_tot = tot - prev_tot;
0606             if (col == 1023) {prev_tot = tot; printf("Total pixels masked on stave %i, chip %i: %i\n", i, chipid, chip_tot);}
0607         }
0608         sum += tot;
0609     }
0610     printf("Total pixels masked: %i\n", sum);
0611     write_mask_file_tb1.close();
0612     write_mask_file_tb2.close();
0613     printf("New pixel mask has been created\n");
0614 
0615     return sum;
0616 }
0617 
0618 //============================================================//
0619 
0620 int analysis()
0621 {
0622     //TH1D* hhitrate_chip;
0623     //TH1D* hnhit_chip_norm[NSTAVE];
0624     //TH2D* h2d_chip_norm[NSTAVE];
0625     //TH1F* hdiffrow_chip_norm[NSTAVE];
0626     //TH1F* hdiffcol_chip_norm[NSTAVE];
0627     char name[500];
0628 
0629     //-- Check if we've already initialized the histograms
0630     //   if not, initialize
0631     if ( !hhitrate_chip )
0632     {
0633         hhitrate_chip = (TH1F*) hchip->Clone("hhitrate_chip");
0634         hwarnrate_chip = (TH1F*) hwarn->Clone("hwarnrate_chip");
0635         herrrate_chip = (TH1F*) herr->Clone("herrrate_chip");
0636 
0637         for (int ichip = 0; ichip < NSTAVE; ichip++)
0638         {
0639             sprintf(name, "hnhit_chip_norm_%i", ichip);
0640             hnhit_chip_norm[ichip] = (TH1F*) hnhit_chip[ichip]->Clone(name);
0641 
0642             sprintf(name, "h2d_chip_norm_%i", ichip);
0643             h2d_chip_norm[ichip] = (TH2F*) h2d_chip[ichip]->Clone(name);
0644 
0645             sprintf(name, "hdiffrow_chip_norm_%i", ichip);
0646             hdiffrow_chip_norm[ichip] = (TH1F*) hdiffrow_chip[ichip]->Clone(name);
0647 
0648             sprintf(name, "hdiffcol_chip_norm_%i", ichip);
0649             hdiffcol_chip_norm[ichip] = (TH1F*) hdiffcol_chip[ichip]->Clone(name);
0650         } // ichip
0651     }
0652 
0653     //-- current number of events
0654     double nevents = hnevnt->Integral();
0655 
0656     if ( nevents <= 0 )
0657         return 1;
0658 
0659 
0660     //-- hitrate vs chip index
0661     hhitrate_chip->Reset();
0662     hwarnrate_chip->Reset();
0663     herrrate_chip->Reset();
0664     for (int ib = 1; ib <= hchip->GetNbinsX(); ib++)
0665     {
0666         // hits
0667         float bc = hchip->GetBinContent(ib);
0668         hhitrate_chip->SetBinContent(ib, bc / nevents);
0669         hhitrate_chip->SetBinError(ib, sqrt(bc) / nevents);
0670 
0671         // warnings
0672         bc = hwarn->GetBinContent(ib);
0673         hwarnrate_chip->SetBinContent(ib, bc / nevents);
0674         hwarnrate_chip->SetBinError(ib, sqrt(bc) / nevents);
0675 
0676         // errors
0677         bc = herr->GetBinContent(ib);
0678         herrrate_chip->SetBinContent(ib, bc / nevents);
0679         herrrate_chip->SetBinError(ib, sqrt(bc) / nevents);
0680     } // ib
0681 
0682     if (show_beam_fit) { // create canvas for beamcenter plots
0683       cBeamCenter = new TCanvas("cBeamCenter", "cBeamCenter", 1350, 900);
0684       cBeamCenter->Divide(3,4);
0685     }
0686 
0687     //-- per chip histograms
0688     for (int ichip = 0; ichip < NSTAVE; ichip++)
0689     {
0690         //-- nhit distribution for each chip
0691         hnhit_chip_norm[ichip]->Reset();
0692         for (int ib = 1; ib <= hnhit_chip[ichip]->GetNbinsX(); ib++)
0693         {
0694             float bc = hnhit_chip[ichip]->GetBinContent(ib);
0695             hnhit_chip_norm[ichip]->SetBinContent(ib, bc / nevents);
0696             hnhit_chip_norm[ichip]->SetBinError(ib, sqrt(bc) / nevents);
0697         } // ib
0698 
0699         //-- 2D distribution for each chip
0700         h2d_chip_norm[ichip]->Reset();
0701         for (int ix = 1; ix <= h2d_chip[ichip]->GetNbinsX(); ix++)
0702             for (int iy = 1; iy <= h2d_chip[ichip]->GetNbinsY(); iy++)
0703             {
0704                 float bc = h2d_chip[ichip]->GetBinContent(ix, iy);
0705                 h2d_chip_norm[ichip]->SetBinContent(ix, iy, bc / nevents);
0706             }
0707 
0708         //-- mean pixel index for each chip
0709         double m, s;
0710         TH1D* h;
0711 
0712         if ( cBeamCenter )
0713         {
0714           cBeamCenter->cd(3*ichip+1);
0715           gPad->SetLogz();
0716           TH2D* h2d_beam = (TH2D*)h2d_chip_norm[ichip]->Clone("h2d_beam");
0717           h2d_beam->GetXaxis()->SetTitle("Stave_Cols");
0718           h2d_beam->GetYaxis()->SetTitle("Stave_Rows");
0719           h2d_beam->Draw("colz");
0720           TLatex lt;
0721           lt.SetNDC();
0722           lt.SetTextAlign(22);
0723           lt.SetTextSize(1.5*lt.GetTextSize());
0724           lt.SetTextColor(1);
0725           lt.DrawText(0.5,0.96, stave_name[ichip].c_str());
0726 
0727           switch ( ichip )
0728           {
0729             case 0 :
0730               dead_chip_forward[0]->Draw();
0731               dead_chip_backward[0]->Draw();
0732               break;
0733 
0734             case 2 :
0735               dead_chip_forward[0]->Draw(); dead_chip_backward[0]->Draw();
0736               dead_chip_forward[3]->Draw(); dead_chip_backward[3]->Draw();
0737               dead_chip_forward[6]->Draw(); dead_chip_backward[6]->Draw();
0738               dead_chip_forward[7]->Draw(); dead_chip_backward[7]->Draw();
0739               dead_chip_forward[8]->Draw(); dead_chip_backward[8]->Draw();
0740               break;
0741 
0742             case 3 :
0743               dead_chip_forward[6]->Draw(); dead_chip_backward[6]->Draw();
0744               break;
0745 
0746             default :
0747              break;
0748           }
0749         }
0750 
0751         // cols projection
0752         h = (TH1D*) h2d_chip_norm[ichip]->ProjectionX();
0753         double m_col = h->GetMean();
0754         double r_col = h->GetRMS();
0755         fg->SetParameters(h->GetBinContent(h->GetMaximumBin()), m_col, r_col);
0756         //YCM: initial range to full stave range
0757         const int lastBinX = 9 * 1024 - 1;
0758         const int lastBinY = 511;
0759         fg->SetRange(0, lastBinX);
0760         h->Fit(fg, "RQ0N");
0761         m = fg->GetParameter(1);
0762         s = fg->GetParameter(2);
0763         //cout << " x1 mean:" << m << " sig:" << s << endl;
0764         m = m - 2.5*s > lastBinX ? lastBinX : m;
0765         fg->SetRange(m - 2.5*s, m + 2.5*s);
0766         fg->SetParameters(fg->GetParameter(0), m, s);
0767         h->Fit(fg, "RQ0N");
0768         m = fg->GetParameter(1);
0769         s = fg->GetParameter(2);
0770         //cout << " x2 mean:" << m << " sig:" << s << endl;
0771         if ( m > lastBinX )
0772         {
0773             s = s - (m - lastBinX);
0774             m = lastBinX;
0775         }
0776         if ( m < 0 )
0777         {
0778             s = s - m;
0779             m = 0;
0780         }
0781         //cout << " (m:" << m << " s:" << s << ")"<< endl;
0782 
0783         // save
0784         m_col = m;
0785         r_col = s;
0786         if ( cBeamCenter )
0787         {
0788           cBeamCenter->cd( 3 * ichip + 2 );
0789           h->GetXaxis()->SetTitle("Stave_Cols");
0790           h->GetYaxis()->SetTitle("nClusters");
0791           h->Draw();
0792           TF1* f_col = (TF1*)fg->Clone("fx");
0793           f_col->Draw("same");
0794           TLatex lt;
0795           lt.SetNDC();
0796           lt.SetTextAlign(22);
0797           lt.SetTextSize(1.5*lt.GetTextSize());
0798           lt.SetTextColor(2);
0799           lt.DrawLatex(.3, .7, Form("#mu_{COL} = %.2f #pm %.2f ", m_col, r_col));
0800           lt.DrawLatex(.3, .6, Form("offset = %.2f", m_col));
0801         }
0802         else
0803           delete h;
0804 
0805         // row projection
0806         h = (TH1D*) h2d_chip_norm[ichip]->ProjectionY();
0807         double m_row = h->GetMean();
0808         double r_row = h->GetRMS();
0809         fg->SetParameters(h->GetBinContent(h->GetMaximumBin()), m_row, r_row);
0810         fg->SetRange(0, 1023);
0811         h->Fit(fg, "RQ0N");
0812         m = fg->GetParameter(1);
0813         s = fg->GetParameter(2);
0814         //cout << " y2 mean:" << m << " sig:" << s << endl;
0815         fg->SetRange(m - 2.5*s, m + 2.5*s);
0816         fg->SetParameters(fg->GetParameter(0), m, s);
0817         h->Fit(fg, "RQ0N");
0818         m = fg->GetParameter(1);
0819         s = fg->GetParameter(2);
0820         //cout << " y2 mean:" << m << " sig:" << s << endl;
0821         if ( m > lastBinY )
0822         {
0823             s = s - (m - lastBinY);
0824             m = lastBinY;
0825         }
0826         if ( m < 0 )
0827         {
0828             s = s - m;
0829             m = 0;
0830         }
0831         //cout << " (m:" << m << " s:" << s << ")"<< endl;
0832 
0833         // save
0834         m_row = m;
0835         r_row = s;
0836         if ( cBeamCenter )
0837         {
0838           cBeamCenter->cd( 3 * ichip + 3 );
0839           h->GetXaxis()->SetTitle("Stave_Cols");
0840           h->GetYaxis()->SetTitle("nClusters");
0841           h->Draw();
0842           TF1* f_row = (TF1*)fg->Clone("fy");
0843           f_row->Draw("same");
0844           TLatex lt;
0845           lt.SetNDC();
0846           lt.SetTextAlign(22);
0847           lt.SetTextSize(1.5*lt.GetTextSize());
0848           lt.SetTextColor(2);
0849           lt.DrawLatex(.3, .7, Form("#mu_{ROW}  = %.2f #pm %.2f ", m_row, r_row));
0850           lt.DrawLatex(.3, .6, Form("offset = %.2f", lastBinY - m_row));
0851         }
0852         else
0853           delete h;
0854 
0855         gmpos_chip[ichip]->SetPoint(0, m_col, m_row);
0856         gmpos_chip[ichip]->SetPointError(0, r_col, r_col, r_row, r_row);
0857 
0858         //-- difference in mean pixel index per event for each chip
0859         hdiffrow_chip_norm[ichip]->Reset();
0860         hdiffcol_chip_norm[ichip]->Reset();
0861         for (int ib = 1; ib <= hdiffrow_chip[ichip]->GetNbinsX(); ib++)
0862         {
0863             hdiffrow_chip_norm[ichip]->SetBinContent(ib, hdiffrow_chip[ichip]->GetBinContent(ib));
0864             hdiffcol_chip_norm[ichip]->SetBinContent(ib, hdiffcol_chip[ichip]->GetBinContent(ib));
0865         }
0866 
0867         float introw = hdiffrow_chip_norm[ichip]->Integral();
0868         if ( introw > 0 )
0869             hdiffrow_chip_norm[ichip]->Scale(1./introw);
0870 
0871         float intcol = hdiffcol_chip_norm[ichip]->Integral();
0872         if ( intcol > 0 )
0873             hdiffcol_chip_norm[ichip]->Scale(1./intcol);
0874     } // ichip
0875 
0876     return 0;
0877 }
0878 
0879 //============================================================//
0880 
0881 int OM()
0882 {
0883     //-- run analysis
0884     analysis();
0885 
0886 
0887 
0888     gStyle->SetOptStat(0);
0889     char name[500];
0890     //-- setup canvas
0891     static bool initialized = false;
0892 
0893     if ( !initialized )
0894     {
0895 
0896         //com3 = new TCanvas("com","MVTX online monitoring", 1600, 800);
0897         com = new TCanvas("com","MVTX online monitoring", 2000, 800);
0898         com->SetMargin(0, 0, 0, 0);
0899 
0900         // hitrate vs chip index
0901         //phitrate = new TPad("phitrate", "hit rate", 0.0, 0.45, 0.29, 0.9);
0902         phitrate = new TPad("phitrate", "hit rate", 0.0, 0.60, 0.29, 0.89);
0903         phitrate->SetMargin(0.15, 0.05, 0.15, 0.05);
0904         phitrate->SetTicks(1, 1);
0905         phitrate->Draw();
0906 /*
0907         // nhit distributions
0908         com->cd();
0909         pnhit = new TPad("pnhit", "n hit", 0.29, 0.60, 0.50, 0.90);
0910         pnhit->SetMargin(0.12, 0.02, 0.15, 0.05);
0911         pnhit->SetTicks(1, 1);
0912         pnhit->Draw();
0913 
0914         // mean hit positions
0915         com->cd();
0916         pmean = new TPad("pmean", "mean", 0.00, 0.00, 0.29, 0.45);
0917         pmean->SetMargin(0.15, 0.05, 0.15, 0.10);
0918         pmean->SetTicks(1, 1);
0919         pmean->Draw();
0920 */
0921         // diff col index
0922         com->cd();
0923         //pdiffcol = new TPad("pdiffcol", "pdiffcol", 0.29, 0.30, 0.50, 0.60);
0924         pdiffcol = new TPad("pdiffcol", "pdiffcol", 0.00, 0.30, 0.29, 0.59);
0925         //pdiffcol->SetMargin(0.12, 0.02, 0.12, 0.10);
0926         pdiffcol->SetMargin(0.15, 0.05, 0.15, 0.05);
0927         pdiffcol->SetTicks();
0928         pdiffcol->Draw();
0929 
0930         // diff row index
0931         com->cd();
0932         //pdiffrow = new TPad("pdiffrow", "pdiffrow", 0.29, 0.00, 0.50, 0.30);
0933         pdiffrow = new TPad("pdiffrow", "pdiffrow", 0.00, 0.00, 0.29, 0.29);
0934         //pdiffrow->SetMargin(0.12, 0.02, 0.12, 0.10);
0935         pdiffrow->SetMargin(0.15, 0.05, 0.15, 0.05);
0936         pdiffrow->SetTicks();
0937         pdiffrow->Draw();
0938 
0939         // 2D hit distributions
0940         for (int i = 0; i < NSTAVE; i++)
0941         {
0942             com->cd();
0943             sprintf(name, "p2d_%i", i);
0944             //p2d[i] = new TPad(name, name,
0945             //        0.50, 0.9-0.22*i,
0946             //        1.00, 0.9-0.22*(i+1) );
0947             p2d[i] = new TPad(name, name,
0948                     0.30, 0.9-0.22*i,
0949                     1.00, 0.9-0.22*(i+1) );
0950             //p2d[i]->SetMargin(0.15, 0.15, 0.15, 0.15);
0951             p2d[i]->SetMargin(0.01, 0.02, 0.25, 0.25);
0952             p2d[i]->SetTicks(1, 1);
0953             p2d[i]->Draw();
0954         } // i
0955 
0956         //-- info
0957         com->cd();
0958         sprintf(name, "Number of Events: %.0f", hnevnt->Integral());
0959         lnevents = new TLatex(0.5, 0.95, name);
0960         lnevents->SetNDC();
0961         lnevents->SetTextAlign(22);
0962         lnevents->Draw();
0963 
0964         TLatex lt;
0965         lt.SetNDC();
0966         lt.SetTextAlign(22);
0967 
0968         //-- hitrate
0969         phitrate->cd();
0970         //gPad->SetLogy();
0971         haxis_chip->GetYaxis()->SetRangeUser(0,5);
0972         haxis_chip->Draw();
0973         hhitrate_chip->Draw("hist e same");
0974         hwarnrate_chip->Draw("hist e same");
0975         herrrate_chip->Draw("hist e same");
0976 
0977         //-- nhit distributions
0978         //pnhit->cd();
0979         //gPad->SetLogy();
0980         haxis_nhit->GetYaxis()->SetRangeUser(1e-4, 1);
0981         //haxis_nhit->DrawCopy();
0982         for (int i = 0; i < NSTAVE; i++)
0983         {
0984             hnhit_chip_norm[i]->Draw("hist e same");
0985 
0986             sprintf(name, "chip %i = %.2f", i, hnhit_chip_norm[i]->GetMean());
0987             lnhitmean[i] = new TLatex(0.8, 0.85 - 0.05*i, name);
0988             lnhitmean[i]->SetNDC();
0989             lnhitmean[i]->SetTextColor(chipColor[i]);
0990             lnhitmean[i]->Draw("same");
0991         }
0992 
0993         //-- mean positions
0994         //pmean->cd();
0995         //haxis_2d->DrawCopy();
0996         //gPad->Update();
0997 /*        if (flip_yaxis)
0998         {
0999             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 
1143 int print_canvas()
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 
1155 int print_status()
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 
1208 void reset_histos()
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 
1225 void get_alignment()
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 }