#include <TAxis.h>
#include <TLinearFitter.h>
#include <TMath.h>
+#include <TStyle.h>
+#include <TCanvas.h>
#include <TTreeStream.h>
-
+#include <TGraphErrors.h>
+#include <TDirectory.h>
+#include <TTreeStream.h>
+#include <TF1.h>
//header file
#include "AliTRDCalibraVdriftLinearFit.h"
AliTRDCalibraVdriftLinearFit::AliTRDCalibraVdriftLinearFit() : /*FOLD00*/
TObject(),
fVersion(0),
+ fNameCalibUsed(""),
fLinearFitterHistoArray(540),
fLinearFitterPArray(540),
- fLinearFitterEArray(540)
+ fLinearFitterEArray(540),
+ fRobustFit(kTRUE),
+ fDebugStreamer(0x0),
+ fDebugLevel(0),
+ fSeeDetector(0)
{
//
// default constructor
AliTRDCalibraVdriftLinearFit::AliTRDCalibraVdriftLinearFit(const AliTRDCalibraVdriftLinearFit &ped) : /*FOLD00*/
TObject(ped),
fVersion(ped.fVersion),
+ fNameCalibUsed(ped.fNameCalibUsed),
fLinearFitterHistoArray(540),
fLinearFitterPArray(540),
- fLinearFitterEArray(540)
+ fLinearFitterEArray(540),
+ fRobustFit(kTRUE),
+ fDebugStreamer(0x0),
+ fDebugLevel(0),
+ fSeeDetector(0)
{
//
// copy constructor
AliTRDCalibraVdriftLinearFit::AliTRDCalibraVdriftLinearFit(const TObjArray &obja) : /*FOLD00*/
TObject(),
fVersion(0),
+ fNameCalibUsed(""),
fLinearFitterHistoArray(540),
fLinearFitterPArray(540),
- fLinearFitterEArray(540)
+ fLinearFitterEArray(540),
+ fRobustFit(kTRUE),
+ fDebugStreamer(0x0),
+ fDebugLevel(0),
+ fSeeDetector(0)
{
//
// constructor from a TObjArray
fLinearFitterPArray.Delete();
fLinearFitterEArray.Delete();
+ if ( fDebugStreamer ) delete fDebugStreamer;
+
}
//_____________________________________________________________________________
void AliTRDCalibraVdriftLinearFit::Copy(TObject &c) const
TAxis *yaxis = linearfitterhisto->GetYaxis();
TLinearFitter linearfitter = TLinearFitter(2,"pol1");
//printf("test\n");
+ Double_t integral = linearfitterhisto->Integral();
+ //printf("Integral is %f\n",integral);
+ Bool_t securitybreaking = kFALSE;
+ if(TMath::Abs(integral-1199) < 0.00001) securitybreaking = kTRUE;
for(Int_t ibinx = 0; ibinx < linearfitterhisto->GetNbinsX(); ibinx++){
for(Int_t ibiny = 0; ibiny < linearfitterhisto->GetNbinsY(); ibiny++){
if(linearfitterhisto->GetBinContent(ibinx+1,ibiny+1)>0){
Double_t x = xaxis->GetBinCenter(ibinx+1);
Double_t y = yaxis->GetBinCenter(ibiny+1);
+
for(Int_t k = 0; k < (Int_t)linearfitterhisto->GetBinContent(ibinx+1,ibiny+1); k++){
- linearfitter.AddPoint(&x,y);
- arrayI[cb]++;
+ if(!securitybreaking){
+ linearfitter.AddPoint(&x,y);
+ arrayI[cb]++;
+ }
+ else {
+ if(arrayI[cb]< 1198){
+ linearfitter.AddPoint(&x,y);
+ arrayI[cb]++;
+ }
+ }
}
+
}
}
}
TVectorD *par = new TVectorD(2);
TVectorD pare = TVectorD(2);
TVectorD *parE = new TVectorD(3);
- if((linearfitter.EvalRobust(0.8)==0)) {
- //if((linearfitter.Eval()==0)) {
-
+ //printf("Fit\n");
+ if(((fRobustFit) && (linearfitter.EvalRobust(0.8)==0)) || ((!fRobustFit) && (linearfitter.Eval()==0))) {
+ //if((linearfitter.Eval()==0)) {
+ //printf("Take the param\n");
linearfitter.GetParameters(*par);
+ //printf("Done\n");
//linearfitter.GetErrors(pare);
//Float_t ppointError = TMath::Sqrt(TMath::Abs(linearfitter.GetChisquare())/arrayI[cb]);
//(*parE)[0] = pare[0]*ppointError;
//par->Print();
//parE->Print();
}
+ //printf("Finish\n");
}
//delete linearfitterhisto;
}
+//____________Functions fit Online CH2d________________________________________
+void AliTRDCalibraVdriftLinearFit::FillPEArray2()
+{
+ //
+ // Fill fFitterPArray and fFitterEArray from inside
+ //
+
+ // Loop over histos
+ TF1 *f1 = new TF1("f1","[0]+[1]*x",-1,1);
+
+ for(Int_t cb = 0; cb < 540; cb++){
+ //const TH2S *fitterhisto = (TH2S*)fLinearFitterHistoArray.UncheckedAt(cb);
+ TH2S *fitterhisto = (TH2S*)fLinearFitterHistoArray.UncheckedAt(cb);
+ //printf("Processing the detector cb %d we find %d\n",cb, (Bool_t) fitterhisto);
+
+ if ( fitterhisto != 0 ){
+
+ Int_t nEntries=0;
+ TGraphErrors *gg=DrawMS(fitterhisto,nEntries);
+ if (!gg) continue;
+ // Number of points of the TGraphErrors
+ if(gg->GetN() < 20) {
+ if(gg) delete gg;
+ continue;
+ }
+ //printf("det %d, number of points %d, nEntries %d\n",cb,gg->GetN(),nEntries);
+ gg->Fit(f1,"Q0");
+
+ TVectorD *par = new TVectorD(2);
+ TVectorD *parE = new TVectorD(3);
+ (*parE)[0] = 0.0;
+ (*parE)[1] = 0.0;
+ (*parE)[2] = (Double_t) nEntries;
+ (*par)[0] = f1->GetParameter(0);
+ (*par)[1] = f1->GetParameter(1);
+ fLinearFitterPArray.AddAt(par,cb);
+ fLinearFitterEArray.AddAt(parE,cb);
+
+ if(fDebugLevel==0) {
+ if(gg) delete gg;
+ }
+ else {
+ if(cb==fSeeDetector) {
+ gStyle->SetPalette(1);
+ gStyle->SetOptStat(1111);
+ gStyle->SetPadBorderMode(0);
+ gStyle->SetCanvasColor(10);
+ gStyle->SetPadLeftMargin(0.13);
+ gStyle->SetPadRightMargin(0.01);
+ TCanvas *stat = new TCanvas("stat","",50,50,600,800);
+ stat->cd(1);
+ fitterhisto->Draw("colz");
+ gg->Draw("P");
+ f1->Draw("same");
+ break;
+ }
+ }
+ }
+ }
+ if(fDebugLevel==0) delete f1;
+}
+
+//_________Helper function__________________________________________________
+TGraphErrors* AliTRDCalibraVdriftLinearFit::DrawMS(const TH2 *const h2, Int_t &nEntries)
+{
+ //
+ // Debug function
+ //
+ TF1 fg("fg", "gaus", -10., 30.);
+ TGraphErrors *gp = new TGraphErrors();
+
+ TAxis *ax(h2->GetXaxis());
+ TAxis *ay(h2->GetYaxis());
+ TH1D *h1(NULL);
+ for(Int_t ipt(0), jpt(1), ig(0); ipt<ax->GetNbins(); ipt++, jpt++){
+ h1 = h2->ProjectionY("py", jpt, jpt);
+ fg.SetParameter(1, h1->GetMean());
+ //Float_t x(ax->GetBinCenter(jpt));
+ Int_t n=Int_t(h1->Integral(1, h1->GetNbinsX()));
+ nEntries+=n;
+ if(n < 15){
+ //Warning("drawMS()", Form("reject x[%d]=%f on n=%d", jpt, x, n));
+ continue;
+ }
+ h1->Fit(&fg, "WWQ0");
+ if(fg.GetNDF()<2){
+ //Warning("drawMS()", Form("reject x[%d]=%f on NDF=%d", jpt, x, fg.GetNDF()));
+ continue;
+ }
+ if(((fg.GetParameter(1)+fg.GetParameter(2)/2)>ay->GetXmax()) || ((fg.GetParameter(1)-fg.GetParameter(2)/2)<ay->GetXmin()) || (TMath::Abs(fg.GetParameter(0))< 0.00001)) continue;
+ gp->SetPoint(ig, ax->GetBinCenter(jpt), fg.GetParameter(1));
+ gp->SetPointError(ig, 0, TMath::Sqrt(pow(fg.GetParError(1),2) + (1/pow(fg.GetParameter(0),2))));
+ ig++;
+ }
+ delete h1;
+ return gp;
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff