1 /**************************************************************************
2 * This file is property of and copyright by the Experimental Nuclear *
3 * Physics Group, Dep. of Physics *
4 * University of Oslo, Norway, 2007 *
6 * Author: Per Thomas Hille <perthi@fys.uio.no> for the ALICE HLT Project.*
7 * Contributors are mentioned in the code where appropriate. *
8 * Please report bugs to perthi@fys.uio.no *
10 * Permission to use, copy, modify and distribute this software and its *
11 * documentation strictly for non-commercial purposes is hereby granted *
12 * without fee, provided that the above copyright notice appears in all *
13 * copies and that both the copyright notice and this permission notice *
14 * appear in the supporting documentation. The authors make no claims *
15 * about the suitability of this software for any purpose. It is *
16 * provided "as is" without express or implied warranty. *
17 **************************************************************************/
18 #include "AliCaloRawAnalyzerComparison.h"
20 #include "AliCaloRawAnalyzerCrude.h"
21 #include "AliCaloRawAnalyzerKStandard.h"
22 #include "AliCaloRawAnalyzerFastFit.h"
23 #include "AliCaloRawAnalyzerNN.h"
24 #include "AliCaloRawAnalyzerPeakFinder.h"
26 #include "AliCaloFitResults.h"
35 AliCaloRawAnalyzerComparison::AliCaloRawAnalyzerComparison() : fMod(0),
43 fReferenceAnalyzer = new AliCaloRawAnalyzerKStandard();
45 AliCaloRawAnalyzerCrude *crude = new AliCaloRawAnalyzerCrude();
46 AliCaloRawAnalyzerKStandard *lms = new AliCaloRawAnalyzerKStandard();
47 AliCaloRawAnalyzerFastFit *fastfit = new AliCaloRawAnalyzerFastFit();
48 AliCaloRawAnalyzerNN *neuralnet = new AliCaloRawAnalyzerNN();
49 AliCaloRawAnalyzerPeakFinder *peakfinder = new AliCaloRawAnalyzerPeakFinder();
51 fReferenceAnalyzer->SetIsZeroSuppressed();
52 peakfinder ->SetIsZeroSuppressed();
53 lms ->SetIsZeroSuppressed();
54 fastfit ->SetIsZeroSuppressed();
55 neuralnet ->SetIsZeroSuppressed();
56 crude ->SetIsZeroSuppressed();
58 fRawAnalyzers.push_back( (AliCaloRawAnalyzer*)crude);
59 fRawAnalyzers.push_back( (AliCaloRawAnalyzer*)lms);
60 fRawAnalyzers.push_back( (AliCaloRawAnalyzer*)fastfit);
61 fRawAnalyzers.push_back( (AliCaloRawAnalyzer*)neuralnet);
63 fRawAnalyzers.push_back( (AliCaloRawAnalyzer*)peakfinder);
64 InitHistograms( fRawAnalyzers, fReferenceAnalyzer );
68 AliCaloRawAnalyzerComparison::Evaluate( const vector<AliCaloBunchInfo> &bunchvector, const UInt_t altrocfg1,
69 const UInt_t altrocfg2, const int event, const int col, const int row )
71 // cout << __FILE__ << __LINE__ << endl;
73 AliCaloFitResults ref = fReferenceAnalyzer->Evaluate( bunchvector, altrocfg1, altrocfg2 );
74 // AliCaloFitResults ref;
76 for(int i=0; i < fRawAnalyzers.size(); i++ )
78 AliCaloFitResults an = fRawAnalyzers.at(i)->Evaluate( bunchvector, altrocfg1, altrocfg2 );
80 // cout << __FILE__ << __LINE__ << ":" << fRawAnalyzers.at(i)->GetAlgoAbbr() << " col = " << col << " row " << row << " amp= " << an.GetAmp() << endl ;
82 // fAmplitudeVsEvent[i]->Fill( event, an.GetAmp() );
83 // fTofVsEvent[i]->Fill( event, an.GetTof() );
85 fRefAmpVsAnalyzers[i]->Fill(ref.GetAmp(), an.GetAmp() );
86 fRefTofVsAnalyzers[i]->Fill(ref.GetTof(), an.GetTof() );
87 fAmpDiff[i]->Fill(ref.GetAmp() - an.GetAmp() );
88 fTofDiff[i]->Fill(ref.GetTof() - an.GetTof() );
90 if(col == fMonCol1 && fMonRow1 == row )
92 fAmplitudeVsEvent[i]->Fill( event, an.GetAmp() );
93 fTofVsEvent[i]->Fill( event, an.GetTof() );
97 if(col == fMonCol2 && row == fMonRow2 )
102 if(row >= 0 && col >= 0 )
104 fAmpHistograms[i][col][row]->Fill( an.GetAmp() );
110 AliCaloRawAnalyzerComparison::EventChanged()
112 for(int i=0; i < fRawAnalyzers.size(); i++ )
115 if( ( fMon1[i].GetAmp() > 50 && fMon2[i].GetAmp() > 50 ) && ( fMon1[i].GetAmp() < 1023 && fMon2[i].GetAmp() < 1023 ) )
117 cout << __FILE__ << __LINE__ << fRawAnalyzers.at(i)->GetAlgoAbbr() << "\tamp1=" << fMon1[i].GetAmp() << "\tamp2 = " << fMon2[i].GetAmp() << endl;
121 if( fMon1[i].GetAmp() > 50 && fMon2[i].GetAmp() > 50 )
123 // if( fMon1[i].GetTof() != fMon2[i].GetTof())
125 // fTofResDifferential[i]->Fill( ( fMon1[i].GetTof() - fMon2[i].GetTof())/TMath::Sqrt(2) );
126 fTofResDifferential[i]->Fill( ( fMon1[i].GetTof() - fMon2[i].GetTof()) );
127 fTofResAbsolute[i]->Fill( fMon1[i].GetTof() );
135 AliCaloRawAnalyzerComparison::WriteHistograms()
136 { // write histograms
137 TFile *f = new TFile("comparison2.root", "recreate");
140 for(int col=0; col < NZCOLSSMOD; col ++ )
142 for(int row=0; row < NXROWSSMOD; row ++ )
144 fAmpHistograms[col][row]->Write();
149 for(int i=0; i < fRawAnalyzers.size(); i++ )
151 for(int col=0; col < NZCOLSSMOD; col ++ )
153 for(int row=0; row < NXROWSSMOD; row ++ )
155 fAmpHistograms[i][col][row]->Write();
159 fAmplitudeVsEvent[i]->Write();
160 fTofVsEvent[i]->Write();
161 fRefAmpVsAnalyzers[i]->Write();
162 fRefTofVsAnalyzers[i]->Write();
163 fAmpDiff[i]->Write();
164 fTofDiff[i]->Write();
165 fTofResDifferential[i]->Write();
166 fTofResAbsolute[i]->Write();
175 AliCaloRawAnalyzerComparison::InitHistograms( vector <AliCaloRawAnalyzer*> analyzers, AliCaloRawAnalyzer* ref )
180 for(int col=0; col < NZCOLSSMOD; col ++ )
182 for(int row=0; row < NXROWSSMOD; row ++ )
184 sprintf(tmpname, "z(col)%d_x(row)%d_amplitude;Counts;Amplitude/ADC counts", col, row);
185 fAmpHistograms[col][row] = new TH1D(tmpname, tmpname, 1024, 0, 1023 );
190 // TH1D *fAmpHistograms[NZCOLSSMOD][NXROWSSMOD];
192 for(int i=0; i < analyzers.size(); i++ )
194 for(int col=0; col < NZCOLSSMOD; col ++ )
196 for(int row=0; row < NXROWSSMOD; row ++ )
198 sprintf(tmpname, "z(col)%d_x(row)%d_amplitude_%s;Amplitude/ADC counts;Counts", col, row, analyzers.at(i)->GetAlgoAbbr() );
199 fAmpHistograms[i][col][row] = new TH1D(tmpname, tmpname, 1024, 0, 1023 );
203 sprintf(tmpname, "%s_amplitude_vs_event_row%d_col%d;Event;Amplitude/ADC counts", analyzers.at(i)->GetAlgoAbbr(), fMonRow1, fMonCol1 );
204 fAmplitudeVsEvent[i]=new TH2D(tmpname, tmpname, 8000, 0, 7999, 1024, 0, 1023 );
206 sprintf(tmpname, "%s_tof_vs_event_row%d_col%d;Event;Amplitude/ADC counts", analyzers.at(i)->GetAlgoAbbr(), fMonRow1, fMonCol1 );
208 // fTofVsEvent[i]=new TH2D(tmpname, tmpname, 8000, 0, 7999, 2000, 3000, 5000);
209 fTofVsEvent[i]=new TH2D(tmpname, tmpname, 8000, 0, 7999, 2000, 1000, 4000);
212 sprintf(tmpname, "%s_vs_%s_ampltude; Amplitude_{%s}/ADC counts; Amplitude_{%s}/ADC counts", ref->GetAlgoAbbr(),
213 analyzers.at(i)->GetAlgoAbbr(), ref->GetAlgoAbbr(), analyzers.at(i)->GetAlgoAbbr() );
214 fRefAmpVsAnalyzers[i] = new TH2D(tmpname, tmpname, 1024, 0, 1023, 1024, 0, 1023);
215 sprintf(tmpname, "%s_vs_%s_tof; tof_{%s}/ns; tof_{%s}/ns", ref->GetAlgoAbbr(), analyzers.at(i)->GetAlgoAbbr(), ref->GetAlgoAbbr(), analyzers.at(i)->GetAlgoAbbr() );
216 fRefTofVsAnalyzers[i] = new TH2D(tmpname, tmpname, 500, 2000, 4999, 500, 2000, 4999 );
217 sprintf( tmpname, "%s-%s amplitude;counts;A_{%s} - A_{%s}", ref->GetAlgoAbbr(),
218 analyzers.at(i)->GetAlgoAbbr(), ref->GetAlgoAbbr(), analyzers.at(i)->GetAlgoAbbr());
219 fAmpDiff[i] = new TH1D(tmpname, tmpname, 100, -10, 10 );
220 sprintf( tmpname, "%s-%s tof;counts;A_{%s} - A_{%s}", ref->GetAlgoAbbr(),
221 analyzers.at(i)->GetAlgoAbbr(), ref->GetAlgoAbbr(), analyzers.at(i)->GetAlgoAbbr());
222 fTofDiff[i] = new TH1D(tmpname, tmpname, 1000, -5000, 5000 );
223 sprintf( tmpname, "%s Differential tof resolution (%d, %d) vs (%d, %d);#sigma_{tof}^{%s}/ns;Counts", analyzers.at(i)->GetAlgoAbbr(), fMonCol1, fMonRow1, fMonCol2, fMonRow2,
224 analyzers.at(i)->GetAlgoAbbr() );
225 fTofResDifferential[i] = new TH1D(tmpname, tmpname, 1000, -250, 250 );
226 sprintf( tmpname, "%s Absolute tof distribution (%d, %d);#sigma_{tof}^{%s}/ns;Counts", analyzers.at(i)->GetAlgoAbbr(), fMonCol1, fMonRow1, analyzers.at(i)->GetAlgoAbbr() );
228 fTofResAbsolute[i] = new TH1D(tmpname, tmpname, 2000 , -1000, 7000 );