BRAT 2.4.5
Class index
Full class index
brahmlib
BRAHMS
ROOT page 
//____________________________________________________________________
// 
// Class BrTofAdcGainCalModule
//
// Calibration of adc gains
// 
//____________________________________________________________________
//
// $Id: BrTofAdcGainCalModule.cxx,v 1.11 2002/09/10 16:19:37 videbaek Exp $
// $Author: videbaek $
// $Date: 2002/09/10 16:19:37 $
// $Copyright: (C) 2001 BRAHMS Collaboration <brahmlib@rhic.bnl.gov>
//

#ifndef WIN32
#include <iostream>
#include <iomanip>
#include <fstream>
#else
#include <iostream.h>
#include <iomanip.h>
#include <fstream.h>
#endif
#ifndef BRAT_BrTofAdcGainCalModule
#include "BrTofAdcGainCalModule.h"
#endif
#ifndef ROOT_TF1
#include "TF1.h"
#endif
#ifndef ROOT_TH2
#include "TH2.h"
#endif
#ifndef ROOT_TNtuple
#include "TNtuple.h"
#endif
#ifndef ROOT_TString
#include "TString.h"
#endif
#ifndef BRAT_BrDataTable
#include "BrDataTable.h"
#endif
#ifndef BRAT_BrTofDig
#include "BrTofDig.h"
#endif
#ifndef BRAT_BrRunInfoManager
#include "BrRunInfoManager.h"
#endif
#ifndef BRAT_BrTofTrackMatch
#include "BrTofTrackMatch.h"
#endif
#ifndef BRAT_BrCalibration
#include "BrCalibration.h"
#endif

//____________________________________________________________________
ClassImp(BrTofAdcGainCalModule);

//____________________________________________________________________
 BrTofAdcGainCalModule::BrTofAdcGainCalModule()
  : BrTofCalModule()
{
  // Default constructor. DO NOT USE
  SetState(kSetup);

  SetDefaultParameters();
}

//____________________________________________________________________
 BrTofAdcGainCalModule::BrTofAdcGainCalModule(const Char_t* name, 
						 const Char_t* title)
  : BrTofCalModule(name, title)
{
  // Named Constructor
  SetState(kSetup);
      
  SetDefaultParameters();
}

//____________________________________________________________________
 void BrTofAdcGainCalModule::SetDefaultParameters()
{
  SetMinGain();  // default is 200
  SetMaxGain();  // default is 1800
  SetFitStart(); // default is 300 chan before histo_max
  SetFitStop();  // default is 250 chan after histo_max
}

//____________________________________________________________________
 void BrTofAdcGainCalModule::DefineHistograms()
{
  // Define histograms. They are:
  // <fill in here>

  if (GetState() != kInit) {
    Stop("DefineHistograms", "Must be called after Init"); 
    return;  
  }

  if (fLoadAscii || fCommitAscii) {
    if (Verbose() > 5)
      Warning("DefineHistograms", "No need to declare histos "
              "since we only load calibration from ascii file"); 
    return;  
  }

  const Int_t nSlats = fParamsTof->GetNoSlats();

  TDirectory* saveDir = gDirectory; 

  fHistDir = gDirectory->mkdir(Form("%sAdcGain",GetName())); 
  fHistDir->cd();
  
  // Make histograms here
  fHistDir->mkdir("Adc_Bot");
  fHistDir->mkdir("Adc_Top");
  fHistDir->mkdir("Cal_Bot");
  fHistDir->mkdir("Cal_Top");
  
  // histos for calib:
  fHistDir->cd("Adc_Top");
  
  fTopAdc = new TH1F* [nSlats];
  for (Int_t i = 1; i <= nSlats; i++) 
    fTopAdc[i-1] = new TH1F(Form("top%03d", i), "", 150, 0.5, 4000.5);
  
  fHistDir->cd("Adc_Bot");
  
  fBotAdc = new TH1F* [nSlats];
  for (Int_t i = 1; i <= nSlats; i++) 
    fBotAdc[i-1] = new TH1F(Form("bot%03d", i), "", 150, 0.5, 4000.5);
  
  
  fHistDir->cd();
  fTSum = new TH1F("tgains", 
		   Form("%s top tubes: ADC Gains", GetName()),
		   nSlats, 0.5, nSlats + 0.5); 
  fBSum = new TH1F("bgains", 
		   Form("%s bottom tubes: ADC Gains", GetName()),
		   nSlats, 0.5, nSlats + 0.5); 

  gDirectory = saveDir;
}

