GANIL/PSAWatchers.C

Go to the documentation of this file.

/***************************************************************************
 *   Copyright (C) 2010 by Olivier Stezowski                               *
 *   stezow(AT)ipnl.in2p3.fr                                               *
 *                                                                         *
 *   This program is free software; you can redistribute it and/or modify  *
 *   it under the terms of the GNU General Public License as published by  *
 *   the Free Software Foundation; either version 2 of the License, or     *
 *   (at your option) any later version.                                   *
 *                                                                         *
 *   This program is distributed in the hope that it will be useful,       *
 *   but WITHOUT ANY WARRANTY; without even the implied warranty of        *
 *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         *
 *   GNU General Public License for more details.                          *
 *                                                                         *
 *   You should have received a copy of the GNU General Public License     *
 *   along with this program; if not, write to the                         *
 *   Free Software Foundation, Inc.,                                       *
 *   59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.             *
 ***************************************************************************/

#include "PSAWatchers.h"

#include "TH3.h"
#include "TH2.h"
#include "TH1.h"
#include "TCanvas.h"
#include "TPad.h"
#include "Hits.h"
#include "TPaveStats.h"
#include "TLegend.h"
#include "CrystalWatchers.h"
#include "TPaveText.h"
#include "TStyle.h"

#include "LoupeOnPad.h"

using namespace ADF;

Bool_t PSAWatcher::SetTrigger(ADF::DFTrigger* trigger ) 
{
        return Watcher::GetFromTrigger(trigger, "data:psa", fFrame);
}

Float_t PSACrystal3D::fgMetric = 10.;
Float_t PSACrystal3D::fRadius = 5. ;
Float_t PSACrystal3D::fDepht = 10.;
Int_t PSACrystal3D::fNBinX = 50;
Int_t PSACrystal3D::fNBinZ = 50;
Float_t PSACrystal3D::fDefaultECutMin = 0.;
Float_t PSACrystal3D::fDefaultECutMax = 100000000.;


PSACrystal3D::PSACrystal3D( const char * name, const char * title, TDirectory *sp_dir, TDirectory *tag_dir):
        PSAWatcher(name,title,sp_dir,tag_dir),
        fH3DHitPosition(0x0),
        fE0(0x0),       
        fE1(0x0),       
        fT0(0x0),       
        fT1(0x0),
        fTE0(0x0),
        fTE1(0x0),      
        fEnergySeg(0x0),        
        fNbHits(0x0),   
        fCurrentECutMin(fDefaultECutMin),
        fCurrentECutMax(fDefaultECutMax)       
{
        fH3DHitPosition = MakeTH3<TH3F>("H3DHitPosition", "Distribution of Hits in the crystal; x (cm); y (cm); z (cm)",
                                                           fNBinX, -fRadius, fRadius,
                                                           fNBinX, -fRadius, fRadius,
                                                           fNBinZ, 0, fDepht );
        fE0 = MakeTH1<TH1F>("ECore0", "Core energy 0, (high gain)",6000, 0, 6000);
    fE1 = MakeTH1<TH1F>("ECore1", "Core energy 1, (low gain) ",1000, 0, 10000);
    fE1->SetLineColor(kRed);
        fT0 = MakeTH1<TH1F>("TCore0", "Core time 0",300, 0, 100);
        fT1 = MakeTH1<TH1F>("TCore1", "Core time 1",300, 0, 100);
    fT1->SetLineColor(kRed);
    fTE0 = MakeTH2<TH2F>("TECore0", "Time energy 0",600,0,6000,300, 0, 100);
    fTE1 = MakeTH2<TH2F>("TECore1", "Time energy 1",600, 0,6000,300, 0, 100);

        fEnergySeg  = MakeTH1<TH1F>("ESum", "Sum of segment",6000, 0, 6000);
        fNbHits     = MakeTH1<TH1F>("NbHits", "Nb hits in crystal",20, 0, 20);
}

PSACrystal3D::~PSACrystal3D()
{
                ;
}


