LLaBr3Watchers.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 "LLaBr3Watchers.h"
#include "TMath.h"
#include <fstream>
#include "MetaWatchers.h"
#include "TrackedWatchers.h"

Short_t BaseLLaBr3Watcher::gLUT[BaseLLaBr3Watcher::gNumberOfItems][BaseLLaBr3Watcher::gMaxNumberOfModules];

ClassImp(BaseLLaBr3Watcher);

void MyLLaBr3::Calibrate () 
{
        Double_t tmp = E;

        // not very useful .... keep it in case as an exemple 
        if ( E < threshold_adc_min || E > threshold_adc_max ) {
                has_fired = false;
        }
        else 
                E = A + B*tmp + C*TMath::Power(tmp,2) ;
}

BaseLLaBr3Watcher::BaseLLaBr3Watcher(const char *name, const char *title, TDirectory *sp_dir, TDirectory *tag_dir): 
        RancLegnaroWatcher(name,title,sp_dir,tag_dir)
{
        // to be sure the LUT is initialized
        InitLUT();
        fNumberOfModules = 0;

        if ( !LoadCalibration("LLaBr3.conf") )
                fLog << warning << "LLaBr3.conf not found in the current directory ... LoadCalibration should be called at least once, please use LoadCalibration !! " << dolog;
} 

BaseLLaBr3Watcher::~BaseLLaBr3Watcher()
{       
}

// should add nmore test to be sure the file is consistant ...
Bool_t BaseLLaBr3Watcher::LoadCalibration(const char *calib_file)
{
        std::ifstream in(calib_file); std::string tmp;
        
        if ( in.is_open() == false ) 
                return false;
        
        getline(in,tmp);
        tmp += " ";
        while ( in.good() && !in.eof() ) { // read input stream line by line

                if ( tmp[0] == '#' ) {
                        getline(in,tmp); 
                        continue; 
                } // this line is a comment 
                {
                        Float_t r, theta, phy, a, b, c; Int_t tmin, tmax;
                        Int_t which; 
                        std::istringstream decode(tmp);
                        
                        decode >> which;
                        if ( which == -1 ) { // this is LUT
                                for (Short_t i = 0; i < fNumberOfModules; i++ )  
                                        for (Short_t j = 0; j < gNumberOfItems; j++ )  
                                                decode >> gLUT[j][i] ;

/*                              cout << " LUT " << endl;
                                for (Short_t i = 0; i < fNumberOfModules; i++ )  
                                        for (Short_t j = 0; j < gNumberOfItems; j++ )  
                                                cout << j << " " <<  i << " " <<  gLUT[j][i] << endl ; */
                        }
                        else { // this a module
                          if ( which < gMaxNumberOfModules ) {

                                  std::cout << " Reading LaBr3 Module # " << fNumberOfModules << std::endl;
                                decode >> r >> theta >> phy >> tmin >> tmax >> a >> b >> c ;

                                fNumberOfModules++;
                                
                                fLaBr3[which].P.SetMagThetaPhi(r,theta*TMath::Pi()/180.0,phy*TMath::Pi()/180.0);

                                fLaBr3[which].threshold_adc_min = tmin;
                                fLaBr3[which].threshold_adc_max = tmax;

                                fLaBr3[which].A = a;
                                fLaBr3[which].B = b;
                                fLaBr3[which].C = c;
        
                           }
                           else {;}
                         }
                }
                // add null character at the end since getline
                // does not store \n and it may give troubles when decoding
                getline(in,tmp);
                tmp += " ";
        }       
        in.close();
  
        fLog << " LLaBr3 has been configured with " << fNumberOfModules << " modules " << dolog;

        return true;    
}

void BaseLLaBr3Watcher::InitLUT()
{
        // module 1
        gLUT[0][0] = 96;
        gLUT[1][0] = 97;

        // module 2
        gLUT[0][1] = 98;
        gLUT[1][1] = 99;

        // module 3
        gLUT[0][2] = 100;
        gLUT[1][2] = 101;

        // module 2
        gLUT[0][3] = 0;
        gLUT[1][3] = 0;

        // module 3
        gLUT[0][4] = 0;
        gLUT[1][4] = 0;
}

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

        // fill the array of floats
        RancLegnaroWatcher::Exec(option);
        
        // fill up the different modules : calibrates and check the threshold
        // to know the number of detector fired in this event
        for (Int_t i = 0u; i < fNumberOfModules; i++) {
                
                // fill the data part
                fLaBr3[i].E = rawBuf[gLUT[0][i]];
                fLaBr3[i].T = rawBuf[gLUT[1][i]];

                // if out of range, 
                fLaBr3[i].has_fired = true; fLaBr3[i].Calibrate(); 
        }
}


ClassImp(DisplayLLaBr3);