//____________________________________________________________________
 void BrTofAdcGainCalModule::Init()
{
  // Job-level initialisation
  SetState(kInit);
  
  //---------------
  // base class initialization (register calibration parameters)
  BrTofCalModule::Init();
  
  Int_t nSlats = fParamsTof->GetNoSlats();
  BrCalibration::EAccessMode mode = BrCalibration::kRead;
  
  // check if we want to load adc gain cal from ascii file
  if (fLoadAscii) {
    mode = BrCalibration::kTransitional;
    fCalibration->Use("topAdcGain", mode, nSlats);
    fCalibration->Use("botAdcGain", mode, nSlats);
  }
  
  else {
    // check needed calibration, if already loaded 
    // if not, try to get them from DB
    if (!fCalibration->GetAccessMode("topPedestal")  ||
	!fCalibration->GetAccessMode("botPedestal")) {
      mode = BrCalibration::kRead;
      if (DebugLevel() > 3)
	Warning("Init", "Pedestals will be read from the SQL database");
      fCalibration->Use("topPedestal", mode, nSlats);
      fCalibration->Use("botPedestal", mode, nSlats);
      fCalibration->Use("topPedestalWidth", mode, nSlats);
      fCalibration->Use("botPedestalWidth", mode, nSlats);
    }
    
    // if we want to save them (ascii or DB)
    if (fSaveAscii || fCommitAscii) {
      mode = BrCalibration::kWrite;
      fCalibration->Use("topAdcGain", mode, nSlats);
      fCalibration->Use("botAdcGain", mode, nSlats);
    }
  }
}


//____________________________________________________________________
 void BrTofAdcGainCalModule::Begin()
{
  // Run-level initialisation
  SetState(kBegin);
  
  if (fLoadAscii) {
    ReadAscii();
    return;
  }
  
  if (!HistBooked())
    if (!fCommitAscii) {
      Abort("Begin", "Need histos for calibration");
      return;
    }
  
  // check if we got parameter revisions:
  if (!fCalibration->RevisionExists("*")) {
    Abort("Begin", "Some calibration revision are missing!");
    return;
  }
  
  for (Int_t s = 1; s <= fParamsTof->GetNoSlats(); s++) {
    fCalibration->SetTopAdcGain(s, BrTofCalibration::kCalException);
    fCalibration->SetBotAdcGain(s, BrTofCalibration::kCalException);
  }
}

//____________________________________________________________________
 void BrTofAdcGainCalModule::Event(BrEventNode* inNode, BrEventNode* outNode)
{
  // Per event method
  SetState(kEvent);

  if (fCommitAscii || fLoadAscii)
    return;
  
  // get the the matched tracks and hits table
  BrDataTable* match = 
    inNode->GetDataTable(fMatchName.Data());
  if (!match)
    match = outNode->GetDataTable(fMatchName.Data());
  
  if (!match) {
    if (Verbose() > 10) Warning("Event", "No tracks matching hits");
    return;
  }
  
  // get hit table
  BrDataTable* hits = inNode->GetDataTable(Form("DigTof %s", GetName()));
  if (!hits)
    hits = outNode->GetDataTable(Form("DigTof %s", GetName()));
  
  if (!hits) {
    if (Verbose() > 3) Warning("Event", "No tof hits at all");
    return;
  }

  // --------------------------------------------------------

  // now select right hit thanks to the match object
  // we are now sure that the hit will correspond to a valid track  
  
  Int_t nhit = hits->GetEntries(); // number of tof hits in event
  Int_t nmch = match->GetEntries();
  BrTofDig* tofHit[nmch]; // array of good hits
  
  for (Int_t h = 0; h < nhit; h++) {
    BrTofDig* hit = (BrTofDig*)hits->At(h);
    
    for (Int_t m = 0; m < nmch; m++) {
      BrTofTrackMatch* pair = (BrTofTrackMatch*)match->At(m);
      if (hit->GetSlatno() != pair->GetHitId())
	continue;
      tofHit[m] = hit;
    }
  }
  
  //------------------------------------------   
  //    can now proceed with calibration
  //------------------------------------------
  
  for(Int_t h = 0; h < nmch; h++) {
    // check calibrated slats
    Int_t   slat  = tofHit[h]->GetSlatno();
    Float_t tped  = fCalibration->GetTopPedestal(slat);
    Float_t bped  = fCalibration->GetBotPedestal(slat);

    if (!IsValid(tped) || !IsValid(bped))
      fValidSlat[slat - 1] = kFALSE;
    else
      fValidSlat[slat - 1] = kTRUE;
    
    if (!fValidSlat[slat - 1])
      continue;
    
    // ------- fill calibration histograms here:
    fTopAdc[slat-1]->Fill(tofHit[h]->GetAdcUp()   - tped);
    fBotAdc[slat-1]->Fill(tofHit[h]->GetAdcDown() - bped);
  }
}