void PSACrystal3D::DoCanvas(TCanvas * c, Option_t *o)
{
        TString option = o;
    if ( option == "3D" ) {
                fH3DHitPosition->Draw(o);       
    }
    else {
        c->Divide(3,2);
        c->cd(1);
        fE0->Draw();
        c->cd(2);
        fE1->Draw();
        c->cd(3);
        fEnergySeg->Draw();
        c->cd(4);
        fT0->Draw();
        fT1->Draw("same");
        c->cd(5);
        fTE0->Draw("col");
        c->cd(6);
        fTE1->Draw("col");
    }
}

void PSACrystal3D::Exec(Option_t * /*option*/)
{
        // be sure the frame has been set properly
        if ( fFrame == 0x0 || !fFrame->IsValid() ) 
                { SetLastExecStatus(-1); return; }
        
        const PSAInterface *data = GetCstDataPointer<PSAInterface>(fFrame); 
        
        Double_t e = 0.;
        for( UShort_t i=0u; i < data->GetNbHits(); i++ ) {
                e +=  data->GetHit(i)->GetE();
                if( data->GetHit(i)->GetE() >= fCurrentECutMin && data->GetHit(i)->GetE() <= fCurrentECutMax )
                        fH3DHitPosition->Fill(data->GetHit(i)->GetX()/fgMetric,
                                                data->GetHit(i)->GetY()/fgMetric,
                                                data->GetHit(i)->GetZ()/fgMetric);
        }
        fNbHits->Fill(data->GetNbHits());
        fE0->Fill(data->GetE(0u));
        fE1->Fill(data->GetE(1u));
        fT0->Fill(data->GetT(0));
        fT1->Fill(data->GetT(1));
        fTE0->Fill(data->GetE(0u),data->GetT(0));
        fTE1->Fill(data->GetE(1u),data->GetT(1));
        fEnergySeg->Fill(e);
}

void PSACrystal3D::SetEnergyCut( Float_t EcutMin,Float_t EcutMax ){
        if( EcutMax < EcutMin ) {
                std::cerr << "ERROR: PSACrystal3D::SetEnergyCut( Float_t EcutMin,Float_t EcutMax)"
                << std::endl
                << "ECutMin greater than ECutMax, the cut has not been changed : "
                << "ECutMin = " << fCurrentECutMin
                << "  EcutMax = " << fCurrentECutMax << std::endl;
                return;
        }
        
        if( EcutMin != fCurrentECutMin && EcutMax != fCurrentECutMax ) {
                fCurrentECutMin = EcutMin;
                fCurrentECutMax = EcutMax;
                fH3DHitPosition->Reset();
                TString newtitle = Form("Distribution of Hits in the crystal E #in [%lf,%lf]; x (cm); y (cm); z (cm)",
                                                                EcutMin,
                                                                EcutMax);
                fH3DHitPosition->SetTitle(newtitle.Data());
        }
}

void PSACrystal3D::SetCurrentECutMin( Float_t EcutMin )
{
        fCurrentECutMin = EcutMin;
}
void PSACrystal3D::SetCurrentECutMax( Float_t EcutMax )
{
        fCurrentECutMax = EcutMax;
}
Float_t PSACrystal3D::GetCurrentECutMin() const
{
        return(fCurrentECutMin);
}
Float_t PSACrystal3D::GetCurrentECutMax() const
{
        return(fCurrentECutMax);
}



Float_t PSACrystal::fgMetric = 10.;
Float_t PSACrystal::fRadius = 5. ;
Float_t PSACrystal::fDepht = 10.;
Int_t PSACrystal::fNBinX = 50;
Int_t PSACrystal::fNBinZ = 50;
Float_t PSACrystal::fDefaultECutMin = 0.;
Float_t PSACrystal::fDefaultECutMax = 100000000.;

