+const TF1* GetPeriodRPFit(const char* histName)
+{
+ TString name = Form("f1RPFit%s", histName);
+
+ // if created earlier, find in list
+ static THashList* list = new THashList;
+ TF1* func = (TF1*) list->FindObject(name.Data());
+ if( func ) return func;
+
+ TFile* file = TFile::Open(fullMergeFileName, "read");
+ TList* hlist = (TList*) file->Get("PHOSPi0Flow_kCentral/PHOSPi0Flow_kCentralCoutput1");
+ TH2F* hist2d = (TH2F*) hlist->FindObject(histName);
+ TH1D* hist = hist2d->ProjectionX();
+ int nBins = hist->GetNbinsX();
+ const Double_t* array = & hist->GetArray()[1];
+
+ //TCanvas* canv = new TCanvas;
+ func = new TF1(name.Data(), "[0]*(1. + [1]*sin([2]*x + [3]))", 0, TMath::Pi());
+ double mean = TMath::Mean(nBins, array);
+ double rms = TMath::RMS(nBins, array) / mean;
+ func->SetParameters(mean, rms, 0.6*TMath::Pi(), -0.2*TMath::Pi());
+ func->SetParLimits(0, 0, 2*mean);
+ func->SetParLimits(1, 0, 1);
+ func->SetParLimits(2, 0, 100);
+ func->SetParLimits(3, -TMath::Pi(), TMath::Pi());
+ func->SetParNames("s_{0}", "s_{1}", "#omega", "#psi");
+ hist->GetXaxis()->SetTitle("#phi");
+ int error = 0;
+ TCanvas* canv = new TCanvas(name.Data());
+ gStyle->SetOptStat(0);
+ gStyle->SetOptFit(1);
+ int error = hist->Fit(func, "Q");
+ if( ! name.Contains("flat") ) {
+ gStyle->SetOptFit(1);
+ } else {
+ gStyle->SetOptFit(1);
+ }
+ hist->Draw("E");
+ canv->SaveAs(Form("imgs/%s.pdf", name.Data()));
+ if( error )
+ Printf("GetPeriodRPFit::ERROR, could not fit %s, error:%d", histName, error);
+
+ //hist->Draw();
+ //func->Draw("same");
+ //canv->SaveAs("rpfit.pdf");
+
+ // Fit
+ list->Add(func);
+ return func;
+}
+
+
+//-----------------------------------------------------------------------------
+void AddNoise(TH1D* hist, double noise)
+{
+ int nBins = hist->GetNbinsX();
+ double mean = TMath::Mean(nBins, &((hist->GetArray())[1]) );
+ //Printf("Adding noise to: %s", hist->GetName());
+ for(int bin=1; bin<=nBins; ++bin) {
+ hist->SetBinContent( bin, hist->GetBinContent(bin) + gRandom->Gaus(0, noise) );
+ double err = hist->GetBinError(bin);
+ hist->SetBinError( bin, TMath::Sqrt( err**2 + noise**2 ) );
+ }
+}
+