Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 12:45:42

0001 //#include "compareGeometries.h"
0002 #include <string>
0003 #include <sstream>
0004 
0005 #include "TProfile.h"
0006 #include "TList.h"
0007 #include "TNtuple.h"
0008 #include "TString.h"
0009 
0010 double arrowSize = 0.0095;
0011 float y_,x_,z_,phi_,r_,dphi_,dr_,dx_,dz_,dy_;
0012 int level_,sublevel_;
0013 char outputDir_[192];
0014 
0015 void Plot10Mu(const char* text,float X, float Y, float size)
0016 {
0017   TPaveText* atext = new TPaveText(X,Y,X+size,Y+size);
0018   atext->AddText(text);
0019   atext->SetLineColor(0);
0020   atext->SetFillColor(0);
0021   atext->SetTextFont(42);
0022   atext->SetTextSize(0.04);
0023   atext->Draw();
0024 }
0025 
0026 
0027 void normArrow(float x, float y, float norm)
0028 {
0029   // draw 100 mu m arrow if norm = 1
0030   TArrow* normArrow = new TArrow(x,y,x+norm,y,arrowSize,">");
0031   // normArrow->SetLineWidth(2);
0032   normArrow->Draw();
0033 }
0034 
0035 void DrawRPhiLegend(double xLim, double yLim, double barrelRPhiRescale)
0036 {
0037   float xTest = 0.9*xLim;
0038   float yTest = yLim/2;
0039   float testBlockSize = 0.2*xLim; //5cm in axis unit
0040   float disty = 0;
0041   float dYTest =0.1*xLim;
0042   float dZTest =2;
0043   float xLegObj = 20;
0044   float yLegObj = 18;
0045 
0046   Plot10Mu("#Delta r:",xTest,yTest,testBlockSize);
0047   Plot10Mu("500 #mum",xTest,yTest-3*dYTest,testBlockSize);
0048   normArrow(xTest+dYTest,yTest-4*dYTest-disty,500./10000*barrelRPhiRescale);
0049 }
0050 
0051 
0052 int makeRPhiArrowPlot( TTree* data, const char* name, double xLim, double yLim, double level, double sublevel, double zMin, double zMax, double rMin, double rMax, double barrelRPhiRescale){
0053 
0054 
0055   TCanvas* OBPCanvas = new TCanvas(name,name,1050,875);
0056   OBPCanvas->DrawFrame(-xLim, -yLim, 1.2*xLim, yLim, ";module position x [cm];module position y [cm]");
0057   OBPCanvas->SetFillColor(0);
0058   OBPCanvas->SetFrameBorderMode(0);
0059 
0060   TFrame* aFrame = OBPCanvas->GetFrame();
0061   aFrame->SetFillColor(0);
0062 
0063   int passcut = 0;
0064   for(int entry = 0;entry<data->GetEntries(); entry++)
0065     {
0066       data->GetEntry(entry);
0067       if ((level_ == level)&&(((sublevel_ == sublevel)&&(sublevel != 0))||(sublevel == 0))){
0068         if ((z_ <= zMax)&&(z_ > zMin)&&(r_ <= rMax)&&(r_ > rMin)){
0069           TArrow* aArraw = new TArrow( x_, y_ , x_ + barrelRPhiRescale*dx_, y_+barrelRPhiRescale*dy_,0.0075,">");
0070           aArraw->Draw();
0071           passcut++;
0072         }
0073       }
0074     }
0075   DrawRPhiLegend( xLim, yLim, barrelRPhiRescale );
0076 
0077   char sliceLeg[192];
0078   sprintf( sliceLeg, "%s: %f < z <= %f", name, zMin, zMax );
0079   //Plot10Mu( name, xLim/2, yLim, 0.2*xLim );
0080   TPaveText* atext = new TPaveText(0.2*xLim,0.85*yLim,0.66*xLim,0.99*yLim);
0081   atext->AddText(sliceLeg);
0082   atext->SetLineColor(0);
0083   atext->SetFillColor(0);
0084   atext->SetTextFont(42);
0085   atext->SetTextSize(0.04);
0086   atext->Draw();
0087 
0088 
0089 
0090   char outfile[192];
0091   sprintf( outfile, "%s/%s.png", outputDir_, name );
0092   OBPCanvas->Print( outfile );
0093   sprintf(outfile, "%s/%s.pdf", outputDir_, name);
0094   OBPCanvas->SaveAs(outfile);
0095 
0096   return passcut;
0097 }
0098 
0099 int makeZPhiArrowPlot( TTree* data, const char* name, double zLim, double phiLim, double level, double sublevel, double zMin, double zMax, double rMin, double rMax, double barrelRPhiRescale){
0100 
0101 
0102   TCanvas* OBPCanvas = new TCanvas(name,name,1050,875);
0103   OBPCanvas->DrawFrame(-zLim, -phiLim, 1.2*zLim, phiLim, ";module position z [cm];module position r*phi [cm]");
0104   OBPCanvas->SetFillColor(0);
0105   OBPCanvas->SetFrameBorderMode(0);
0106 
0107   TFrame* aFrame = OBPCanvas->GetFrame();
0108   aFrame->SetFillColor(0);
0109 
0110   int passcut = 0;
0111   for(int entry = 0;entry<data->GetEntries(); entry++)
0112     {
0113       data->GetEntry(entry);
0114       if ((level_ == level)&&(((sublevel_ == sublevel)&&(sublevel != 0))||(sublevel == 0))){
0115         if ((z_ <= zMax)&&(z_ > zMin)&&(r_ <= rMax)&&(r_ > rMin)){
0116           TArrow* aArraw = new TArrow( z_, r_*phi_ , z_ + barrelRPhiRescale*dz_, r_*phi_+barrelRPhiRescale*r_*dphi_,0.0075,">");
0117           aArraw->Draw();
0118           passcut++;
0119         }
0120       }
0121     }
0122   DrawRPhiLegend( zLim, phiLim, barrelRPhiRescale );
0123 
0124   char sliceLeg[192];
0125   sprintf( sliceLeg, "%s: %f < r <= %f", name, rMin, rMax );
0126   //Plot10Mu( name, xLim/2, yLim, 0.2*xLim );
0127   TPaveText* atext = new TPaveText(0.2*zLim,0.85*phiLim,0.66*zLim,0.99*phiLim);
0128   atext->AddText(sliceLeg);
0129   atext->SetLineColor(0);
0130   atext->SetFillColor(0);
0131   atext->SetTextFont(42);
0132   atext->SetTextSize(0.04);
0133   atext->Draw();
0134 
0135 
0136 
0137   char outfile[192];
0138   sprintf( outfile, "%s/%s.png", outputDir_, name );
0139   OBPCanvas->Print( outfile );
0140   sprintf(outfile, "%s/%s.pdf", outputDir_, name);
0141   OBPCanvas->SaveAs(outfile);
0142 
0143   return passcut;
0144 }
0145 
0146 /*
0147   int makeRZArrowPlot( TTree* data, char* name, double zLim, double zLimMax, double rLim, double rLimMax, double level, double sublevel, double zMin, double zMax, double rMin, double rMax, double barrelRPhiRescale){
0148 
0149 
0150   TCanvas* OBPCanvas = new TCanvas(name,name,1050,875);
0151   OBPCanvas->DrawFrame(zLim, rLim, zLimMax, rLimMax, ";module position z [cm];module position r [cm]");
0152   OBPCanvas->SetFillColor(0);
0153   OBPCanvas->SetFrameBorderMode(0);
0154 
0155   TFrame* aFrame = OBPCanvas->GetFrame();
0156   aFrame->SetFillColor(0);
0157 
0158   int passcut = 0;
0159   for(int entry = 0;entry<data->GetEntries(); entry++)
0160   {
0161   data->GetEntry(entry);
0162   if ((level_ == level)&&(sublevel_ == sublevel)){
0163   if ((z_ <= zMax)&&(z_ > zMin)&&(r_ <= rMax)&&(r_ > rMin)){
0164   TArrow* aArraw = new TArrow( z_, r_ , z_ + barrelRPhiRescale*dz_, r_+barrelRPhiRescale*dr_,0.0075,">");
0165   aArraw->Draw();
0166   passcut++;
0167   }
0168   }
0169   }
0170   // legend
0171   double xpos = 0.8*(zLimMax-zLim) + zLim;
0172   double ypos = 0.7*(rLimMax-rLim) + rLim;
0173   double sizer = 0.05*(zLimMax-zLim);
0174   Plot10Mu("#Delta r:",xpos,ypos,sizer);
0175   Plot10Mu("500 #mum",xpos,ypos-1.5*sizer,sizer);
0176   normArrow(xpos,ypos-2*sizer,500./10000*barrelRPhiRescale);
0177 
0178 
0179   char sliceLeg[192];
0180   sprintf( sliceLeg, "%s: %d < z <= %d", name, zMin, zMax );
0181   //Plot10Mu( name, xLim/2, yLim, 0.2*xLim );
0182   TPaveText* atext = new TPaveText(0.4*(zLimMax-zLim) + zLim,0.85*(rLimMax-rLim)+rLim,0.86*(zLimMax-zLim) + zLim,0.99*(rLimMax-rLim)+rLim);
0183   atext->AddText(sliceLeg);
0184   atext->SetLineColor(0);
0185   atext->SetFillColor(0);
0186   atext->SetTextFont(42);
0187   atext->SetTextSize(0.04);
0188   atext->Draw();
0189 
0190 
0191 
0192   char outfile[192];
0193   sprintf( outfile, "%s/%s.png", outputDir_, name );
0194   OBPCanvas->Print( outfile );
0195 
0196   return passcut;
0197   }
0198 */
0199 
0200 void makeArrowPlots(const char* filename, const char* outputDir)
0201 {
0202 
0203   TFile fin(filename);
0204   fin.cd();
0205 
0206   bool plotPXB = true;
0207   bool plotTIB = true;
0208   bool plotTOB = true;
0209   bool plotPXF = true;
0210   bool plotTID = true;
0211   bool plotTEC = true;
0212 
0213   TString outputfile("OUTPUT_");
0214   outputfile.Append(filename);
0215   TFile output(outputfile, "recreate");
0216 
0217   sprintf( outputDir_, "%s", outputDir );
0218   gSystem->mkdir(outputDir_, true);
0219 
0220   TTree* data = static_cast<TTree*>(fin.Get("alignTree"));
0221   data->SetBranchAddress("sublevel",&sublevel_);
0222   data->SetBranchAddress("level",&level_);
0223   data->SetBranchAddress("x",&x_);
0224   data->SetBranchAddress("y",&y_);
0225   data->SetBranchAddress("z",&z_);
0226   data->SetBranchAddress("dx",&dx_);
0227   data->SetBranchAddress("dy",&dy_);
0228   data->SetBranchAddress("dz",&dz_);
0229   data->SetBranchAddress("dphi",&dphi_);
0230   data->SetBranchAddress("dr",&dr_);
0231   data->SetBranchAddress("phi",&phi_);
0232   data->SetBranchAddress("r",&r_);
0233 
0234   // args are: tree, title, xLim, yLim
0235   // cuts are: level, sublevel, zMin, zMax, rMin, rMax, scaleFactor
0236   // PXB slices
0237   int totalPXB_modules = 0;
0238   int totalPXB_modules_zphi = 0;
0239   int dummy = 0;
0240   if (plotPXB){
0241     double pxbScale = 30.0;
0242     totalPXB_modules += makeRPhiArrowPlot( data, "PXB_BarrelXY-4", 20, 20, 1, 1, -26, -20, 0, 19, pxbScale);
0243     totalPXB_modules += makeRPhiArrowPlot( data, "PXB_BarrelXY-3", 20, 20, 1, 1, -19, -14, 0, 19, pxbScale);
0244     totalPXB_modules += makeRPhiArrowPlot( data, "PXB_BarrelXY-2", 20, 20, 1, 1, -14, -6.5, 0, 19, pxbScale);
0245     totalPXB_modules += makeRPhiArrowPlot( data, "PXB_BarrelXY-1", 20, 20, 1, 1, -6.5, 0, 0, 19, pxbScale);
0246     totalPXB_modules += makeRPhiArrowPlot( data, "PXB_BarrelXY+1", 20, 20, 1, 1, 0, 6.5, 0, 19, pxbScale);
0247     totalPXB_modules += makeRPhiArrowPlot( data, "PXB_BarrelXY+2", 20, 20, 1, 1, 6.5, 14, 0, 19, pxbScale);
0248     totalPXB_modules += makeRPhiArrowPlot( data, "PXB_BarrelXY+3", 20, 20, 1, 1, 14, 19, 0, 19, pxbScale);
0249     totalPXB_modules += makeRPhiArrowPlot( data, "PXB_BarrelXY+4", 20, 20, 1, 1, 19, 26, 0, 19, pxbScale);
0250     double pxbScale_zphi = 40.0;
0251     totalPXB_modules_zphi += makeZPhiArrowPlot( data, "PXB_BarrelZPhi_1", 35, 15, 1, 1, -300, 300, 0, 5, pxbScale_zphi);
0252     totalPXB_modules_zphi += makeZPhiArrowPlot( data, "PXB_BarrelZPhi_2", 35, 30, 1, 1, -300, 300, 5, 9, pxbScale_zphi);
0253     totalPXB_modules_zphi += makeZPhiArrowPlot( data, "PXB_BarrelZPhi_3", 35, 45, 1, 1, -300, 300, 9, 14, pxbScale_zphi);
0254     totalPXB_modules_zphi += makeZPhiArrowPlot( data, "PXB_BarrelZPhi_4", 35, 65, 1, 1, -300, 300, 14, 19, pxbScale_zphi);
0255   }
0256 
0257   // TIB slices
0258   int totalTIB_modules = 0;
0259   int totalTIB_modules_zphi = 0;
0260   if (plotTIB){
0261     double tibScale = 30.0;
0262     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY-6", 80, 80, 1, 3, -70, -56, 0, 120, tibScale);
0263     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY-5", 80, 80, 1, 3, -56, -42, 0, 120, tibScale);
0264     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY-4", 80, 80, 1, 3, -42, -32, 0, 120, tibScale);
0265     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY-3", 80, 80, 1, 3, -32, -20, 0, 120, tibScale);
0266     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY-2", 80, 80, 1, 3, -20, -10, 0, 120, tibScale);
0267     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY-1", 80, 80, 1, 3, -10, 0, 0, 120, tibScale);
0268     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY+1", 80, 80, 1, 3, 0, 10, 0, 120, tibScale);
0269     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY+2", 80, 80, 1, 3, 10, 20, 0, 120, tibScale);
0270     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY+3", 80, 80, 1, 3, 20, 32, 0, 120, tibScale);
0271     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY+4", 80, 80, 1, 3, 32, 42, 0, 120, tibScale);
0272     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY+5", 80, 80, 1, 3, 42, 56, 0, 120, tibScale);
0273     totalTIB_modules += makeRPhiArrowPlot( data, "TIB_BarrelXY+6", 80, 80, 1, 3, 56, 70, 0, 120, tibScale);
0274     double tibScale_zphi = 40.0;
0275     totalTIB_modules_zphi += makeZPhiArrowPlot( data, "TIB_BarrelZPhi_1", 80, 120, 1, 3, -300, 300, 20.0, 29.0, tibScale_zphi);
0276     totalTIB_modules_zphi += makeZPhiArrowPlot( data, "TIB_BarrelZPhi_2", 80, 140, 1, 3, -300, 300, 29.0, 37.5, tibScale_zphi);
0277     totalTIB_modules_zphi += makeZPhiArrowPlot( data, "TIB_BarrelZPhi_3", 80, 170, 1, 3, -300, 300, 37.5, 47.5, tibScale_zphi);
0278     totalTIB_modules_zphi += makeZPhiArrowPlot( data, "TIB_BarrelZPhi_4", 80, 200, 1, 3, -300, 300, 47.5, 60.0, tibScale_zphi);
0279 
0280   }
0281 
0282   // TOB slices
0283   int totalTOB_modules = 0;
0284   int totalTOB_modules_zphi = 0;
0285   if (plotTOB){
0286     double tobScale = 100.0;
0287     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY-6", 150, 150, 1, 5, -108, -90, 0, 120, tobScale);
0288     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY-5", 150, 150, 1, 5, -90, -72, 0, 120, tobScale);
0289     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY-4", 150, 150, 1, 5, -72, -54, 0, 120, tobScale);
0290     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY-3", 150, 150, 1, 5, -54, -36, 0, 120, tobScale);
0291     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY-2", 150, 150, 1, 5, -36, -18, 0, 120, tobScale);
0292     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY-1", 150, 150, 1, 5, -18, 0, 0, 120, tobScale);
0293     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY+1", 150, 150, 1, 5, 0, 18, 0, 120, tobScale);
0294     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY+2", 150, 150, 1, 5, 18, 36, 0, 120, tobScale);
0295     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY+3", 150, 150, 1, 5, 36, 54, 0, 120, tobScale);
0296     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY+4", 150, 150, 1, 5, 54, 72, 0, 120, tobScale);
0297     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY+5", 150, 150, 1, 5, 72, 90, 0, 120, tobScale);
0298     totalTOB_modules += makeRPhiArrowPlot( data, "TOB_BarrelXY+6", 150, 150, 1, 5, 90, 108, 0, 120, tobScale);
0299     double tobScale_zphi = 40.0;
0300     totalTOB_modules_zphi += makeZPhiArrowPlot( data, "TOB_BarrelZPhi_1", 130, 250, 1, 5, -300, 300, 55.0, 65.0, tobScale_zphi);
0301     totalTOB_modules_zphi += makeZPhiArrowPlot( data, "TOB_BarrelZPhi_2", 130, 280, 1, 5, -300, 300, 65.0, 75.0, tobScale_zphi);
0302     totalTOB_modules_zphi += makeZPhiArrowPlot( data, "TOB_BarrelZPhi_3", 130, 320, 1, 5, -300, 300, 75.0, 85.0, tobScale_zphi);
0303     totalTOB_modules_zphi += makeZPhiArrowPlot( data, "TOB_BarrelZPhi_4", 130, 350, 1, 5, -300, 300, 85.0, 93.0, tobScale_zphi);
0304     totalTOB_modules_zphi += makeZPhiArrowPlot( data, "TOB_BarrelZPhi_5", 130, 380, 1, 5, -300, 300, 93.0, 101.0, tobScale_zphi);
0305     totalTOB_modules_zphi += makeZPhiArrowPlot( data, "TOB_BarrelZPhi_6", 130, 410, 1, 5, -300, 300, 101.0, 110.0, tobScale_zphi);
0306   }
0307 
0308   // PXF slices
0309   int totalPXF_modules = 0;
0310   int totalPXF_modules_rz = 0;
0311   if (plotPXF){
0312     double pxfScale = 30.0;
0313     totalPXF_modules += makeRPhiArrowPlot( data, "PXF_DiskXY+1", 20, 20, 1, 2, 25, 36, 0, 120, pxfScale);
0314     totalPXF_modules += makeRPhiArrowPlot( data, "PXF_DiskXY+2", 20, 20, 1, 2, 36, 44, 0, 120, pxfScale);
0315     totalPXF_modules += makeRPhiArrowPlot( data, "PXF_DiskXY+3", 20, 20, 1, 2, 44, 55, 0, 120, pxfScale);
0316     totalPXF_modules += makeRPhiArrowPlot( data, "PXF_DiskXY-1", 20, 20, 1, 2, -36, -25, 0, 120, pxfScale);
0317     totalPXF_modules += makeRPhiArrowPlot( data, "PXF_DiskXY-2", 20, 20, 1, 2, -44, -36, 0, 120, pxfScale);
0318     totalPXF_modules += makeRPhiArrowPlot( data, "PXF_DiskXY-3", 20, 20, 1, 2, -55, -44, 0, 120, pxfScale);
0319     /*
0320       double pxfScale_zphi = 10.0;
0321       totalPXF_modules_rz += makeRZArrowPlot( data, "PXF_DiskRZ-1", -38, -30, 5, 17, 1, 2, -40, -25, 0, 120.0, pxfScale_zphi);
0322       totalPXF_modules_rz += makeRZArrowPlot( data, "PXF_DiskRZ-2", -52, -44, 5, 17, 1, 2, -55, -40, 0, 120.0, pxfScale_zphi);
0323       totalPXF_modules_rz += makeRZArrowPlot( data, "PXF_DiskRZ+1", 32, 40, 5, 17, 1, 2, 25, 40, 0, 120.0, pxfScale_zphi);
0324       totalPXF_modules_rz += makeRZArrowPlot( data, "PXF_DiskRZ+2", 46, 54, 5, 17, 1, 2, 40, 55, 0, 120.0, pxfScale_zphi);
0325     */
0326   }
0327 
0328   // TID slices
0329   int totalTID_modules = 0;
0330   int totalTID_modules_rz = 0;
0331   if (plotTID){
0332     double tidScale = 50.0;
0333     totalTID_modules += makeRPhiArrowPlot( data, "TID_DiskXY+1", 70, 70, 1, 4, 70, 80, 0, 120, tidScale);
0334     totalTID_modules += makeRPhiArrowPlot( data, "TID_DiskXY+2", 70, 70, 1, 4, 80, 95, 0, 120, tidScale);
0335     totalTID_modules += makeRPhiArrowPlot( data, "TID_DiskXY+3", 70, 70, 1, 4, 95, 110, 0, 120, tidScale);
0336     totalTID_modules += makeRPhiArrowPlot( data, "TID_DiskXY-1", 70, 70, 1, 4, -80, -70, 0, 120, tidScale);
0337     totalTID_modules += makeRPhiArrowPlot( data, "TID_DiskXY-2", 70, 70, 1, 4, -95, -80, 0, 120, tidScale);
0338     totalTID_modules += makeRPhiArrowPlot( data, "TID_DiskXY-3", 70, 70, 1, 4, -110, -95, 0, 120, tidScale);
0339     /*
0340       double tidScale_zphi = 10.0;
0341       totalTID_modules_rz += makeRZArrowPlot( data, "TID_DiskRZ-1", -80, -70, -3, 55, 1, 4, -80, -70, 0, 120.0, tidScale_zphi);
0342       totalTID_modules_rz += makeRZArrowPlot( data, "TID_DiskRZ-2", -95, -80, 20, 55, 1, 4, -95, -80, 0, 120.0, tidScale_zphi);
0343       totalTID_modules_rz += makeRZArrowPlot( data, "TID_DiskRZ-3", -110, -95, 20, 55, 1, 4, -110, -95, 0, 120.0, tidScale_zphi);
0344       totalTID_modules_rz += makeRZArrowPlot( data, "TID_DiskRZ+1", 70, 80, 20, 55, 1, 4, 70, 80, 0, 120.0, tidScale_zphi);
0345       totalTID_modules_rz += makeRZArrowPlot( data, "TID_DiskRZ+2", 80, 95, 20, 55, 1, 4, 80, 95, 0, 120.0, tidScale_zphi);
0346       totalTID_modules_rz += makeRZArrowPlot( data, "TID_DiskRZ+3", 95, 110, 20, 55, 1, 4, 95, 110, 0, 120.0, tidScale_zphi);
0347     */
0348   }
0349 
0350 
0351   // TEC slices
0352   int totalTEC_modules = 0;
0353   if (plotTEC){
0354     double tecScale = 100.0;
0355     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY+1", 150, 150, 1, 6, 120, 130, 0, 120, tecScale);
0356     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY+2", 150, 150, 1, 6, 130, 145, 0, 120, tecScale);
0357     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY+3", 150, 150, 1, 6, 145, 160, 0, 120, tecScale);
0358     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY+4", 150, 150, 1, 6, 160, 175, 0, 120, tecScale);
0359     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY+5", 150, 150, 1, 6, 175, 190, 0, 120, tecScale);
0360     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY+6", 150, 150, 1, 6, 190, 215, 0, 120, tecScale);
0361     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY+7", 150, 150, 1, 6, 215, 235, 0, 120, tecScale);
0362     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY+8", 150, 150, 1, 6, 235, 260, 0, 120, tecScale);
0363     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY+9", 150, 150, 1, 6, 260, 280, 0, 120, tecScale);
0364     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY-1", 150, 150, 1, 6, -130, -120, 0, 120, tecScale);
0365     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY-2", 150, 150, 1, 6, -145, -130, 0, 120, tecScale);
0366     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY-3", 150, 150, 1, 6, -160, -145, 0, 120, tecScale);
0367     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY-4", 150, 150, 1, 6, -175, -160, 0, 120, tecScale);
0368     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY-5", 150, 150, 1, 6, -190, -175, 0, 120, tecScale);
0369     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY-6", 150, 150, 1, 6, -215, -190, 0, 120, tecScale);
0370     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY-7", 150, 150, 1, 6, -235, -215, 0, 120, tecScale);
0371     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY-8", 150, 150, 1, 6, -260, -235, 0, 120, tecScale);
0372     totalTEC_modules += makeRPhiArrowPlot( data, "TEC_DiskXY-9", 150, 150, 1, 6, -280, -260, 0, 120, tecScale);
0373   }
0374 
0375   std::cout << "Plotted PXB modules: " << totalPXB_modules << std::endl;
0376   std::cout << "Plotted PXB modules (zphi): " << totalPXB_modules_zphi << std::endl;
0377   std::cout << "Plotted TIB modules: " << totalTIB_modules << std::endl;
0378   std::cout << "Plotted TIB modules (zphi): " << totalTIB_modules_zphi << std::endl;
0379   std::cout << "Plotted TOB modules: " << totalTOB_modules << std::endl;
0380   std::cout << "Plotted TOB modules (zphi): " << totalTOB_modules_zphi << std::endl;
0381   std::cout << "Plotted PXF modules: " << totalPXF_modules << std::endl;
0382   //std::cout << "Plotted PXF modules (rz): " << totalPXF_modules_rz << std::endl;
0383   std::cout << "Plotted TID modules: " << totalTID_modules << std::endl;
0384   //std::cout << "Plotted TID modules (rz): " << totalTID_modules_rz << std::endl;
0385   std::cout << "Plotted TEC modules: " << totalTEC_modules << std::endl;
0386 
0387 
0388 
0389 
0390 }