PSACrystal::PSACrystal( const char * name, const char * title, TDirectory *sp_dir, TDirectory *tag_dir):
    PSAWatcher(name,title,sp_dir,tag_dir),
    CanvasVisu(name,title),
    fH3DHitPosition(0x0),
    fE0(0x0),
    fE1(0x0),
    fT0(0x0),
    fT1(0x0),
    fTE0(0x0),
    fTE1(0x0),
    fEnergySeg(0x0),
    fNbHits(0x0),
    fCurrentECutMin(fDefaultECutMin),
    fCurrentECutMax(fDefaultECutMax),
    fSpectra( new TObjArray(CrystalInterface::kNbSegments) ),
    fNevts(0),
    fCrystalID(0)
{
    fH3DHitPosition = MakeTH3<TH3F>("H3DHitPosition", "Distribution of Hits in the crystal; x (cm); y (cm); z (cm)",
                                    fNBinX, -fRadius, fRadius,
                                    fNBinX, -fRadius, fRadius,
                                    fNBinZ, 0, fDepht );
    fE0 = MakeTH1<TH1F>("ECore_HG", "Core energy 0, (high gain)",6000, 0, 6000);
    fE1 = MakeTH1<TH1F>("ECore_LG", "Core energy 1, (low gain) ",6000, 0, 6000);
    for (Int_t i = 0; i < CrystalInterface::kNbSegments; i++ ) {
        PadPos HistPosition = fHistPositionList[i];
        TString Name = (TString)"PSA_" + "CrystalName" + "_Seg_" + HistPosition.SegName + "_Specta";
        TH1F *h = MakeTH1<TH1F>(Name.Data(),Name.Data(), 6000,0,6000);
        h->GetXaxis()->SetTitle("Energy (keV)");
        h->GetYaxis()->SetTitle("Counts");
        fSpectra->AddAt(h,i);
    }
    fT0 = MakeTH1<TH1F>("TCore_HG", "Core time 0 (high gain)",1000, 0, 100);
    fT1 = MakeTH1<TH1F>("TCore_LG", "Core time 1 (low gain)",1000, 0, 100);
    fTE0 = MakeTH2<TH2F>("TECore_HG", "Time energy 0 (high gain)",6000,0,6000,1000, 0, 100);
    fTE1 = MakeTH2<TH2F>("TECore_LG", "Time ebergy 1 (low gain)",6000, 0,6000,1000, 0, 100);

    fEnergySeg  = MakeTH1<TH1F>("ESum", "Sum energy of segments",6000, 0, 6000);
    fNbHits     = MakeTH1<TH1F>("NbHits", "Nb hits in crystal",20, 0, 20);
    
    for(int i=0 ; i<6 ; i++)
    {
        fHitsSlices[i] = MakeTH2<TH2F>(Form("2DHits_%d",i), Form("2DHits_%d",i),fNBinX,-fRadius,fRadius,fNBinX, -fRadius, fRadius);
        fHitsSlices[i]->GetXaxis()->SetTitle("X (cm)");
        fHitsSlices[i]->GetYaxis()->SetTitle("Y (cm)");
    }
}

PSACrystal::~PSACrystal()
{

}

Bool_t PSACrystal::SetTrigger(ADF::DFTrigger* trigger )
{
    return Watcher::GetFromTrigger(trigger, "data:psa", fFrame);
}

void PSACrystal::Exec(Option_t * /*option*/)
{
    // be sure the frame has been set properly
    if ( fFrame == 0x0 || !fFrame->IsValid() )
    { SetLastExecStatus(-1); return; }

    fNevts++;

    const PSAInterface *data = GetCstDataPointer<PSAInterface>(fFrame);

    if(fCrystalID == 0) fCrystalID = data->GetUID();

    Double_t SumSeg = 0.;
    int SegID=0;
    double ESeg=0.;

    for( UShort_t i=0u; i < data->GetNbHits(); i++ )
    {
        SumSeg +=  data->GetHit(i)->GetE();
        ESeg = data->GetHit(i)->GetE();
        SegID = data->GetHit(i)->GetID(0);

        ((TH1*)(fSpectra->At(SegID)))->Fill(ESeg);

        if( data->GetHit(i)->GetE() >= fCurrentECutMin && data->GetHit(i)->GetE() <= fCurrentECutMax )
        {
            fH3DHitPosition->Fill(data->GetHit(i)->GetX()/fgMetric, data->GetHit(i)->GetY()/fgMetric, data->GetHit(i)->GetZ()/fgMetric);

            int slice = data->GetHit(i)->GetID(0)%6;

            fHitsSlices[slice]->Fill(data->GetHit(i)->GetX()/fgMetric, data->GetHit(i)->GetY()/fgMetric);
        }
        
    }
    fNbHits->Fill(data->GetNbHits());
    fE0->Fill(data->GetE(0u));
    fE1->Fill(data->GetE(1u));
    fT0->Fill(data->GetT(0));
    fT1->Fill(data->GetT(1));
    fTE0->Fill(data->GetE(0u),data->GetT(0));
    fTE1->Fill(data->GetE(1u),data->GetT(1));
    fEnergySeg->Fill(SumSeg);
}

