|
//____________________________________________________________________ // // TDC gain calibration module // Use events from TDC calibration run (e.g 2544) // Fill histograms and profile, check linearity of the TDCs //____________________________________________________________________ // // $Id: BrBbTdcGainCalModule.cxx,v 1.6 2002/03/21 15:02:17 ouerdane Exp $ // $Author: ouerdane $ // $Date: 2002/03/21 15:02:17 $ // $Copyright: (C) 2001 BRAHMS Collaboration <brahmlib@rhic.bnl.gov> // #ifndef BRAT_BrBbTdcGainCalModule #include "BrBbTdcGainCalModule.h" #endif #ifndef ROOT_TDirectory #include "TDirectory.h" #endif #ifndef WIN32 #include <iostream> #include <iomanip> #include <fstream> #else #include <iostream.h> #include <iomanip.h> #include <fstream.h> #endif #ifndef BRAT_BrDataTable #include "BrDataTable.h" #endif #ifndef ROOT_TString #include "TString.h" #endif #ifndef BRAT_BrBbDig #include "BrBbDig.h" #endif #ifndef BRAT_BrRunInfoManager #include "BrRunInfoManager.h" #endif #ifndef ROOT_TH2 #include "TH2.h" #endif #ifndef ROOT_TFile #include "TFile.h" #endif //____________________________________________________________________ ClassImp(BrBbTdcGainCalModule); //____________________________________________________________________ BrBbTdcGainCalModule::BrBbTdcGainCalModule() { // Default constructor. DO NOT USE SetState(kSetup); SetTimeStep(); // default is 10 ns SetMinGain(); // default is 0.023 ns / chan SetMaxGain(); // default is 0.027 ns / chan } //____________________________________________________________________ BrBbTdcGainCalModule::BrBbTdcGainCalModule(const Char_t* name, const Char_t* title) : BrBbCalModule(name, title) { // Named Constructor SetState(kSetup); SetTimeStep(); // default is 10 ns SetMinGain(); // default is 0.023 ns / chan SetMaxGain(); // default is 0.027 ns / chan } //____________________________________________________________________ void BrBbTdcGainCalModule::DefineHistograms() { // Define histograms. They are: // <fill in here> if (GetState() != kInit) { Stop("DefineHistograms", "Must be called after Init"); return; } TDirectory* saveDir = gDirectory; // list of histograms TDirectory* histDir = gDirectory->mkdir(Form("%sTdcGains",GetName())); histDir->cd(); // Make histograms here histDir->mkdir("Tubes"); // top tube histos histDir->cd("Tubes"); fTdc = new TH1F* [fNoTubes]; for (Int_t i = 1; i <= fNoTubes; i++) fTdc[i-1] = new TH1F(Form("%sTdc%d", GetName(), i), Form("%s Tdc, tube %d", GetName(), i), 4000, 0.5, 4000.5); fTdcGain = new TH1F* [fNoTubes]; for (Int_t i = 1; i <= fNoTubes; i++) { fTdcGain[i-1] = new TH1F(Form("%s%d", GetName(), i), Form("%s tube %d", GetName(), i), 1200, 0.5, 120.5); fTdcGain[i-1]->SetMarkerStyle(22); fTdcGain[i-1]->SetMarkerSize(1.2); fTdcGain[i-1]->SetMarkerColor(2); } // summaries histDir->cd(); // top tubes fSum = new TH1F("Summary", Form("%s TDC Gains tubes", GetName()), fNoTubes, 0.5, fNoTubes + 0.5); fSum->SetMarkerStyle(22); fSum->SetMarkerSize(1.1); fSum->SetMarkerColor(2); gDirectory = saveDir; } //____________________________________________________________________ void BrBbTdcGainCalModule::Init() { // Job-level initialisation SetState(kInit); // init base class (BrBbCalibration and BrDetectorParamsBB) BrBbCalModule::Init(); fCalibration->Use("tdcGain", BrCalibration::kWrite, fNoTubes); } //____________________________________________________________________ void BrBbTdcGainCalModule::Begin() { // Run-level initialisation SetState(kBegin); // if reading calibration parameters from ascii file to commit them, // don't need to check if histos are booked if (fCommitAscii) return; if (!HistBooked()) { Fatal("Begin", "Histograms must be booked for the" " calibration to work"); return; } } //____________________________________________________________________ void BrBbTdcGainCalModule::Event(BrEventNode* inNode, BrEventNode* outNode) { // Fill histograms with TDC values for top and bot tubes // note that you should use only a calibration run (e.g. 2544) // Per event method SetState(kEvent); BrDataTable* hits = (BrDataTable*)inNode->GetDataTable(fDigName); if (!hits) { if (DebugLevel() > 5) Warning("Event", "No data table for %s", GetName()); return; } if (!hits->GetEntries()) { if (DebugLevel() > 5) Warning("Event", "No hits in datatable for %s", GetName()); return; } for(Int_t i = 0; i < fNoTubes; i++) { BrBbDig* hit = (BrBbDig*)hits->At(i); Int_t tube = hit->GetTubeNo(); if (hit->GetTdc() < 2 && hit->GetTdc() >= 4080) continue; fTdc[tube-1]->Fill(hit->GetTdc()); } } //____________________________________________________________________ void BrBbTdcGainCalModule::End() { // Run-level finalisation SetState(kEnd); } //____________________________________________________________________ void BrBbTdcGainCalModule::Finish() { // Job-level finalisation SetState(kFinish); if (fCommitAscii) { ReadAscii(); return; } // determine maxima (sharp peaks) // between channel 0 and 4000, there are 10 peaks, the time step in // the calibrations are assumed to be 10 nsec for (Int_t i = 1; i <= fNoTubes; i++) { // scan histo, find each peak // result stored in profile for final fit // FIXE ME: have to set the error on tdc channel otherwise the chi2 // is very low... FindMax(fTdc[i-1], fTdcGain[i-1]); fTdcGain[i-1]->Fit("pol1", "Q0"); TF1* pol = (TF1*)fTdcGain[i-1]->GetListOfFunctions()->At(0); if (Verbose() > 3) cout << fTdcGain[i-1]->GetName() << " : " << pol->GetParameter(1) << " +/- " << pol->GetParError(1) << endl; if (1/pol->GetParameter(1) < fMinGain || 1/pol->GetParameter(1) > fMaxGain) fCalibration->SetTdcGain(i, BrBbCalibration::kCalException); else fCalibration->SetTdcGain(i, 1/pol->GetParameter(1)); } // time to summarize: for (Int_t i = 1; i <= fNoTubes; i++) fSum->Fill(i, (Stat_t)fCalibration->GetTdcGain(i)*1000); // save stuff to ascii file if required if (fSaveAscii) SaveAscii(); } //____________________________________________________________________ void BrBbTdcGainCalModule::FindMax(TH1F* h, TH1F* p) { // ---------------- // scan histo h, find each peak // result stored in profile for final fit // ---------------- Int_t bin = 0; // running bin Double_t tdcAtMax = 0; // tdc channel at max Int_t noPeak = 0; while (1) { // check if we are out of bound if (bin >= 4000) break; // increment bin to peek the start of the peak while (h->GetBinContent(bin) == 0) { bin++; if (bin >= 4000) break; } // recheck if we are out of bound if (bin >= 4000) break; // now we are at start of peak // the width is less than 50 channels // ---- get channel with max stat Int_t bRange = 50; Stat_t max = 0; Int_t binMax = 0; for (Int_t j = bin; j < bin + bRange; j++) { if( h->GetBinContent(j)>max) { max = h->GetBinContent(j); binMax = j; } } tdcAtMax = h->GetBinCenter(binMax); if (DebugLevel() > 10) cout << GetName() << " TdcGain cal: bin is: " << binMax << " with content: " << h->GetBinContent(binMax) << endl; // // ---- Find precise mean and RMS // Stat_t stat[4]; TH1F hh(*h); hh.GetXaxis()->SetRange(tdcAtMax-10,tdcAtMax+10); hh.GetStats(stat); Double_t mean = stat[2]/stat[0]; Double_t rms = TMath::Sqrt(stat[3]/stat[0]-mean*mean); Double_t err = rms/TMath::Sqrt(stat[0]); noPeak++; if (Verbose() > 5) cout << " Peak " << setw(3) << noPeak << " (" << setw(5) << stat[0] << " ) " << " at mean chan: " << setw(7) << mean << " +- " << setw(4) << err << endl; // set mean and error Float_t time = noPeak*fTimeStep; // time step in ns Double_t range = p->GetXaxis()->GetXmax() - p->GetXaxis()->GetXmin(); Double_t binX = time*p->GetNbinsX()/range; p->SetBinContent((Int_t)binX, mean); p->SetBinError((Int_t)binX, err); bin += 380; // go to next peak (380 is hardcoded...) } return; } //____________________________________________________________________ void BrBbTdcGainCalModule::SaveAscii() { // save pedestal to ascii file BrRunInfoManager* runMan = BrRunInfoManager::Instance(); const BrRunInfo* run = runMan->GetCurrentRun(); if (run->GetRunNo() == -1) { Failure("SaveAscii", "RunInfo has run number = -1"); return; } BrBbCalModule::SaveAscii(); ofstream file(fCalibFile.Data(), ios::out); file << "****************************************** " << endl; file << "* Calibration for BB counter " << GetName() << endl; file << "* Used events from run " << run->GetRunNo() << endl; file << "****************************************** " <<endl; file << "*" << endl; file << "* tube | TDC gain " << endl; file << "* --------------------" << endl << endl; for (Int_t i = 0; i < fNoTubes; i++) { Int_t tube = i + 1; file << setw(5) << tube << setw(13) << fCalibration->GetTdcGain(tube) << endl; } file << "* ------------------------------------" << endl << endl; cout << GetName() << " TDC gains saved to file " << fCalibFile.Data() << endl << endl; } //____________________________________________________________________ void BrBbTdcGainCalModule::ReadAscii() { BrBbCalModule::ReadAscii(); ifstream file(fCalibFile.Data(), ios::in); if (!file) { Failure("ReadAscii", "File %s was not found", fCalibFile.Data()); return; } Float_t gain; Int_t tube; 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 <= fNoTubes; i++) { file >> tube >> gain; if (DebugLevel() > 5) cout << setw(4) << tube << setw(12) << gain << endl; fCalibration->SetTdcGain(tube, gain); } fCalibration->SetComment("tdcGain", fComment.Data()); } //____________________________________________________________________ void BrBbTdcGainCalModule::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: ouerdane $" << endl << " $Date: 2002/03/21 15:02:17 $" << endl << " $Revision: 1.6 $ " << endl << endl << "-------------------------------------------------" << endl; } //____________________________________________________________________ // // $Log: BrBbTdcGainCalModule.cxx,v $ // Revision 1.6 2002/03/21 15:02:17 ouerdane // added fComment member to base class module and method SetComment so that the user can set comments about the calibration at commit time // // Revision 1.5 2001/12/07 21:21:59 videbaek // Change UserForRead to Use( ,mode,..) // // Revision 1.4 2001/10/15 00:35:08 ouerdane // Updated all Beam-Beam counter calibration modules (like I did for the TOF // some weeks ago): // // Common changes: // added methods Set[Save,Commit,Load]Ascii(Bool_t) // removed methods Set[ReadFrom,LoadFrom,SaveTo]File // // BrBbGainCalModule: // cleaned up a lot the code, added diagnostic histograms showing calibrated // ADC after the 1st MIP peak was found. // Still at a stage where the 1st MIP peak is the only one to be checked. // Later on, will add algorithm to find other peaks. // // Added BrBbSlewingCalModule in Makefile.am, etc // The fit function introduced is dT = a + b/sqrt(dE) + c/dE (where dE = ADC/Gain0) // // Revision 1.3 2001/09/10 18:21:23 videbaek // Added a few more debug statements. // More GetBincenter out of loop. // // Revision 1.2 2001/07/31 09:19:02 ouerdane // Forgot to insert a break statement in a while loop // // Revision 1.1 2001/07/23 11:40:13 ouerdane // Added calibration modules for Beam-Beam counters // BrBbCalModule: base class // BrBbPedCalModule: pedestals // BrBbTdcGaincalModule tdc gains // // Updated Makefile.am, LinkDef.h, Include.h // |
||||||
|
This page automatically generated by script docBrat by Christian Holm |
Copyright ; 2002 BRAHMS Collaboration
<brahmlib@rcf.rhic.bnl.gov>
|