CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_5/src/DPGAnalysis/SiStripTools/bin/OccupancyPlotMacros.cc

Go to the documentation of this file.
00001 #include "TText.h"
00002 #include "TLine.h"
00003 #include "TGaxis.h"
00004 #include "TFile.h"
00005 #include "TDirectory.h"
00006 #include "TH1F.h"
00007 #include "TProfile.h"
00008 #include "TH1D.h"
00009 #include "TList.h"
00010 #include "TBox.h"
00011 #include "TFrame.h"
00012 #include "TStyle.h"
00013 #include "TCanvas.h"
00014 #include "TColor.h"
00015 #include <cstring>
00016 #include <iostream>
00017 #include <math.h>
00018 #include "TROOT.h"
00019 #include "OccupancyPlotMacros.h"
00020 
00021 void PlotOccupancyMap(TFile* ff, const char* module, const float min, const float max, const float mmin, const float mmax, const int color) {
00022 
00023   gROOT->SetStyle("Plain");
00024 
00025   if(color == 1) {
00026     // A not-so-great color version
00027     const Int_t NRGBs = 5;
00028     const Int_t NCont = 255;
00029     Double_t stops[NRGBs] = { 0.00, 0.25, 0.50, 0.75, 1.00 };
00030     Double_t red[NRGBs]   = { 0.00, 0.00, 0.40, 1.00, 1.00 };
00031     Double_t green[NRGBs] = { 0.00, 0.40, 0.70, 0.60, 1.00 };
00032     Double_t blue[NRGBs]  = { 0.30, 0.60, 0.00, 0.00, 0.20 };
00033     TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
00034     gStyle->SetNumberContours(NCont);
00035   }
00036   else if(color==2) {
00037     // Gray scale
00038     const Int_t NRGBs = 3;
00039     const Int_t NCont = 255;
00040     Double_t stops[NRGBs] = { 0.00, 0.50, 1.00 };
00041     Double_t red[NRGBs]   = { 0.90, 0.50, 0.00};
00042     Double_t green[NRGBs] = { 0.90, 0.50, 0.00};
00043     Double_t blue[NRGBs]  = { 0.90, 0.50, 0.00};
00044     TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
00045     gStyle->SetNumberContours(NCont);
00046   }
00047   else if(color==3) {
00048     // used by Kevin in the TRK-11-001 paper
00049     const Int_t NRGBs = 7;
00050     const Int_t NCont = 255;
00051     Double_t stops[NRGBs] = { 0.00, 0.15, 0.30, 0.45, 0.65, 0.85, 1.00 };
00052     Double_t red[NRGBs]   = { 0.60, 0.30, 0.00, 0.00, 0.60, 0.40, 0.00 };
00053     Double_t green[NRGBs] = { 1.00, 0.90, 0.80, 0.75, 0.20, 0.00, 0.00 };
00054     Double_t blue[NRGBs]  = { 1.00, 1.00, 1.00, 0.30, 0.00, 0.00, 0.00 };
00055     TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
00056     gStyle->SetNumberContours(NCont);
00057   }
00058 
00059   int ncol = gStyle->GetNumberOfColors();
00060   std::cout << "Number of colors "  << ncol << std::endl;
00061   
00062   if(ff->cd(module)) {
00063 
00064     TProfile* aveoccu= (TProfile*)gDirectory->Get("aveoccu");
00065     TProfile* avemult= (TProfile*)gDirectory->Get("avemult");
00066     TH1F* nchannels = (TH1F*)gDirectory->Get("nchannels_real");
00067 
00068     TProfile* averadius = (TProfile*)gDirectory->Get("averadius"); 
00069     TProfile* avez = (TProfile*)gDirectory->Get("avez"); 
00070 
00071     std::cout << "pointers " << aveoccu << " " << avemult << " " << nchannels << " " << averadius << " " << avez << std::endl;
00072 
00073     if(aveoccu && avemult && nchannels && averadius && avez) {
00074 
00075       nchannels->Sumw2();
00076       for(int i=1;i<nchannels->GetNbinsX()+1;++i) {
00077         nchannels->SetBinError(i,0.);
00078       }
00079 
00080       TH1D* haveoccu = aveoccu->ProjectionX("haveoccu");
00081       haveoccu->SetDirectory(0);
00082       haveoccu->Divide(nchannels);
00083       TH1D* havemult = avemult->ProjectionX("havemult");
00084       havemult->SetDirectory(0);
00085       havemult->Divide(nchannels);
00086 
00087       TH1D* havewidth = (TH1D*)haveoccu->Clone("havewidth");
00088       havewidth->SetDirectory(0);
00089       havewidth->SetTitle("Average Cluster Size");
00090       havewidth->Divide(havemult);
00091 
00092 
00093       new TCanvas("occupancy","occupancy",1200,500);
00094       gPad->SetLogy(1);
00095       haveoccu->SetStats(0);
00096       haveoccu->DrawCopy();
00097       TLine* l1 = new TLine(1000,0,1000,haveoccu->GetMaximum()); l1->DrawClone(); 
00098       TLine* l2 = new TLine(2000,0,2000,haveoccu->GetMaximum()); l2->DrawClone(); 
00099       TLine* l3 = new TLine(3000,0,3000,haveoccu->GetMaximum()); l3->DrawClone(); 
00100       TLine* l4 = new TLine(4000,0,4000,haveoccu->GetMaximum()); l4->DrawClone(); 
00101       TLine* l5 = new TLine(5000,0,5000,haveoccu->GetMaximum()); l5->DrawClone(); 
00102       TText* tpix = new TText(500,haveoccu->GetMaximum(),"BPIX+FPIX"); tpix->SetTextAlign(22); tpix->DrawClone();
00103       TText* ttib = new TText(1500,haveoccu->GetMaximum(),"TIB"); ttib->SetTextAlign(22); ttib->DrawClone();
00104       TText* ttid = new TText(2500,haveoccu->GetMaximum(),"TID"); ttid->SetTextAlign(22); ttid->DrawClone();
00105       TText* ttob = new TText(3500,haveoccu->GetMaximum(),"TOB"); ttob->SetTextAlign(22); ttob->DrawClone();
00106       TText* ttecm = new TText(4500,haveoccu->GetMaximum(),"TEC-"); ttecm->SetTextAlign(22); ttecm->DrawClone();
00107       TText* ttecp = new TText(5500,haveoccu->GetMaximum(),"TEC+"); ttecp->SetTextAlign(22); ttecp->DrawClone();
00108       
00109       new TCanvas("multiplicity","multiplicity",1200,500);
00110       gPad->SetLogy(1);
00111       havemult->SetStats(0);
00112       havemult->DrawCopy();
00113       tpix->SetY(havemult->GetMaximum()); tpix->DrawClone();
00114       ttib->SetY(havemult->GetMaximum()); ttib->DrawClone();
00115       ttid->SetY(havemult->GetMaximum()); ttid->DrawClone();
00116       ttob->SetY(havemult->GetMaximum()); ttob->DrawClone();
00117       ttecm->SetY(havemult->GetMaximum()); ttecm->DrawClone();
00118       ttecp->SetY(havemult->GetMaximum()); ttecp->DrawClone();
00119       l1->SetY2(havemult->GetMaximum()); l1->DrawClone(); 
00120       l2->SetY2(havemult->GetMaximum()); l2->DrawClone(); 
00121       l3->SetY2(havemult->GetMaximum()); l3->DrawClone(); 
00122       l4->SetY2(havemult->GetMaximum()); l4->DrawClone(); 
00123       l5->SetY2(havemult->GetMaximum()); l5->DrawClone(); 
00124       
00125       new TCanvas("width","width",1200,500);
00126       havewidth->SetStats(0);
00127       havewidth->DrawCopy();
00128       tpix->SetY(havewidth->GetMaximum()); tpix->DrawClone();
00129       ttib->SetY(havewidth->GetMaximum()); ttib->DrawClone();
00130       ttid->SetY(havewidth->GetMaximum()); ttid->DrawClone();
00131       ttob->SetY(havewidth->GetMaximum()); ttob->DrawClone();
00132       ttecm->SetY(havewidth->GetMaximum()); ttecm->DrawClone();
00133       ttecp->SetY(havewidth->GetMaximum()); ttecp->DrawClone();
00134       l1->SetY2(havewidth->GetMaximum()); l1->DrawClone(); 
00135       l2->SetY2(havewidth->GetMaximum()); l2->DrawClone(); 
00136       l3->SetY2(havewidth->GetMaximum()); l3->DrawClone(); 
00137       l4->SetY2(havewidth->GetMaximum()); l4->DrawClone(); 
00138       l5->SetY2(havewidth->GetMaximum()); l5->DrawClone();
00139       
00140       TCanvas * o2 = new TCanvas("occupancy2","occupancy2",1200,800);
00141       o2->Divide(3,2);
00142       o2->cd(1);
00143       haveoccu->SetAxisRange(100,270);
00144       haveoccu->DrawCopy();
00145       tpix->SetY(haveoccu->GetMaximum()); tpix->SetX(185); tpix->DrawClone();
00146       o2->cd(2);
00147       haveoccu->SetAxisRange(1050,1450);
00148       haveoccu->DrawCopy();
00149       ttib->SetY(haveoccu->GetMaximum()); ttib->SetX(1250); ttib->DrawClone();
00150       o2->cd(3);
00151       haveoccu->SetAxisRange(2070,2400);
00152       haveoccu->DrawCopy();
00153       ttid->SetY(haveoccu->GetMaximum()); ttid->SetX(2235); ttid->DrawClone();
00154       o2->cd(4);
00155       haveoccu->SetAxisRange(3000,3700);
00156       haveoccu->DrawCopy();
00157       ttob->SetY(haveoccu->GetMaximum()); ttob->SetX(3350); ttob->DrawClone();
00158       o2->cd(5);
00159       haveoccu->SetAxisRange(4000,4850);
00160       haveoccu->DrawCopy();
00161       ttecm->SetY(haveoccu->GetMaximum()); ttecm->SetX(4425); ttecm->DrawClone();
00162       o2->cd(6);
00163       haveoccu->SetAxisRange(5000,5850);
00164       haveoccu->DrawCopy();
00165       ttecp->SetY(haveoccu->GetMaximum()); ttecp->SetX(5425); ttecp->DrawClone();
00166 
00167       TCanvas * m2 = new TCanvas("multiplicity2","multiplicity2",1200,800);
00168       m2->Divide(3,2);
00169       m2->cd(1);
00170       havemult->SetAxisRange(100,270);
00171       havemult->DrawCopy();
00172       tpix->SetY(havemult->GetMaximum()); tpix->SetX(185); tpix->DrawClone();
00173       m2->cd(2);
00174       havemult->SetAxisRange(1050,1450);
00175       havemult->DrawCopy();
00176       ttib->SetY(havemult->GetMaximum()); ttib->SetX(1250); ttib->DrawClone();
00177       m2->cd(3);
00178       havemult->SetAxisRange(2070,2400);
00179       havemult->DrawCopy();
00180       ttid->SetY(havemult->GetMaximum()); ttid->SetX(2235); ttid->DrawClone();
00181       m2->cd(4);
00182       havemult->SetAxisRange(3000,3700);
00183       havemult->DrawCopy();
00184       ttob->SetY(havemult->GetMaximum()); ttob->SetX(3350); ttob->DrawClone();
00185       m2->cd(5);
00186       havemult->SetAxisRange(4000,4850);
00187       havemult->DrawCopy();
00188       ttecm->SetY(havemult->GetMaximum()); ttecm->SetX(4425); ttecm->DrawClone();
00189       m2->cd(6);
00190       havemult->SetAxisRange(5000,5850);
00191       havemult->DrawCopy();
00192       ttecp->SetY(havemult->GetMaximum()); ttecp->SetX(5425); ttecp->DrawClone();
00193 
00194       TCanvas * w2 = new TCanvas("width2","width2",1200,800);
00195       w2->Divide(3,2);
00196       w2->cd(1);
00197       havewidth->SetAxisRange(100,270);
00198       havewidth->DrawCopy();
00199       tpix->SetY(havewidth->GetMaximum()); tpix->SetX(185); tpix->DrawClone();
00200       w2->cd(2);
00201       havewidth->SetAxisRange(1050,1450);
00202       havewidth->DrawCopy();
00203       ttib->SetY(havewidth->GetMaximum()); ttib->SetX(1250); ttib->DrawClone();
00204       w2->cd(3);
00205       havewidth->SetAxisRange(2070,2400);
00206       havewidth->DrawCopy();
00207       ttid->SetY(havewidth->GetMaximum()); ttid->SetX(2235); ttid->DrawClone();
00208       w2->cd(4);
00209       havewidth->SetAxisRange(3000,3700);
00210       havewidth->DrawCopy();
00211       ttob->SetY(havewidth->GetMaximum()); ttob->SetX(3350); ttob->DrawClone();
00212       w2->cd(5);
00213       havewidth->SetAxisRange(4000,4850);
00214       havewidth->DrawCopy();
00215       ttecm->SetY(havewidth->GetMaximum()); ttecm->SetX(4425); ttecm->DrawClone();
00216       w2->cd(6);
00217       havewidth->SetAxisRange(5000,5850);
00218       havewidth->DrawCopy();
00219       ttecp->SetY(havewidth->GetMaximum()); ttecp->SetX(5425); ttecp->DrawClone();
00220 
00221       // Loop on bins and creation of boxes
00222 
00223       TList modulesoccu;
00224       TList modulesmult;
00225 
00226       for(int i=1;i<haveoccu->GetNbinsX();++i) {
00227 
00228         if(averadius->GetBinEntries(i)*avez->GetBinEntries(i)) {
00229 
00230           double dz = 2.;
00231           double dr = 1.;
00232           // determine module size
00233           
00234           if(i > 100 && i < 200) { dz=3.33;dr=0.4;}
00235 
00236           if(i > 200 && i < 1000 && ( i%10 == 1 || i%10 == 7)) { dz=0.8;dr=0.4;}
00237           if(i > 200 && i < 1000 && !( i%10 == 1 || i%10 == 7)) { dz=0.8;dr=0.8;}
00238 
00239           if(i > 1000 && i < 2000) { dz=5.948;dr=0.4;}
00240 
00241           if(i > 3000 && i < 4000) { dz=9.440;dr=0.4;}
00242 
00243           if(i > 2000 && i < 3000  && (i%1000)/100 == 1) { dz=0.8;dr=5.647;} 
00244           if(i > 2000 && i < 3000  && (i%1000)/100 == 2) { dz=0.8;dr=4.512;} 
00245           if(i > 2000 && i < 3000  && (i%1000)/100 == 3) { dz=0.8;dr=5.637;} 
00246 
00247           if(i > 4000 && i < 6000  && (i%1000)/100 == 1) { dz=0.8;dr=4.362;} 
00248           if(i > 4000 && i < 6000  && (i%1000)/100 == 2) { dz=0.8;dr=4.512;} 
00249           if(i > 4000 && i < 6000  && (i%1000)/100 == 3) { dz=0.8;dr=5.637;} 
00250           if(i > 4000 && i < 6000  && (i%1000)/100 == 4) { dz=0.8;dr=5.862;} 
00251           if(i > 4000 && i < 6000  && (i%1000)/100 == 5) { dz=0.8;dr=7.501;} 
00252           if(i > 4000 && i < 6000  && (i%1000)/100 == 6) { dz=0.8;dr=9.336;} 
00253           if(i > 4000 && i < 6000  && (i%1000)/100 == 7) { dz=0.8;dr=10.373;} 
00254         
00255           {  
00256             TBox* modoccu = new TBox(avez->GetBinContent(i)-dz,averadius->GetBinContent(i)-dr,avez->GetBinContent(i)+dz,averadius->GetBinContent(i)+dr);
00257             modoccu->SetFillStyle(1001);
00258             int icol=int(ncol*(log(haveoccu->GetBinContent(i))-log(min))/(log(max)-log(min)));
00259             if(icol < 0) icol=0;
00260             if(icol > (ncol-1)) icol=(ncol-1);
00261             std::cout << i << " " << icol << " " << haveoccu->GetBinContent(i) << std::endl; 
00262             modoccu->SetFillColor(gStyle->GetColorPalette(icol));
00263             modulesoccu.Add(modoccu);
00264           }
00265           {
00266             TBox* modmult = new TBox(avez->GetBinContent(i)-dz,averadius->GetBinContent(i)-dr,avez->GetBinContent(i)+dz,averadius->GetBinContent(i)+dr);
00267             modmult->SetFillStyle(1001);
00268             int icol=int(ncol*(log(havemult->GetBinContent(i))-log(mmin))/(log(mmax)-log(mmin)));
00269             if(icol < 0) icol=0;
00270             if(icol > (ncol-1)) icol=(ncol-1);
00271             std::cout << i << " " << icol << " " << havemult->GetBinContent(i) << std::endl; 
00272             modmult->SetFillColor(gStyle->GetColorPalette(icol));
00273             modulesmult.Add(modmult);
00274           }
00275 
00276         }
00277 
00278       }
00279       // eta boundaries lines
00280       TList etalines;
00281       TList etalabels;
00282       for(int i=0;i<8;++i) {
00283         double eta = 3.0-i*0.2;
00284         TLine* lin = new TLine(295,2*295/(exp(eta)-exp(-eta)),305,2*305/(exp(eta)-exp(-eta)));
00285         etalines.Add(lin);
00286         char lab[100];
00287         sprintf(lab,"%3.1f",eta);
00288         TText* label = new TText(285,2*285/(exp(eta)-exp(-eta)),lab);
00289         label->SetTextSize(.03);
00290         label->SetTextAlign(22);
00291         etalabels.Add(label);
00292       }
00293       for(int i=0;i<8;++i) {
00294         double eta = -3.0+i*0.2;
00295         TLine* lin = new TLine(-295,-2*295/(exp(eta)-exp(-eta)),-305,-2*305/(exp(eta)-exp(-eta)));
00296         etalines.Add(lin);
00297         char lab[100];
00298         sprintf(lab,"%3.1f",eta);
00299         TText* label = new TText(-285,-2*285/(exp(eta)-exp(-eta)),lab);
00300         label->SetTextSize(.03);
00301         label->SetTextAlign(22);
00302         etalabels.Add(label);
00303       }
00304       for(int i=0;i<15;++i) {
00305         double eta = -1.4+i*0.2;
00306         TLine* lin = new TLine(130.*(exp(eta)-exp(-eta))/2.,130,138.*(exp(eta)-exp(-eta))/2.,138);
00307         etalines.Add(lin);
00308         char lab[100];
00309         sprintf(lab,"%3.1f",eta);
00310         TText* label = new TText(125.*(exp(eta)-exp(-eta))/2.,125,lab);
00311         label->SetTextSize(.03);
00312         label->SetTextAlign(22);
00313         etalabels.Add(label);
00314       }
00315 
00316 
00317       TGaxis *raxis = new TGaxis(-310,0,-310,140,0,140,10,"S");
00318       TGaxis *zaxis = new TGaxis(-310,0,310,0,-310,310,10,"S");
00319       raxis->SetTickSize(.01);      zaxis->SetTickSize(.01);
00320       raxis->SetTitle("R (cm)"); zaxis->SetTitle("Z (cm)");
00321 
00322       TList palette;
00323       TList mpalette;
00324 
00325       for(int i = 0;i< ncol ; ++i) {
00326         TBox* box= new TBox(315,0+140./ncol*i,330,0+140./ncol*(i+1));
00327         box->SetFillStyle(1001);
00328         box->SetFillColor(gStyle->GetColorPalette(i));
00329         palette.Add(box);
00330         mpalette.Add(box);
00331 
00332       }
00333 
00334       TGaxis *paxis = new TGaxis(330,0,330,140,min,max,510,"SLG+");
00335       paxis->SetTickSize(.02);
00336       paxis->SetLabelOffset(paxis->GetLabelOffset()*0.5);
00337       palette.Add(paxis);
00338 
00339       TGaxis *mpaxis = new TGaxis(330,0,330,140,mmin,mmax,510,"SLG+");
00340       mpaxis->SetTickSize(.02);
00341       mpaxis->SetLabelOffset(paxis->GetLabelOffset()*0.5);
00342       mpalette.Add(mpaxis);
00343 
00344       TCanvas* cc1 = new TCanvas("occumap","occumap",1000,500);
00345       cc1->Range(-370.,-20.,390.,150.);
00346       TFrame* fr1 = new TFrame(-310,0,310,140);
00347       fr1->UseCurrentStyle();
00348       fr1->Draw();
00349       raxis->Draw(); zaxis->Draw();
00350       std::cout << modulesoccu.GetSize() << std::endl;
00351       etalines.Draw();
00352       etalabels.Draw();
00353       palette.Draw();
00354       modulesoccu.Draw();
00355 
00356       TCanvas* cc2 = new TCanvas("multmap","multmap",1000,500); 
00357       cc2->Range(-370.,-20.,390.,150.);
00358       TFrame* fr2 = new TFrame(-310,0,310,140);
00359       fr2->UseCurrentStyle();
00360       fr2->Draw();
00361       raxis->Draw(); zaxis->Draw();
00362       std::cout << modulesmult.GetSize() << std::endl;
00363       etalines.Draw();
00364       etalabels.Draw();
00365       mpalette.Draw();
00366       modulesmult.Draw();
00367 
00368     }
00369 
00370 
00371   }
00372 
00373 }