void PSACrystal::SetEnergyCut( Float_t EcutMin,Float_t EcutMax ){
    if( EcutMax < EcutMin ) {
        std::cerr << "ERROR: PSACrystal3D::SetEnergyCut( Float_t EcutMin,Float_t EcutMax)"
                  << std::endl
                  << "ECutMin greater than ECutMax, the cut has not been changed : "
                  << "ECutMin = " << fCurrentECutMin
                  << "  EcutMax = " << fCurrentECutMax << std::endl;
        return;
    }

    if( EcutMin != fCurrentECutMin && EcutMax != fCurrentECutMax ) {
        fCurrentECutMin = EcutMin;
        fCurrentECutMax = EcutMax;
        fH3DHitPosition->Reset();
        TString newtitle = Form("Distribution of Hits in the crystal E #in [%lf,%lf]; x (cm); y (cm); z (cm)",
                                EcutMin,
                                EcutMax);
        fH3DHitPosition->SetTitle(newtitle.Data());
    }
}

void PSACrystal::SetCurrentECutMin( Float_t EcutMin )
{
    fCurrentECutMin = EcutMin;
}
void PSACrystal::SetCurrentECutMax( Float_t EcutMax )
{
    fCurrentECutMax = EcutMax;
}
Float_t PSACrystal::GetCurrentECutMin() const
{
    return(fCurrentECutMin);
}
Float_t PSACrystal::GetCurrentECutMax() const
{
    return(fCurrentECutMax);
}

void PSACrystal::ShowSpectraMap(const char *option)
{
    TString opt = (TString)option;
    TCanvas *c = 0x0;

    const PSAInterface *data = GetCstDataPointer<PSAInterface>(fFrame);
    TString CrystalName = GetCrystalNameFromID(data->GetUID());

    opt.ToLower();

    if(opt == "same") c = (TCanvas*)gPad->GetCanvas();
    else
    {
        TString Name = "PSA_Spectra_" + CrystalName;
        TString Title = "Mapping of the spectra for crystal " + CrystalName + "at the post psa level";

        c = CanvasVisu::BuildSegmentedCanvas(Name,Title);
    }
    if(c == 0x0) return;

    c->cd();

    for(int i=0 ; i<CrystalInterface::kNbSegments ; i++)
    {
        PadPos HistPosition = fHistPositionList[i];

        TPad *pad = (TPad*)c->cd(HistPosition.PadNmbr);

        TH1 *h = 0x0;

        TPaveText *p = (TPaveText*)pad->GetListOfPrimitives()->FindObject("InfoPad");

        h = (TH1*)fSpectra->At(i);

        TString NameTitle = (TString)"Psa_" + CrystalName + "_Seg_" + HistPosition.SegName + "_Specta";

        h->SetNameTitle(NameTitle,NameTitle);

        int col = 0;
        for(int k=0 ; k<pad->GetListOfPrimitives()->GetSize() ; k++)
        {
            if(pad->GetListOfPrimitives()->At(k)->InheritsFrom("TH1"))col++;
        }
        h->SetLineColor(GetColor(col));
        h->SetFillColor(0);
        h->SetStats(0);
        h->Draw(opt);
        if(col==0) h->GetXaxis()->SetRangeUser(1,h->GetXaxis()->GetXmax());

        if(p)p->Draw();

        pad->Modified();
        pad->Update();
    }
    for(int i=0 ; i<CrystalInterface::kNbCores ; i++)
    {
        PadPos HistPosition = fHistPositionList[CrystalInterface::kNbSegments+i];

        TPad *pad = (TPad*)c->cd(HistPosition.PadNmbr);

        TH1 *h = 0x0;

        TPaveText *p = (TPaveText*)pad->GetListOfPrimitives()->FindObject("InfoPad");

        if(i==0) h = fE0;
        else h = fE1;

        TString NameTitle = (TString)"Psa_" + CrystalName + "_Seg_" + HistPosition.SegName + "_Specta";

        h->SetNameTitle(NameTitle,NameTitle);

        int col = 0;
        for(int k=0 ; k<pad->GetListOfPrimitives()->GetSize() ; k++)
        {
            if(pad->GetListOfPrimitives()->At(k)->InheritsFrom("TH1"))col++;
        }
        h->SetLineColor(GetColor(col));
        h->SetFillColor(0);
        h->SetStats(0);
        h->Draw(opt);
        if(col==0) h->GetXaxis()->SetRangeUser(1,h->GetXaxis()->GetXmax());

        if(p)p->Draw();

        pad->Modified();
        pad->Update();
    }

    if(!opt.Contains("same"))
    {
        new LoupeOnPad(c);
        SetLoupe(false);
    }
}