DisplayLLaBr3::DisplayLLaBr3(const char *name, const char *title, TDirectory *sp_dir, TDirectory *tag_dir): 
        BaseLLaBr3Watcher(name,title,sp_dir,tag_dir),
        fVertexBuilder(0x0),
        fCounts(0x0),
        fHitPattern(0x0),
        fE_no_doppler(gMaxNumberOfModules),
        fE_doppler(gMaxNumberOfModules),
        fT(gMaxNumberOfModules)
{
        // Vertex used for doppler correction, last registered
        fVertexBuilder = VertexBuilder::theCurrentVertexBuilder();
        
        fCounts = MakeTH1<TH1F>("NbModules","To get the number of modules fired",gMaxNumberOfModules,0,gMaxNumberOfModules);
        fHitPattern = MakeTH1<TH1F>("HitPattern","Hit pattern",gMaxNumberOfModules,0,gMaxNumberOfModules);

        for (Short_t i = 0u; i < fNumberOfModules; i++) {
                
                TString tmp; TH1F *h;

                tmp = Form("LaBr3_E_%d",i);
                h = MakeTH1<TH1F>(tmp.Data(),"Calibrated energy",1000,0,6000);
                fE_no_doppler.Add ( h ) ;

                tmp = Form("LaBr3_E_Doppler_%d",i);
                h = MakeTH1<TH1F>(tmp.Data(),"Doppler corrected energy",1000,0,6000);
                fE_doppler.Add ( h ) ;

                tmp = Form("LaBr3_T_%d",i);
                h = MakeTH1<TH1F>(tmp.Data(),"Time",1000,0,6000);
                fT.Add ( h ) ;
        }
} 

DisplayLLaBr3::~DisplayLLaBr3()
{ 

}

void DisplayLLaBr3::DoCanvas(TCanvas *c, Option_t * /*o*/)
{
        c->Divide(4, 2, 0.001, 0.001); // 
                
        for (Int_t i = 0u; i < fNumberOfModules; i++) {
                c->cd(i+1);
                TH1 *h = (TH1 *)fE_no_doppler.At(i); h->Draw();
        }
}

void DisplayLLaBr3::Exec(Option_t *option) 
{       
        // fill the array of floats
        BaseLLaBr3Watcher::Exec(option);
        
        // count the nummber of detectors that have fired
        UShort_t nb_fired = 0u;

        for (Int_t i = 0u; i < fNumberOfModules; i++) {

                if ( fLaBr3[i].has_fired ) {
                          fHitPattern->Fill(i);
                          nb_fired++;
                }
                else { continue; }
                //
                TH1 *h;
                
                // calibrated
                h = (TH1 *)fE_no_doppler.At(i);
                h->Fill(fLaBr3[i].E);
                h = (TH1 *)fT.At(i);
                h->Fill(fLaBr3[i].T);

                //cout << fLaBr3[i].Xr << " " << fLaBr3[i].Yr <<  " " << fLaBr3[i].Zr << endl;

                if ( fVertexBuilder->GetVertex()->GetBeta() > 0.01 )  {
                        h = (TH1 *)fE_doppler.At(i);
                        h->Fill(fLaBr3[i].E);
                }

        }
        fCounts->Fill(nb_fired);
}

ClassImp(MyLLaBr3Tree);

MyLLaBr3Tree::MyLLaBr3Tree(const char *name, const char *title, TTree *tree) : 
        BaseLLaBr3Watcher(name,title), 
        TTreeBuilder(name,title,tree)
{       
        if ( GetTree() ) { // just add the branch for that watcher
                SetBranches();
        }
}

void MyLLaBr3Tree::SetBranches()
{
        GetTree()->Branch("NbLaBr", &fNbLaBr3, "NbLaBr/I");
        GetTree()->Branch("WhichLaBr3", fWhichLaBr3, "WhichLaBr3[NbLaBr]/s");

        GetTree()->Branch("ELaBr3", fELaBr3, "E[NbLaBr]/F");
        GetTree()->Branch("TLaBr3", fTLaBr3, "T[NbLaBr]/F");
}

MyLLaBr3Tree::~MyLLaBr3Tree()
{
}

void MyLLaBr3Tree::Exec(Option_t *option) 
{
        if ( gDebug > 3 ) 
                printf("MyLLaBr3Tree::Exec is called \n");
        
        // fill the array of floats and compute the PRISMA physics
        BaseLLaBr3Watcher::Exec(option);

        fNbLaBr3 = 0;
        for (Int_t i = 0u; i < fNumberOfModules; i++) {
/*              if ( !fLaBr3[i].has_fired ) {
                          continue;
                } */

                fWhichLaBr3[fNbLaBr3] = i;
                fELaBr3[fNbLaBr3] = fLaBr3[i].E;
                fTLaBr3[fNbLaBr3] = fLaBr3[i].T;

                fNbLaBr3++;
        }

        // fill the tree only if this is the owner of the tree ... otherwise by TreeMaster
        FillTree(); 
}