WireGain/wiregVsWire_bylayer.cxx

Go to the documentation of this file.

//This code can display wiregains from two files on one pad. One pad for one layer. 
 
//#include "filename.h"
#include "TStyle.h"
#include "TString.h"
 
int wirebind[44] = {0,40,84,132,188,252,324,404,484,560,636,724,812,912,1012,1124,1236,1364,1492,1632,1772,1932,2092,2252,2412,2588,2764,2940,3116,3324,3532,3740,3948,4188,4428,4668,4908,5164,5420,5676,5932,6220,6508,6796};
const int LayNoi = 0;
const int LayNo = 43;
const int WireNo = 6796;
double wiregmin = 0.;
double wiregmax = 3;   //select good hits which are put in to calculate norm and to be drawn
TString LegendName = "wiregen vs wire"; 
const int Ipadx = 1;
const int Ipady = 1;
const int CanvasX = 800;  //canvas size for ps file
const int CanvasY = 600;
bool NormFlag = true;  //true: nomalize wiregain for each layer
const double histmin = 0.7; //Y axis range of histogram when NormFlag is true
const double histmax = 1.3;
const double delta = 0.15; //Y axis range of histogram when NormFlag is false
 
using namespace std;
 
vector<string>file;
void pushNo()
{
    file.push_back("wiregain.root");
    file.push_back("wiregain_new.root");
}
 