void PSACrystal::Show3DHits()
{
    const PSAInterface *data = GetCstDataPointer<PSAInterface>(fFrame);
    TString CrystalName = GetCrystalNameFromID(data->GetUID());

    TString Name = "PSA_3DHits" + CrystalName;
    TString Title = "PSA 3D Hits for crystal " + CrystalName;

    TCanvas *canvas = CanvasVisu::NewCanvas(Name,Title);

    gStyle->SetOptTitle(0);

    canvas->Divide(3,1,0.001,0.001);
    canvas->cd(1)->Divide(1,3,0.001,0.001);
    for(int i=0 ; i<3 ; i++)
    {
        canvas->cd(1)->cd(3-i);
        fHitsSlices[i]->Draw("col");
    }

    canvas->cd(2)->Divide(1,3,0.001,0.001);
    for(int i=3 ; i<6 ; i++)
    {
        canvas->cd(2)->cd(3-(i-3));
        fHitsSlices[i]->Draw("col");
    }

    canvas->cd(3);

    fH3DHitPosition->Draw("");

}

void PSACrystal::ShowCoreVsSumSeg()
{
    const PSAInterface *data = GetCstDataPointer<PSAInterface>(fFrame);
    TString CrystalName = GetCrystalNameFromID(data->GetUID());

    TString Name = "PSA_CoreVsSumSeg" + CrystalName;
    TString Title = "PSA Core versus sum segments spectra " + CrystalName;

    TCanvas *c = CanvasVisu::NewCanvas(Name,Title);

    c->Divide(1,2,0.001,0.001);

    double topmarg = 0.0657599;
    double rightmarg = 0.00726832;
    double leftmarg = 0.0490612;
    double bottommarg = 0.114119;

    TH1 *h[2] = {fNbHits,fEnergySeg};
    TString titles[2] = {"Fold","Energy (keV)"};

    for(int i=0 ; i<2 ; i++)
    {
        c->cd(i+1);

        gPad->SetRightMargin(rightmarg);
        gPad->SetTopMargin(topmarg);
        gPad->SetLeftMargin(leftmarg);
        gPad->SetBottomMargin(bottommarg);

        TH1 *htemp = h[i];

        htemp->SetLineColor(GetColor(0));
        htemp->Draw();

        htemp->GetXaxis()->SetTitle(titles[i]);
        htemp->GetXaxis()->SetTitleSize(0.06);
        htemp->GetXaxis()->SetLabelSize(0.06);
        htemp->GetXaxis()->SetTitleOffset(0.91);

        htemp->GetYaxis()->SetTitle("Counts");
        htemp->GetYaxis()->SetTitleSize(0.06);
        htemp->GetYaxis()->SetLabelSize(0.06);
        htemp->GetYaxis()->SetTitleOffset(0.45);

        if(i==1)
        {
            htemp->SetBit(TH1::kNoTitle);
            fE0->SetLineColor(GetColor(1));
            fE0->Draw("same");
            fE1->SetLineColor(GetColor(2));
            fE1->Draw("same");
        }

        gPad->Modified();
        gPad->Update();

        TPaveStats *st = (TPaveStats*)htemp->FindObject("stats");
        if(st)
        {
            st->SetX2NDC(1-rightmarg);
            st->SetY2NDC(1-topmarg);
        }

        TPaveText *txt = (TPaveText*)gPad->FindObject("title");
        if(txt)
        {
            double width = txt->GetY2NDC()-txt->GetY1NDC();
            txt->SetY2NDC(1-topmarg);
            txt->SetY1NDC(txt->GetY2NDC()-width);
        }

        if(i==1)
        {
            TLegend *leg = new TLegend(st->GetX1NDC(),st->GetY1NDC()-(st->GetY2NDC()-st->GetY1NDC()),st->GetX2NDC(),st->GetY1NDC(),"","NDC");
            leg->AddEntry(fEnergySeg,fEnergySeg->GetTitle(),"l");
            leg->AddEntry(fE0,fE0->GetTitle(),"l");
            leg->AddEntry(fE1,fE1->GetTitle(),"l");
            leg->Draw();
        }
    }
}

