From: mivanov Date: Fri, 5 Sep 2014 09:43:14 +0000 (+0200) Subject: ATO-78 - Technical changes to compare different calibrations X-Git-Url: http://git.uio.no/git/?a=commitdiff_plain;h=0e00d63077b01a5f43bbef21e03fc3bf39b748e5;p=u%2Fmrichter%2FAliRoot.git ATO-78 - Technical changes to compare different calibrations TTree interface: 1.) AliTPCCalPad::MakePadFromTree : Import TTree::Draw("query") --> AliTPCCalPad 2.) AliTPCCalPad::AddFriend(TTree * treePad, const char *friendName, const char *fname) :Export AliTPCCalPad -> TTree Filters: 1.) Bool_t MedianFilter(Int_t deltaRow, Int_t deltaPad); - Replace the pad calibration with median in neibrhood 2.) Bool_t LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type) - Replace the pad calibration with LTM statistic in neiborhood --- diff --git a/TPC/Base/AliTPCCalPad.cxx b/TPC/Base/AliTPCCalPad.cxx index e91e4e8d36c..1b3ef412bc3 100644 --- a/TPC/Base/AliTPCCalPad.cxx +++ b/TPC/Base/AliTPCCalPad.cxx @@ -45,8 +45,8 @@ #include #include #include - - +#include "AliTPCPreprocessorOnline.h" +#include "AliTPCCalibViewer.h" ClassImp(AliTPCCalPad) //_____________________________________________________________________________ @@ -159,7 +159,31 @@ void AliTPCCalPad::SetCalROC(AliTPCCalROC* roc, Int_t sector){ fROC[sector]->SetValue(ichannel, roc->GetValue(ichannel)); } +Bool_t AliTPCCalPad::MedianFilter(Int_t deltaRow, Int_t deltaPad){ + // + // replace constent with median in the neigborhood + // + Bool_t isOK=kTRUE; + for (Int_t isec = 0; isec < kNsec; isec++) { + if (fROC[isec]){ + isOK&=fROC[isec]->MedianFilter(deltaRow,deltaPad); + } + } + return isOK; +} +Bool_t AliTPCCalPad::LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type){ + // + // replace constent with LTM statistic in neigborhood + // + Bool_t isOK=kTRUE; + for (Int_t isec = 0; isec < kNsec; isec++) { + if (fROC[isec]){ + isOK&=fROC[isec]->LTMFilter(deltaRow, deltaPad,fraction,type); + } + } + return isOK; +} //_____________________________________________________________________________ void AliTPCCalPad::Add(Float_t c1) @@ -884,3 +908,36 @@ AliTPCCalPad * AliTPCCalPad::MakeCalPadFromHistoRPHI(TH2 * hisA, TH2* hisC){ } return calPad; } + +AliTPCCalPad *AliTPCCalPad::MakePadFromTree(TTree * treePad, const char *query, const char* name){ + // + // make cal pad from the tree + // + if (treePad->GetEntries()!=kNsec) return 0; + AliTPCCalPad * calPad= new AliTPCCalPad(name,name); + if (name) calPad->SetName(name); + for (Int_t iSec=0; iSec<72; iSec++){ + AliTPCCalROC* calROC = calPad->GetCalROC(iSec); + UInt_t nchannels = (UInt_t)treePad->Draw(query,"1","goff",1,iSec); + if (nchannels!=calROC->GetNchannels()) { + ::Error("AliTPCCalPad::MakePad",TString::Format("Wrong query sector\t%d\t%d",iSec,nchannels).Data()); + break; + } + for (UInt_t index=0; indexSetValue(index,treePad->GetV1()[index]); + } + return calPad; +} + +void AliTPCCalPad::AddFriend(TTree * treePad, const char *friendName, const char *fname){ + // + // + // + TObjArray *fArray = new TObjArray(1); + fArray->AddLast(this); + this->SetName(friendName); + AliTPCCalibViewer::MakeTree(fname, fArray,0); + TFile * f = TFile::Open(fname); + TTree * tree = (TTree*)f->Get("calPads"); + treePad->AddFriend(tree,friendName); + // tree->AddFriend(TString::Format("%s = calPads",friendName).Data(),fname); +} diff --git a/TPC/Base/AliTPCCalPad.h b/TPC/Base/AliTPCCalPad.h index 3fe540763f3..814718e8d8a 100644 --- a/TPC/Base/AliTPCCalPad.h +++ b/TPC/Base/AliTPCCalPad.h @@ -43,6 +43,13 @@ class AliTPCCalPad : public TNamed { AliTPCCalROC *GetCalROC(Int_t sector) const {return fROC[sector]; }; void SetCalROC(AliTPCCalROC* roc, Int_t sector = -1); virtual void Draw(Option_t* option = ""); + // TTree functions + static AliTPCCalPad *MakePadFromTree(TTree * tree, const char *query, const char*name=0); + void AddFriend(TTree * tree, const char *friendName, const char *fname=0); + // + // convolution + Bool_t MedianFilter(Int_t deltaRow, Int_t deltaPad); + Bool_t LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type); // // algebra void Add(Float_t c1); // add constant c1 to all channels of all ROCs diff --git a/TPC/Base/AliTPCCalROC.cxx b/TPC/Base/AliTPCCalROC.cxx index 2aafaa5b255..933b4aed0b6 100644 --- a/TPC/Base/AliTPCCalROC.cxx +++ b/TPC/Base/AliTPCCalROC.cxx @@ -143,6 +143,74 @@ void AliTPCCalROC::Streamer(TBuffer &R__b) } } +Bool_t AliTPCCalROC::MedianFilter(Int_t deltaRow, Int_t deltaPad){ + // + // Modify content of the class + // write median + Float_t *newBuffer=new Float_t[fNChannels] ; + Float_t *cacheBuffer=new Float_t[fNChannels] ; + // + for (UInt_t iRow=0; iRow=jPads) continue; + cacheBuffer[counter++]=GetValue(jRow,jPad); + } + } + newBuffer[fkIndexes[iRow]+iPad] = TMath::Median(counter,cacheBuffer); + } + } + memcpy(fData, newBuffer,GetNchannels()*sizeof(Float_t)); + delete []newBuffer; + delete []cacheBuffer; + return kTRUE; +} + +Bool_t AliTPCCalROC::LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type){ + // + // + // // + // Modify content of the class + // write median + if (fraction<0 || fraction>1) return kFALSE; + Float_t *newBuffer=new Float_t[fNChannels] ; + Double_t *cacheBuffer=new Double_t[fNChannels]; + // + for (UInt_t iRow=0; iRow=jPads) continue; + cacheBuffer[counter++]=GetValue(jRow,jPad); + } + } + Double_t mean,rms; + AliMathBase::EvaluateUni(counter,cacheBuffer,mean,rms,fraction*Double_t(counter)); + newBuffer[fkIndexes[iRow]+iPad] = (type==0)? mean:rms; + } + } + memcpy(fData, newBuffer,GetNchannels()*sizeof(Float_t)); + delete []newBuffer; + delete []cacheBuffer; + return kTRUE; +} + + // algebra fuctions: diff --git a/TPC/Base/AliTPCCalROC.h b/TPC/Base/AliTPCCalROC.h index 43c59e16676..374ba8523a2 100644 --- a/TPC/Base/AliTPCCalROC.h +++ b/TPC/Base/AliTPCCalROC.h @@ -38,6 +38,9 @@ class AliTPCCalROC : public TNamed { void SetValue(UInt_t channel, Float_t vd) {fData[channel]= vd; }; virtual void Draw(Option_t* option = ""); // + Bool_t MedianFilter(Int_t deltaRow, Int_t deltaPad); + Bool_t LTMFilter(Int_t deltaRow, Int_t deltaPad, Float_t fraction, Int_t type); + // // algebra void Add(Float_t c1); // add c1 to each channel of the ROC void Multiply(Float_t c1); // multiply each channel of the ROC with c1 diff --git a/TPC/Base/test/UnitTest.C b/TPC/Base/test/UnitTest.C new file mode 100644 index 00000000000..acdeae6931d --- /dev/null +++ b/TPC/Base/test/UnitTest.C @@ -0,0 +1,62 @@ +/* + Unit test for fucntions used in the calibrations + gSystem->SetIncludePath("-I$ROOTSYS/include -I$ALICE_ROOT/ -I$ALICE_ROOT/include -I$ALICE_ROOT/STEER\ + -I$ALICE_ROOT/TPC -I$ALICE_ROOT/ITS -I$ALICE_ROOT/TRD -I$ALICE_ROOT/TOF -I$ALICE_ROOT/RAW -I$ALICE_ROOT/PWG1 -I$ALICE_ROOT/STAT -I$ALICE_ROOT/TPC/Base -I$ALICE_ROOT/TPC/Calib"); + + .L $ALICE_ROOT/TPC/Base/test/UnitTest.C+ + +*/ + +#include "TF1.h" +#include "TMath.h" +#include "TLinearFitter.h" +#include "TFile.h" +#include "AliSysInfo.h" +#include "TTree.h" +#include "AliLog.h" +#include "THn.h" +#include "TRandom.h" +#include "AliTPCCalPad.h" +#include "AliTPCCalibViewer.h" +#include "AliTPCcalibDButil.h" + +// +// PARAMETERS: +// +TString baseDir="/hera/alice/wiechula/calib/guiTrees"; +// +// + + +void UnitTestAliTPCCalPadTree(){ + // + // + // + TObjArray *fArray = new TObjArray(100); + TTree * treePad=AliTPCcalibDButil::ConnectGainTrees(baseDir); + for (Int_t i=0; i<3; i++){ + AliTPCCalPad * padLx = AliTPCCalPad::MakePadFromTree(treePad,"lx.fElements","QMax"); + AliTPCCalPad * padLy = AliTPCCalPad::MakePadFromTree(treePad,"ly.fElements","QMax"); + AliTPCCalPad * padMax = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MaxCharge.fElements","QMax"); + AliTPCCalPad * padMean = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MeanCharge.fElements","QTot"); + if (i>0) { + padLx->MedianFilter(i,2*i); + padLy->MedianFilter(i,2*i); + padMax->MedianFilter(i,2*i); + padMean->MedianFilter(i,2*i); + } + padLx->SetName(TString::Format("QLx%d",i).Data()); + padLy->SetName(TString::Format("QLy%d",i).Data()); + padMax->SetName(TString::Format("QMax%d",i).Data()); + padMean->SetName(TString::Format("QTot%d",i).Data()); + fArray->AddLast(padLx); + fArray->AddLast(padLy); + fArray->AddLast(padMax); + fArray->AddLast(padMean); + } + AliTPCCalibViewer::MakeTree("QAtest.root", fArray,0); + // + TFile*fout= TFile::Open("QAtest.root"); + TTree * tree = (TTree*)fout->Get("calPads"); + +}