void wiregVsWire_bylayer()
{
    pushNo();
    const int FileNo = file.size();
    TH2F*** h = new TH2F**[LayNo];
    stringstream histname;
    TCanvas*** c = new TCanvas**[LayNo];
    TLegend*** lgd = new TLegend**[LayNo];
    gStyle->SetOptStat(0);
    gStyle->SetCanvasColor(10);
    for(int i=0;i<LayNo;i++)
    {
        h[i] = new TH2F*[FileNo/2];
        c[i] = new TCanvas*[FileNo/2];
        lgd[i] = new TLegend*[FileNo/2];
    }
 
 
    TFile* f1;
    TFile* f2;
    TTree* tree1;
    TTree* tree2;
    double wireg1[WireNo],wireg2[WireNo],Norm1,Norm2;
 
    stringstream ss_file1;  //file name 
    stringstream ss_file2;
    stringstream ss1;  //picture's drawing content
    stringstream ss2;
    stringstream cut;
    stringstream legendheader; 
 
 
    for(int k=0;k<FileNo;k+=2)
    {
        ss_file1.str("");
        ss_file1<<file[k];
        //cout<<ss_file1.str().c_str()<<endl;
        f1=new TFile(ss_file1.str().c_str());
        tree1=(TTree*)f1->Get("wiregcalib");
        tree1->SetBranchAddress("wireg", wireg1);
        //tree1->SetBranchAddress("mean", wireg1);
        tree1->GetEntry(0);
        tree1->SetMarkerStyle(22);
        tree1->SetMarkerColor(4);
        tree1->SetLineColor(4);
 
        ss_file2.str("");
        ss_file2<<file[k+1];
        //cout<<ss_file2.str().c_str()<<endl;
        f2=new TFile(ss_file2.str().c_str());
        tree2=(TTree*)f2->Get("wiregcalib");
        tree2->SetBranchAddress("wireg", wireg2);
        //tree2->SetBranchAddress("mean", wireg2);
        tree2->GetEntry(0);
        tree2->SetMarkerStyle(20);
        tree2->SetMarkerColor(2);
        tree2->SetLineColor(2);
 
 
        for(int i=LayNoi;i<LayNo;i++)
        {
            cout<<"layer= "<<i<<endl;
            Norm1=0;
            Norm2=0;
            int k1=0;
            int k2=0;
            for(int j=wirebind[i];j<wirebind[i+1];j++)
            {
                //cout<<"wireg1["<<j<<"]= "<<wireg1[j]<<endl;
                if(wireg1[j]>wiregmin && wireg1[j]<wiregmax) {Norm1 += wireg1[j];k1++;}
                if(wireg2[j]>wiregmin && wireg2[j]<wiregmax) {Norm2 += wireg2[j];k2++;}
                //if(NormFlag==true) {Norm1 += wireg1[j];k1++;}
                //if(NormFlag==true) {Norm2 += wireg2[j];k2++;}
            }
            Norm1/= k1;
        Norm2/= k2;
            cout<<"Norm1= "<<Norm1<<" Norm2= "<<Norm2<<endl;
 
            legendheader.str("");
            legendheader<<"lay"<<i<<": "<<LegendName;
            lgd[i][k/2] = new TLegend(.3, .75, .89, .89);
            lgd[i][k/2] ->SetFillColor(10);
            lgd[i][k/2] ->SetBorderSize(1);
            lgd[i][k/2] ->SetHeader(legendheader.str().c_str());
            histname.str("");
            histname<<"lay"<<i<<"_"<<file[k]<<"-"<<file[k+1];
            if(NormFlag) h[i][k/2] = new TH2F(histname.str().c_str(),"",100,wirebind[i],wirebind[i+1],100,histmin,histmax);
            else h[i][k/2] = new TH2F(histname.str().c_str(),"",100,wirebind[i],wirebind[i+1],100,(Norm1+Norm2)/2-delta,(Norm1+Norm2)/2+delta);
            h[i][k/2] -> GetYaxis()->SetTitle("wiregain");
            h[i][k/2] -> GetYaxis()->SetTitleOffset(1.);
            h[i][k/2] -> GetXaxis()->SetTitle("wire");
            c[i][k/2] = new TCanvas(histname.str().c_str(),histname.str().c_str(),800,600);
            h[i][k/2] -> Draw();
 
            cut.str("");
            //cut<<"wire>="<<wirebind[i]<<" && wire<"<<wirebind[i+1]<<" && wireg>"<<wiregmin<<" && wireg<"<<wiregmax;
            cut<<"wire>="<<wirebind[i]<<" && wire<"<<wirebind[i+1]<<" && wireg/"<<Norm1<<" >"<<histmin<<" && wireg/"<<Norm1<<" <"<<histmax;
            cout<<"cut: "<<cut.str().c_str()<<endl;
 
            ss1.str("");
            //if(NormFlag) ss1<<"wireg/"<<Norm1<<":wire";
            if(NormFlag) ss1<<"wireg/"<<Norm1<<":wire";
            else ss1<<"wireg"<<":wire";
            tree1->Draw(ss1.str().c_str(),cut.str().c_str(),"sameLP");
            lgd[i][k/2] ->AddEntry(tree1, ss_file1.str().c_str() , "p");
            //f1->Close();
 
            cut.str("");
            cut<<"wire>="<<wirebind[i]<<" && wire<"<<wirebind[i+1]<<" && wireg/"<<Norm2<<" >"<<histmin<<" && wireg/"<<Norm2<<" <"<<histmax;
            ss2.str("");
            //if(NormFlag) ss2<<"wireg/"<<Norm2<<":wire";
            if(NormFlag) ss2<<"wireg/"<<Norm2<<":wire";
            else ss2<<"wireg"<<":wire";
            tree2->Draw(ss2.str().c_str(),cut.str().c_str(),"sameLP");
            lgd[i][k/2] ->AddEntry(tree2, ss_file2.str().c_str() , "p");
            //f2->Close();
            lgd[i][k/2] ->Draw();
        }
    }
 
    TFile* g = new TFile("diff.root","RECREATE");
    for(int i=LayNoi;i<LayNo;i++)
    {
        for(int k=0;k<FileNo;k+=2)
        {
            c[i][k/2] -> Write();  
        }
    }
    g->Close();
 
 
    TCanvas* c1 = new TCanvas("c1", "canvas", CanvasX, CanvasY);
    c1->Divide(Ipadx,Ipady);
    c1->Print("diff.ps[");
    for(int i=LayNoi;i<LayNo;i++)
    {
        for(int k=0;k<FileNo;k+=2)
        {
            c1->cd(k/2+1);
            gPad->Clear();
            c[i][k/2]->DrawClonePad();
            c[i][k/2] -> Close();
            lgd[i][k/2] ->Delete();
        }
        c1->Update();
        c1->Print("diff.ps");
    }
    c1->Print("diff.ps]");
    c1->Close();
 
    delete c1;
    for(int i=LayNoi;i<LayNo;i++)
    {   
        for(int k=0;k<FileNo;k+=2)
        {
            delete h[i][k/2];
            delete c[i][k/2];
        }
    }   
 
}