void PSACrystal::ShowTimeVsE()
{
    const PSAInterface *data = GetCstDataPointer<PSAInterface>(fFrame);
    TString CrystalName = GetCrystalNameFromID(data->GetUID());

    TString Name = "PSA_TimeVsE" + CrystalName;
    TString Title = "PSA time versus energy of the core" + CrystalName;

    TCanvas *c = CanvasVisu::NewCanvas(Name,Title);

    c->Divide(2,2,0.001,0.001);

    double topmarg = 0.0560089;
    double rightmarg = 0.0234869;
    double leftmarg = 0.0747608;
    double bottommarg = 0.104442;

    TH1 *h[4] = {fT0,fT1,fTE0,fTE1};

    for(int i=0 ; i<4 ; i++)
    {
        c->cd(i+1);

        gPad->SetRightMargin(rightmarg);
        gPad->SetTopMargin(topmarg);
        gPad->SetLeftMargin(leftmarg);
        gPad->SetBottomMargin(bottommarg);

        TH1 *htemp = h[i];

        htemp->SetLineColor(GetColor(0));

        if(i<2)htemp->Draw();
        else
        {
            htemp->Draw("col");
            gPad->SetLogz();
        }

        if(i<2) htemp->GetXaxis()->SetTitle("Time (ns)");
        else htemp->GetXaxis()->SetTitle("Energy (keV)");
        htemp->GetXaxis()->SetTitleSize(0.05);
        htemp->GetXaxis()->SetLabelSize(0.05);


        if(i<2) htemp->GetYaxis()->SetTitle("Counts");
        else htemp->GetYaxis()->SetTitle("Time (ns)");
        htemp->GetYaxis()->SetTitleSize(0.05);
        htemp->GetYaxis()->SetLabelSize(0.05);
        htemp->GetYaxis()->SetTitleOffset(0.8);

        htemp->SetBit(TH1::kNoTitle);

        gPad->Modified();
        gPad->Update();

        TPaveStats *st = (TPaveStats*)htemp->FindObject("stats");
        if(st)
        {
            st->SetX2NDC(1-rightmarg);
            st->SetY2NDC(1-topmarg);
        }
    }
}

void PSACrystal::Print(const char * /*option*/ ) const
{
    std::cout<<std::endl;
    int size = 3 + 3 + ((TString)GetName()).Length() + 5 + ((TString)GetTitle()).Length() + 3 + 3;
    for(int i=0; i<size ; i++) std::cout<<"*";
    std::cout<<std::endl;
    std::cout<<"***   "<<GetName()<<" --- "<<GetTitle()<<"   ***"<<std::endl;
    for(int i=0; i<size ; i++) std::cout<<"*";
    std::cout<<std::endl;

    std::cout<<"Number of treated events : "<<fNevts<<std::endl;
    std::cout<<"Crystal ID               : "<<fCrystalID<<std::endl;
}


ClassImp(PSAWatcher);
ClassImp(PSACrystal3D);
ClassImp(PSACrystal);