//____________________________________________________________________
 void BrTofAdcGainCalModule::Finish()
{
  // Job-level finalisation
  SetState(kFinish);

  // if load ascii mode
  if (fLoadAscii) 
    return;
  
  // if commit mode
  if (fCommitAscii) {
    ReadAscii();
    return;
  }
  
  // --------------------------------------------------------------------
  // fit MIP peak with a gaus from peak - fFitStart to peak + fFitStop
  // --------------------------------------------------------------------

  Double_t topGain[fParamsTof->GetNoSlats()];
  Double_t botGain[fParamsTof->GetNoSlats()];
  
  Stat_t max  = 0;
  Axis_t peak = 0;
  
  for (Int_t s = 1; s <= fParamsTof->GetNoSlats(); s++) {
    if (!fValidSlat[s - 1]) {
      topGain[s - 1] = BrTofCalibration::kCalException;
      continue;
    }
    
    max = 0;
    
    // get histos
    for (Int_t i = 1; i <= fTopAdc[s-1]->GetNbinsX(); i++)
      max = TMath::Max(max, fTopAdc[s-1]->GetBinContent(i));
    
    for (Int_t i = 1; i <= fTopAdc[s-1]->GetNbinsX(); i++)
      if (fTopAdc[s-1]->GetBinContent(i) == max) {
	peak = fTopAdc[s-1]->GetBinCenter(i);
	break;
      }
    
    fTopAdc[s-1]->Fit("gaus", "Q0", "", peak - fFitStart, peak + fFitStop);
    topGain[s - 1] = ((TF1*)fTopAdc[s-1]->GetListOfFunctions()->At(0))->
      GetParameter(1);
    
    // check width and remove bad fits
    if (topGain[s - 1] < fMinGain || topGain[s - 1] > fMaxGain)
      topGain[s - 1] = BrTofCalibration::kCalException;
  }
  

  // -------------- bottom tubes

  for (Int_t s = 1; s <= fParamsTof->GetNoSlats(); s++) {
    if (!fValidSlat[s - 1]) {
      botGain[s - 1] = BrTofCalibration::kCalException;
      continue;
    }
    
    max = 0;
    
    // get histos
    for (Int_t i = 1; i <= fBotAdc[s-1]->GetNbinsX(); i++)
      max = TMath::Max(max, fBotAdc[s-1]->GetBinContent(i));
    
    for (Int_t i = 1; i <= fBotAdc[s-1]->GetNbinsX(); i++)
      if (fBotAdc[s-1]->GetBinContent(i) == max) {
	peak = fBotAdc[s-1]->GetBinCenter(i);
	break;
      }
    
    fBotAdc[s-1]->Fit("gaus", "Q0", "", peak - fFitStart, peak + fFitStop);
    
    botGain[s - 1] = ((TF1*)fBotAdc[s-1]->GetListOfFunctions()->At(0))->
      GetParameter(1);
    
    // check width and remove bad fits
    if (botGain[s - 1] < fMinGain || botGain[s - 1] > fMaxGain)
      botGain[s - 1] = BrTofCalibration::kCalException;
    
  }
  

  // fill histo and calibration parameters
  for (Int_t s = 1; s <= fParamsTof->GetNoSlats(); s++) {
    // set offset in calibration parameter element
    fCalibration->SetTopAdcGain(s, topGain[s - 1]);
    fCalibration->SetBotAdcGain(s, botGain[s - 1]);
    
    // fill summaries
    if(TMath::IsNaN(topGain[s - 1]))
      fTSum ->Fill(s, -111);
    else
      fTSum ->Fill(s, topGain[s - 1]);

    if(TMath::IsNaN(botGain[s - 1]))
      fBSum ->Fill(s, -111);
    else
      fBSum ->Fill(s, botGain[s - 1]);
    
    if (Verbose() > 2) 
      cout << " Slat " << s 
	   << " : top -> " << topGain[s - 1] 
	   << "   bot -> " << botGain[s - 1] << endl;
  }
  
  // check results by normalizing the histos by the gain 
  NormalizeHistos();

  // save calibration to ascii file
  if (fSaveAscii)
    SaveAscii();
}

//____________________________________________________________________
 void BrTofAdcGainCalModule::NormalizeHistos()
{
  // private

  // --------------
  // check results
  // --------------

  TDirectory* saveDir = gDirectory;
  
  // bot tubes
  fHistDir->cd("Cal_Bot");

  for (Int_t s = 1; s <= fParamsTof->GetNoSlats(); s++) {
    if (!IsValid(fCalibration->GetBotAdcGain(s)))
      continue;
    
    Int_t minBin = fBotAdc[s-1]->GetXaxis()->GetFirst();
    Int_t maxBin = fBotAdc[s-1]->GetXaxis()->GetLast();
    Axis_t xmin  = fBotAdc[s-1]->GetXaxis()->GetXmin()/fCalibration->GetBotAdcGain(s);  
    Axis_t xmax  = fBotAdc[s-1]->GetXaxis()->GetXmax()/fCalibration->GetBotAdcGain(s);  
    Float_t factor = (xmax - xmin)/Float_t(fBotAdc[s-1]->GetNbinsX());
    
    TH1F* h = new TH1F(Form("calBot%03d", s),
		       Form("Calibrated Adc, bot slat %d", s),
		       maxBin - minBin, xmin, xmax);

    for (Int_t i = minBin; i <= maxBin; i++) {
      Axis_t x = i*factor;
      h->Fill(x, fBotAdc[s-1]->GetBinContent(i));
    }
  }
  // top tubes
  fHistDir->cd("Cal_Top");
  
  for (Int_t s = 1; s <= fParamsTof->GetNoSlats(); s++) {
    if (!IsValid(fCalibration->GetTopAdcGain(s)))
      continue;
    
    Int_t minBin = fTopAdc[s-1]->GetXaxis()->GetFirst();
    Int_t maxBin = fTopAdc[s-1]->GetXaxis()->GetLast();
    Axis_t xmin  = fTopAdc[s-1]->GetXaxis()->GetXmin()/fCalibration->GetTopAdcGain(s);  
    Axis_t xmax  = fTopAdc[s-1]->GetXaxis()->GetXmax()/fCalibration->GetTopAdcGain(s);
    Float_t factor = (xmax - xmin)/Float_t(fTopAdc[s-1]->GetNbinsX());  
    
    TH1F* h = new TH1F(Form("calTop%03d", s),
		       Form("Calibrated Adc, top slat %d", s),
		       maxBin - minBin, xmin, xmax);
    
    for (Int_t i = minBin; i <= maxBin; i++) {
      Axis_t x = i*factor;
      h->Fill(x, fTopAdc[s-1]->GetBinContent(i));
    }
  }
}

//____________________________________________________________________
 void BrTofAdcGainCalModule::SaveAscii() 
{

  // save pedestal to ascii file
  
  BrRunInfoManager* runMan = BrRunInfoManager::Instance();
  Int_t* runs = runMan->GetRunNumbers();
  Int_t  nrun = runMan->GetNumberOfRuns();
  
  BrTofCalModule::SaveAscii();
  
  ofstream file(fCalibFile.Data(), ios::out);
  
  file << "****************************************** " << endl;
  file << "*  Calibration for Tof detector " << GetName() << endl;
  file << "*  Adc Gain calibration         " << endl;
  file << "*     Used events from run(s) ";
  for (Int_t r = 0; r < nrun; r++) file << runs[r] << " ";
  file << endl;
  file << "****************************************** " <<endl;
  file << "*" << endl;    
  file << "* slat  |     top    |    bot     " << endl;
  file << "* --------------------------------" << endl << endl;
  
  for (Int_t i = 0; i < fParamsTof->GetNoSlats(); i++) {
    Int_t slat = i + 1;
    file << setw(4) << slat 
	 << setw(15) << fCalibration->GetTopAdcGain(slat) 
	 << setw(15) << fCalibration->GetBotAdcGain(slat) << endl;
  }
  
  file << "* ------------------------------------" << endl << endl;
  
}


//____________________________________________________________________
 void BrTofAdcGainCalModule::ReadAscii() 
{
  // read calibration from file created by this module

  BrTofCalModule::ReadAscii();
  
  ifstream file(fCalibFile.Data(), ios::in);
  
  if (!file) {
    Stop("ReadAscii", "File %s was not found", fCalibFile.Data());
    return;
  }
  
  Float_t top, bot;
  Int_t slat;
  Char_t comment[256];
  
  file.getline(comment, 256);
  while(comment[0] == '*') {
    file.getline(comment, 256);
    if (DebugLevel() > 5)
      cout << comment << endl;
  } 

  for (Int_t i = 1; i <= fParamsTof->GetNoSlats(); i++) {
    file >> slat >> top >> bot;
    if (DebugLevel() > 5) 
      cout << setw(4)  << slat 
	   << setw(12) << top 
	   << setw(12) << bot << endl;
        
    fCalibration->SetTopAdcGain(slat, top);
    fCalibration->SetBotAdcGain(slat, bot);
  }

  fCalibration->SetComment("topAdcGain", fComment.Data());
  fCalibration->SetComment("botAdcGain", fComment.Data());

}

//____________________________________________________________________
 void BrTofAdcGainCalModule::Print(Option_t* option) const
{
  // Print module information
  // See BrModule::Print for options.
  // In addition this module defines the Option:
  // <fill in here>

  TString opt(option);
  opt.ToLower(); 
  
  BrModule::Print(option); 
  if (opt.Contains("d")) 
   cout << endl 
         << "  Original author: Djamel Ouerdane" << endl
         << "  Last Modifications: " << endl 
         << "    $Author: videbaek $" << endl  
         << "    $Date: 2002/09/10 16:19:37 $"   << endl 
         << "    $Revision: 1.11 $ " << endl  
         << endl 
         << "-------------------------------------------------" << endl;
}

//____________________________________________________________________
//
// $Log: BrTofAdcGainCalModule.cxx,v $
// Revision 1.11  2002/09/10 16:19:37  videbaek
// Fix problem with NaN being stored in summary histograms.
//
// Revision 1.10  2002/03/21 15:04:25  ouerdane
// added fComment member to base class module and method SetComment so that the user can set comments about the calibration at commit time. Removed mean momentum stuff in slewing cal module
//
// Revision 1.9  2001/12/20 15:28:14  ouerdane
// minor changes
//
// Revision 1.8  2001/11/05 06:59:43  ouerdane
// changed to deal with new track classes and fixed some bugs
//
// Revision 1.7  2001/10/19 15:27:55  ouerdane
// Changed some verbosity level and fixed small things in Init
//
// Revision 1.6  2001/10/08 11:27:46  cholm
// Changed to use new DB access classes
//
// Revision 1.5  2001/10/07 13:19:38  ouerdane
// Changed fit parameters and function to gaus
//
// Revision 1.4  2001/10/02 01:53:28  ouerdane
// Added SetSaveAscii, SetLoadAscii, SetCommitAscii and updated the way parameters are tagged for Use
//
// Revision 1.3  2001/08/04 16:04:16  ouerdane
// minor update in the Begin method
//
// Revision 1.2  2001/07/31 09:22:20  ouerdane
// Removed all references to a TList member, replaced
// by histogram members, declare ntuple if fNtuple is true (set
// via SetNtuple)
//
// Revision 1.1  2001/06/29 13:52:36  cholm
// Imported TOF classes from Djamel
//
//

This page automatically generated by script docBrat by Christian Holm

Copyright ; 2002 BRAHMS Collaboration <brahmlib@rcf.rhic.bnl.gov>
Last Update on by

Validate HTML
Validate CSS