diff --git a/Alignment/LaserAlignment/plugins/TkLasBeamFitter.cc b/Alignment/LaserAlignment/plugins/TkLasBeamFitter.cc
index 6d5cf707005d9..ef83ea530dcd2 100644
--- a/Alignment/LaserAlignment/plugins/TkLasBeamFitter.cc
+++ b/Alignment/LaserAlignment/plugins/TkLasBeamFitter.cc
@@ -676,19 +676,19 @@ void TkLasBeamFitter::fitter(TkFittedLasBeam &beam, AlgebraicSymMatrix &covMatri
     TF1 tecPlus("tecPlus", tecPlusFunction, zMin, zMax, nFitParams );
     tecPlus.SetParameter( 1, 0 ); // slope
     tecPlus.SetParameter( nFitParams - 1, 0 ); // BS 
-    lasData->Fit("tecPlus", "R"); // "R", "RV" or "RQ"
+    lasData->Fit(&tecPlus, "R"); // "R", "RV" or "RQ"
   }
   else if(beam.isTecInternal(-1)){
     TF1 tecMinus("tecMinus", tecMinusFunction, zMin, zMax, nFitParams );
     tecMinus.SetParameter( 1, 0 ); // slope
     tecMinus.SetParameter( nFitParams - 1, 0 ); // BS 
-    lasData->Fit("tecMinus", "R");
+    lasData->Fit(&tecMinus, "R");
   }
   else{
     TF1 at("at", atFunction, zMin, zMax, nFitParams );
     at.SetParameter( 1, 0 ); // slope
     at.SetParameter( nFitParams - 1, 0 ); // BS 
-    lasData->Fit("at","R");
+    lasData->Fit(&at,"R");
   }
   
   // get values and errors for offset and slope
diff --git a/Alignment/MuonAlignmentAlgorithms/src/MuonResidualsFitter.cc b/Alignment/MuonAlignmentAlgorithms/src/MuonResidualsFitter.cc
index a6409786319c1..6565c871fb728 100644
--- a/Alignment/MuonAlignmentAlgorithms/src/MuonResidualsFitter.cc
+++ b/Alignment/MuonAlignmentAlgorithms/src/MuonResidualsFitter.cc
@@ -503,7 +503,7 @@ void MuonResidualsFitter::histogramChi2GaussianFit(int which, double &fit_mean,
   f1->SetParameter(0, hist->GetEntries());
   f1->SetParameter(1, 0);
   f1->SetParameter(2, hist->GetRMS());
-  hist->Fit("f1","RQ");
+  hist->Fit(f1,"RQ");
   
   fit_mean  = f1->GetParameter(1);
   fit_sigma = f1->GetParameter(2);
diff --git a/CalibCalorimetry/EcalLaserAnalyzer/plugins/EcalLaserAnalyzerYousi.cc b/CalibCalorimetry/EcalLaserAnalyzer/plugins/EcalLaserAnalyzerYousi.cc
index eb3f14dcaa920..e98481b43dfa1 100644
--- a/CalibCalorimetry/EcalLaserAnalyzer/plugins/EcalLaserAnalyzerYousi.cc
+++ b/CalibCalorimetry/EcalLaserAnalyzer/plugins/EcalLaserAnalyzerYousi.cc
@@ -451,7 +451,7 @@ EcalLaserAnalyzerYousi::endJob() {
       min_r = C_APD[j]->GetBinCenter(C_APD[j]->GetMaximumBin()) - K*RMS;
       max_r = C_APD[j]->GetBinCenter(C_APD[j]->GetMaximumBin()) + K*RMS;
       gs1 = new TF1("gs1", "gaus", min_r, max_r);
-      C_APD[j]->Fit("gs1", "RQI");
+      C_APD[j]->Fit(gs1, "RQI");
       Sigma = gs1->GetParameter(2);
       K = K*1.5;
       iCount++;
@@ -471,7 +471,7 @@ EcalLaserAnalyzerYousi::endJob() {
       min_r = C_APDPN[j]->GetBinCenter(C_APDPN[j]->GetMaximumBin()) - K*RMS;
       max_r = C_APDPN[j]->GetBinCenter(C_APDPN[j]->GetMaximumBin()) + K*RMS;
       gs2 = new TF1("gs2", "gaus", min_r, max_r);
-      C_APDPN[j]->Fit("gs2", "RQI");
+      C_APDPN[j]->Fit(gs2, "RQI");
       Sigma = gs2->GetParameter(2);
       K = K*1.5;
       iCount++;
diff --git a/CalibCalorimetry/EcalPedestalOffsets/src/TPedValues.cc b/CalibCalorimetry/EcalPedestalOffsets/src/TPedValues.cc
index e5385a55ab249..c6049a65b4601 100644
--- a/CalibCalorimetry/EcalPedestalOffsets/src/TPedValues.cc
+++ b/CalibCalorimetry/EcalPedestalOffsets/src/TPedValues.cc
@@ -240,7 +240,7 @@ int TPedValues::makePlots (TFile * rootFile, const std::string & dirName,
           sprintf (name,"XTL%04d_GAIN%02d",endcapCrystalNumbers[xtl],gainHuman) ;
           graph.GetXaxis()->SetTitle("DAC value");
           graph.GetYaxis()->SetTitle("Average pedestal ADC");
-          graph.Fit("fitFunction","RWQ");
+          graph.Fit(&fitFunction,"RWQ");
           graph.Write (name);
           
           double slope1 = fitFunction.GetParameter(0);
diff --git a/CalibMuon/DTCalibration/src/DTTimeBoxFitter.cc b/CalibMuon/DTCalibration/src/DTTimeBoxFitter.cc
index 88f21ac623a93..2355ecdf5a9aa 100644
--- a/CalibMuon/DTCalibration/src/DTTimeBoxFitter.cc
+++ b/CalibMuon/DTCalibration/src/DTTimeBoxFitter.cc
@@ -84,7 +84,7 @@ pair<double, double> DTTimeBoxFitter::fitTimeBox(TH1F *hTimeBox) {
   if(theVerbosityLevel >= 2)
     option = "";
 
-  hTimeBox->Fit("IntGauss", option.c_str(), "",xFitMin, xFitMax);
+  hTimeBox->Fit(fIntGaus, option.c_str(), "",xFitMin, xFitMax);
 
   // Get fitted parameters
   double mean =  fIntGaus->GetParameter("Mean");
diff --git a/CalibTracker/SiPixelGainCalibration/plugins/SiPixelGainCalibrationAnalysis.cc b/CalibTracker/SiPixelGainCalibration/plugins/SiPixelGainCalibrationAnalysis.cc
index 2df992882bc17..c78f11d4f049c 100644
--- a/CalibTracker/SiPixelGainCalibration/plugins/SiPixelGainCalibrationAnalysis.cc
+++ b/CalibTracker/SiPixelGainCalibration/plugins/SiPixelGainCalibrationAnalysis.cc
@@ -343,7 +343,7 @@ SiPixelGainCalibrationAnalysis::doFits(uint32_t detid, std::vector<SiPixelCalibD
     graph_->SetPoint(ipointtemp,xvals[ipointtemp],yvals[ipointtemp]);
     graph_->SetPointError(ipointtemp,0,yerrvals[ipointtemp]);
   }
-  Int_t tempresult = graph_->Fit("func","FQ0N");
+  Int_t tempresult = graph_->Fit(func_,"FQ0N");
   slope=func_->GetParameter(1);
   slopeerror=func_->GetParError(1);
   intercept=func_->GetParameter(0);
@@ -359,7 +359,7 @@ SiPixelGainCalibrationAnalysis::doFits(uint32_t detid, std::vector<SiPixelCalibD
     for(int ii=0; ii<npoints; ++ii){
       edm::LogWarning("SiPixelGainCalibrationAnalysis")<< "vcal " << xvals[ii] << " response: " << yvals[ii] << "+/-" << yerrvals[ii] << std::endl; 
     } 
-    tempresult = graph_->Fit("func","FQ0NW");
+    tempresult = graph_->Fit(func_,"FQ0NW");
     slope=func_->GetParameter(1);
     slopeerror=func_->GetParError(1);
     intercept = func_->GetParameter(0);
diff --git a/CalibTracker/SiStripChannelGain/plugins/SiStripGainCosmicCalculator.cc b/CalibTracker/SiStripChannelGain/plugins/SiStripGainCosmicCalculator.cc
index f55fd370c6d63..e473dfebda3a3 100644
--- a/CalibTracker/SiStripChannelGain/plugins/SiStripGainCosmicCalculator.cc
+++ b/CalibTracker/SiStripChannelGain/plugins/SiStripGainCosmicCalculator.cc
@@ -222,8 +222,9 @@ std::pair<double,double> SiStripGainCosmicCalculator::getPeakOfLandau( TH1F * in
 //  delete landaufit;
 //
   // perform fit with standard landau
-  inputHisto->Fit("landau","0Q");
-  TF1 * fitfunction = (TF1*) inputHisto->GetListOfFunctions()->First();
+  // make our own copy to avoid problems with threads
+  std::unique_ptr<TF1> fitfunction( new TF1("landaufit","landau") );
+  inputHisto->Fit(fitfunction.get(),"0Q");
   adcs = fitfunction->GetParameter("MPV");
   error = fitfunction->GetParError(1); // MPV is parameter 1 (0=constant, 1=MPV, 2=Sigma)
   double chi2 = fitfunction->GetChisquare();
@@ -231,15 +232,14 @@ std::pair<double,double> SiStripGainCosmicCalculator::getPeakOfLandau( TH1F * in
   double chi2overndf = chi2 / ndf;
   // in case things went wrong, try to refit in smaller range
   if(adcs< 2. || (error/adcs)>1.8 ){
-     inputHisto->Fit("landau","0Q",0,0.,400.);
-     TF1 * fitfunction2 = (TF1*) inputHisto->GetListOfFunctions()->First();
+     inputHisto->Fit(fitfunction.get(),"0Q",0,0.,400.);
      std::cout<<"refitting landau for histogram "<<inputHisto->GetTitle()<<std::endl;
      std::cout<<"initial error/adcs ="<<error<<" / "<<adcs<<std::endl;
-     std::cout<<"new     error/adcs ="<<fitfunction2->GetParError(1)<<" / "<<fitfunction2->GetParameter("MPV")<<std::endl;
-     adcs = fitfunction2->GetParameter("MPV");
-     error = fitfunction2->GetParError(1); // MPV is parameter 1 (0=constant, 1=MPV, 2=Sigma)
-     chi2 = fitfunction2->GetChisquare();
-     ndf = fitfunction2->GetNDF();
+     std::cout<<"new     error/adcs ="<<fitfunction->GetParError(1)<<" / "<<fitfunction->GetParameter("MPV")<<std::endl;
+     adcs = fitfunction->GetParameter("MPV");
+     error = fitfunction->GetParError(1); // MPV is parameter 1 (0=constant, 1=MPV, 2=Sigma)
+     chi2 = fitfunction->GetChisquare();
+     ndf = fitfunction->GetNDF();
      chi2overndf = chi2 / ndf;
    }
    // if still wrong, give up
diff --git a/CalibTracker/SiStripChannelGain/plugins/SiStripGainFromCalibTree.cc b/CalibTracker/SiStripChannelGain/plugins/SiStripGainFromCalibTree.cc
index 49c8f4c5c4856..ef0b9d515bc1a 100644
--- a/CalibTracker/SiStripChannelGain/plugins/SiStripGainFromCalibTree.cc
+++ b/CalibTracker/SiStripChannelGain/plugins/SiStripGainFromCalibTree.cc
@@ -376,7 +376,7 @@ void SiStripGainFromCalibTree::getPeakOfLandau(TH1* InputHisto, double* FitResul
    // perform fit with standard landau
    TF1* MyLandau = new TF1("MyLandau","landau",LowRange, HighRange);
    MyLandau->SetParameter(1,300);
-   InputHisto->Fit("MyLandau","0QR WW");
+   InputHisto->Fit(MyLandau,"0QR WW");
 
    // MPV is parameter 1 (0=constant, 1=MPV, 2=Sigma)
    FitResults[0]         = MyLandau->GetParameter(1);  //MPV
diff --git a/CalibTracker/SiStripChannelGain/plugins/SiStripGainFromData.cc b/CalibTracker/SiStripChannelGain/plugins/SiStripGainFromData.cc
index ce78b9e14668f..04e9fbb048d4f 100644
--- a/CalibTracker/SiStripChannelGain/plugins/SiStripGainFromData.cc
+++ b/CalibTracker/SiStripChannelGain/plugins/SiStripGainFromData.cc
@@ -1255,15 +1255,14 @@ void SiStripGainFromData::getPeakOfLandau(TH1* InputHisto, double* FitResults, d
    TF1* MyLandau = new TF1("MyLandau","landau",LowRange, HighRange);
    MyLandau->SetParameter("MPV",300);
 
-   InputHisto->Fit("MyLandau","QR WW");
-   TF1 * fitfunction = (TF1*) InputHisto->GetListOfFunctions()->First();
+   InputHisto->Fit(MyLandau,"QR WW");
 
    // MPV is parameter 1 (0=constant, 1=MPV, 2=Sigma)
-    adcs        = fitfunction->GetParameter("MPV");
-    adcs_err    = fitfunction->GetParError(1);
-    width       = fitfunction->GetParameter(2);
-    width_err   = fitfunction->GetParError(2);
-    chi2overndf = fitfunction->GetChisquare() / fitfunction->GetNDF();
+    adcs        = MyLandau->GetParameter("MPV");
+    adcs_err    = MyLandau->GetParError(1);
+    width       = MyLandau->GetParameter(2);
+    width_err   = MyLandau->GetParError(2);
+    chi2overndf = MyLandau->GetChisquare() / MyLandau->GetNDF();
 
     // if still wrong, give up
     if(adcs<2. || chi2overndf>MaxChi2OverNDF){
diff --git a/CalibTracker/SiStripLorentzAngle/src/LA_Results.cc b/CalibTracker/SiStripLorentzAngle/src/LA_Results.cc
index 9b697336e6da9..1d4c8c3126f81 100644
--- a/CalibTracker/SiStripLorentzAngle/src/LA_Results.cc
+++ b/CalibTracker/SiStripLorentzAngle/src/LA_Results.cc
@@ -95,9 +95,11 @@ summarize_ensembles(Book& book) const {
       book.fill( p.measured.second, name+"_merr",         500, 0, 0.01);
       book.fill( (p.measured.first-p.reco.first)/p.measured.second, name+"_pull", 500, -10,10);
     }
-    book[name+"_measure"]->Fit("gaus","LLQ");
-    book[name+"_merr"]->Fit("gaus","LLQ");
-    book[name+"_pull"]->Fit("gaus","LLQ");
+    //Need our own copy for thread safety
+    TF1 gaus("mygaus","gaus");
+    book[name+"_measure"]->Fit(&gaus,"LLQ");
+    book[name+"_merr"]->Fit(&gaus,"LLQ");
+    book[name+"_pull"]->Fit(&gaus,"LLQ");
   }
 }
 
@@ -137,7 +139,9 @@ offset_slope(const std::vector<LA_Filler_Fitter::EnsembleSummary>& ensembles) {
       yerr.push_back(ensemble.meanMeasured.second);
     }
     TGraphErrors graph(x.size(),&(x[0]),&(y[0]),&(xerr[0]),&(yerr[0]));
-    graph.Fit("pol1");
+    //Need our own copy for thread safety
+    TF1 pol1("mypol1","pol1");
+    graph.Fit(&pol1);
     const TF1*const fit = graph.GetFunction("pol1");
     
     return std::make_pair( std::make_pair(fit->GetParameter(0), fit->GetParError(0)),
diff --git a/CalibTracker/SiStripLorentzAngle/src/SymmetryFit.cc b/CalibTracker/SiStripLorentzAngle/src/SymmetryFit.cc
index a2bc4001c95e6..bfb473415a501 100644
--- a/CalibTracker/SiStripLorentzAngle/src/SymmetryFit.cc
+++ b/CalibTracker/SiStripLorentzAngle/src/SymmetryFit.cc
@@ -1,6 +1,7 @@
 #include "CalibTracker/SiStripLorentzAngle/interface/SymmetryFit.h"
 #include <cmath>
 #include <cassert>
+#include <memory>
 #include "boost/foreach.hpp"
 
 TH1* SymmetryFit::symmetryChi2(std::string basename, const std::vector<TH1*>& candidates, const std::pair<unsigned,unsigned> range)
@@ -128,12 +129,12 @@ int SymmetryFit::fit() {
       p[0] > chi2_->GetBinCenter(chi2_->GetNbinsX()))
     return 7;
 
-  TF1* f = fitfunction();
+  std::unique_ptr<TF1> f( fitfunction() );
   f->SetParameter(0, p[0]);  f->SetParLimits(0, p[0], p[0]);
   f->SetParameter(1, p[1]);  f->SetParLimits(1, p[1], p[1]);
   f->SetParameter(2, p[2]);  f->SetParLimits(2, p[2], p[2]);
   f->SetParameter(3, ndf_);  f->SetParLimits(3, ndf_,ndf_); //Fixed
-  chi2_->Fit(f,"WQ");
+  chi2_->Fit(f.get(),"WQ");
   return 0;
 }
 
@@ -151,12 +152,14 @@ TF1* SymmetryFit::fitfunction()
 std::vector<double> SymmetryFit::pol2_from_pol2(TH1* hist) {
   std::vector<double> v;
 
-  int status = hist->Fit("pol2","WQ");
+  //Need our own copy for thread safety
+  TF1 func("mypol2","pol2"); 
+  int status = hist->Fit(&func,"WQ");
   if(!status) {
     std::vector<double> p;
-    p.push_back(hist->GetFunction("pol2")->GetParameter(0));
-    p.push_back(hist->GetFunction("pol2")->GetParameter(1));
-    p.push_back(hist->GetFunction("pol2")->GetParameter(2));
+    p.push_back(func.GetParameter(0));
+    p.push_back(func.GetParameter(1));
+    p.push_back(func.GetParameter(2));
     if(p[2]>0) {
       v.push_back( -0.5*p[1]/p[2] );
       v.push_back( 1./sqrt(p[2]) );
@@ -169,13 +172,14 @@ std::vector<double> SymmetryFit::pol2_from_pol2(TH1* hist) {
 std::vector<double> SymmetryFit::pol2_from_pol3(TH1* hist) {
   std::vector<double> v;
 
-  int status = hist->Fit("pol3","WQ");
+  std::unique_ptr<TF1> func( new TF1("mypol3","pol3")); 
+  int status = hist->Fit(func.get(),"WQ");
   if(!status) {
     std::vector<double> p;
-    p.push_back(hist->GetFunction("pol3")->GetParameter(0));
-    p.push_back(hist->GetFunction("pol3")->GetParameter(1));
-    p.push_back(hist->GetFunction("pol3")->GetParameter(2));
-    p.push_back(hist->GetFunction("pol3")->GetParameter(3));
+    p.push_back(func->GetParameter(0));
+    p.push_back(func->GetParameter(1));
+    p.push_back(func->GetParameter(2));
+    p.push_back(func->GetParameter(3));
     double radical = p[2]*p[2] - 3*p[1]*p[3] ;
     if(radical>0) {
       double x0 = ( -p[2] + sqrt(radical) ) / ( 3*p[3] ) ;
diff --git a/Calibration/EcalCalibAlgos/src/PhiSymmetryCalibration.cc b/Calibration/EcalCalibAlgos/src/PhiSymmetryCalibration.cc
index 2ced788edc95a..9af4b73a1a802 100644
--- a/Calibration/EcalCalibAlgos/src/PhiSymmetryCalibration.cc
+++ b/Calibration/EcalCalibAlgos/src/PhiSymmetryCalibration.cc
@@ -464,6 +464,8 @@ void PhiSymmetryCalibration::getKfactors()
   std::vector<TGraph*>  k_barl_graph(kBarlRings);
   std::vector<TCanvas*> k_barl_plot(kBarlRings);
 
+  //Create our own TF1 to avoid threading problems
+  TF1 mypol1("mypol1","pol1");
   for (int ieta=0; ieta<kBarlRings; ieta++) {
     for (int imiscal=0; imiscal<kNMiscalBinsEB; imiscal++) {
       int middlebin =  int (kNMiscalBinsEB/2);
@@ -471,8 +473,7 @@ void PhiSymmetryCalibration::getKfactors()
       epsilon_M_eb[imiscal] = miscalEB_[imiscal] - 1.;
     }
     k_barl_graph[ieta] = new TGraph (kNMiscalBinsEB,epsilon_M_eb,epsilon_T_eb);
-    k_barl_graph[ieta]->Fit("pol1");
-
+    k_barl_graph[ieta]->Fit(&mypol1);
 
     ostringstream t;
     t<< "k_barl_" << ieta+1; 
@@ -504,7 +505,7 @@ void PhiSymmetryCalibration::getKfactors()
       epsilon_M_ee[imiscal] = miscalEE_[imiscal] - 1.;
     }
     k_endc_graph[ring] = new TGraph (kNMiscalBinsEE,epsilon_M_ee,epsilon_T_ee);
-    k_endc_graph[ring]->Fit("pol1");
+    k_endc_graph[ring]->Fit(&mypol1);
 
     ostringstream t;
     t<< "k_endc_"<< ring+1;
diff --git a/Calibration/Tools/src/ZIterativeAlgorithmWithFit.cc b/Calibration/Tools/src/ZIterativeAlgorithmWithFit.cc
index aea0038a3dd36..f55f0a0e6fb46 100644
--- a/Calibration/Tools/src/ZIterativeAlgorithmWithFit.cc
+++ b/Calibration/Tools/src/ZIterativeAlgorithmWithFit.cc
@@ -365,11 +365,11 @@ ZIterativeAlgorithmWithFit::~ZIterativeAlgorithmWithFit()
 }
 
 void ZIterativeAlgorithmWithFit::gausfit(TH1F * histoou,double* par,double* errpar,float nsigmalow, float nsigmaup, double* myChi2, int* iterations) {
-  TF1 *gausa = new TF1("gausa","gaus",histoou->GetMean()-3*histoou->GetRMS(),histoou->GetMean()+3*histoou->GetRMS());
+  std::unique_ptr<TF1> gausa{ new TF1("gausa","gaus",histoou->GetMean()-3*histoou->GetRMS(),histoou->GetMean()+3*histoou->GetRMS()) };
   
   gausa->SetParameters(histoou->GetMaximum(),histoou->GetMean(),histoou->GetRMS());
   
-  histoou->Fit("gausa","qR0N");
+  histoou->Fit(gausa.get(),"qR0N");
   
   double p1    = gausa->GetParameter(1);
   double sigma = gausa->GetParameter(2);
@@ -383,7 +383,6 @@ void ZIterativeAlgorithmWithFit::gausfit(TH1F * histoou,double* par,double* errp
   double xmax_fit=p1+nsigmaup*sigma;
   
   int iter=0;
-  TF1* fitFunc;
 
   while ((chi2>1. && iter<5) || iter<2 )
     {
@@ -394,12 +393,12 @@ void ZIterativeAlgorithmWithFit::gausfit(TH1F * histoou,double* par,double* errp
       
       char suffix[20];
       sprintf (suffix,"_iter_%d",iter); 
-      fitFunc = new TF1("FitFunc"+TString(suffix),"gaus",xmin_fit,xmax_fit);
+      std::unique_ptr<TF1> fitFunc{ new TF1("FitFunc"+TString(suffix),"gaus",xmin_fit,xmax_fit) };
       fitFunc->SetParameters(nor,p1,sigma);
       fitFunc->SetLineColor((int)(iter+1));
       fitFunc->SetLineWidth(1);
       //histoou->Fit("FitFunc","lR+","");
-      histoou->Fit("FitFunc"+TString(suffix),"qR0+","");
+      histoou->Fit(fitFunc.get(),"qR0+","");
       
       histoou->GetXaxis()->SetRangeUser(xmi,xma);
       histoou->GetXaxis()->SetLabelSize(0.055);
diff --git a/DQM/BeamMonitor/plugins/BeamMonitor.cc b/DQM/BeamMonitor/plugins/BeamMonitor.cc
index d672842e3a6c5..6d5f179aea59d 100644
--- a/DQM/BeamMonitor/plugins/BeamMonitor.cc
+++ b/DQM/BeamMonitor/plugins/BeamMonitor.cc
@@ -33,6 +33,7 @@ V00-03-25
 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
 #include <numeric>
 #include <math.h>
+#include <memory>
 #include <TMath.h>
 #include <iostream>
 #include <TStyle.h>
@@ -854,10 +855,10 @@ void BeamMonitor::FitAndFill(const LuminosityBlock& lumiSeg,int &lastlumi,int &n
       sprintf(tmpTitle,"%s %i %s %i","Fitted Primary Vertex (cm) of LS: ",beginLumiOfPVFit_," to ",endLumiOfPVFit_);
       pvResults->setAxisTitle(tmpTitle,1);
 
-      TF1 *fgaus = new TF1("fgaus","gaus");
+      std::unique_ptr<TF1> fgaus{ new TF1("fgaus","gaus") };
       double mean,width,meanErr,widthErr;
       fgaus->SetLineColor(4);
-      h_PVx[0]->getTH1()->Fit("fgaus","QLM0");
+      h_PVx[0]->getTH1()->Fit(fgaus.get(),"QLM0");
       mean = fgaus->GetParameter(1);
       width = fgaus->GetParameter(2);
       meanErr = fgaus->GetParError(1);
@@ -892,10 +893,10 @@ void BeamMonitor::FitAndFill(const LuminosityBlock& lumiSeg,int &lastlumi,int &n
       tmphisto = static_cast<TH1D *>((h_PVx[0]->getTH1())->Clone("tmphisto"));
       h_PVx[1]->getTH1()->SetBins(tmphisto->GetNbinsX(),tmphisto->GetXaxis()->GetXmin(),tmphisto->GetXaxis()->GetXmax());
       h_PVx[1] = dbe_->book1D(tmpfile,h_PVx[0]->getTH1F());
-      h_PVx[1]->getTH1()->Fit("fgaus","QLM");
+      h_PVx[1]->getTH1()->Fit(fgaus.get(),"QLM");
 
 
-      h_PVy[0]->getTH1()->Fit("fgaus","QLM0");
+      h_PVy[0]->getTH1()->Fit(fgaus.get(),"QLM0");
       mean = fgaus->GetParameter(1);
       width = fgaus->GetParameter(2);
       meanErr = fgaus->GetParError(1);
@@ -921,10 +922,10 @@ void BeamMonitor::FitAndFill(const LuminosityBlock& lumiSeg,int &lastlumi,int &n
       h_PVy[1]->getTH1()->SetBins(tmphisto->GetNbinsX(),tmphisto->GetXaxis()->GetXmin(),tmphisto->GetXaxis()->GetXmax());
       h_PVy[1]->update();
       h_PVy[1] = dbe_->book1D(tmpfile,h_PVy[0]->getTH1F());
-      h_PVy[1]->getTH1()->Fit("fgaus","QLM");
+      h_PVy[1]->getTH1()->Fit(fgaus.get(),"QLM");
 
 
-      h_PVz[0]->getTH1()->Fit("fgaus","QLM0");
+      h_PVz[0]->getTH1()->Fit(fgaus.get(),"QLM0");
       mean = fgaus->GetParameter(1);
       width = fgaus->GetParameter(2);
       meanErr = fgaus->GetParError(1);
@@ -950,7 +951,7 @@ void BeamMonitor::FitAndFill(const LuminosityBlock& lumiSeg,int &lastlumi,int &n
       h_PVz[1]->getTH1()->SetBins(tmphisto->GetNbinsX(),tmphisto->GetXaxis()->GetXmin(),tmphisto->GetXaxis()->GetXmax());
       h_PVz[1]->update();
       h_PVz[1] = dbe_->book1D(tmpfile,h_PVz[0]->getTH1F());
-      h_PVz[1]->getTH1()->Fit("fgaus","QLM");
+      h_PVz[1]->getTH1()->Fit(fgaus.get(),"QLM");
 
     }//check if found min Vertices
   }//do PVfit
@@ -1165,10 +1166,10 @@ void BeamMonitor::FitAndFill(const LuminosityBlock& lumiSeg,int &lastlumi,int &n
 
 	double mean = bs.z0();
 	double width = bs.sigmaZ();
-	TF1 *fgaus = new TF1("fgaus","gaus");
+	std::unique_ptr<TF1> fgaus{ new TF1("fgaus","gaus") };
 	fgaus->SetParameters(mean,width);
 	fgaus->SetLineColor(4);
-	h_trk_z0->getTH1()->Fit("fgaus","QLRM","",mean-3*width,mean+3*width);
+	h_trk_z0->getTH1()->Fit(fgaus.get(),"QLRM","",mean-3*width,mean+3*width);
       }
 
       fitResults->Reset();
diff --git a/DQM/BeamMonitor/plugins/Vx3DHLTAnalyzer.cc b/DQM/BeamMonitor/plugins/Vx3DHLTAnalyzer.cc
index cfc4d455dc889..0021b9ac82f03 100644
--- a/DQM/BeamMonitor/plugins/Vx3DHLTAnalyzer.cc
+++ b/DQM/BeamMonitor/plugins/Vx3DHLTAnalyzer.cc
@@ -996,47 +996,47 @@ void Vx3DHLTAnalyzer::endLuminosityBlock(const LuminosityBlock& lumiBlock,
       mXlumi->ShiftFillLast(vals[0], std::sqrt(vals[8]), nLumiReset);
       myLinFit->SetParameter(0, mXlumi->getTH1()->GetMean(2));
       myLinFit->SetParameter(1, 0.0);
-      mXlumi->getTH1()->Fit("myLinFit","QR");
+      mXlumi->getTH1()->Fit(myLinFit,"QR");
 
       mYlumi->ShiftFillLast(vals[1], std::sqrt(vals[9]), nLumiReset);
       myLinFit->SetParameter(0, mYlumi->getTH1()->GetMean(2));
       myLinFit->SetParameter(1, 0.0);
-      mYlumi->getTH1()->Fit("myLinFit","QR");
+      mYlumi->getTH1()->Fit(myLinFit,"QR");
 
       mZlumi->ShiftFillLast(vals[2], std::sqrt(vals[10]), nLumiReset);
       myLinFit->SetParameter(0, mZlumi->getTH1()->GetMean(2));
       myLinFit->SetParameter(1, 0.0);
-      mZlumi->getTH1()->Fit("myLinFit","QR");
+      mZlumi->getTH1()->Fit(myLinFit,"QR");
 
       sXlumi->ShiftFillLast(vals[6], std::sqrt(vals[14]), nLumiReset);
       myLinFit->SetParameter(0, sXlumi->getTH1()->GetMean(2));
       myLinFit->SetParameter(1, 0.0);
-      sXlumi->getTH1()->Fit("myLinFit","QR");
+      sXlumi->getTH1()->Fit(myLinFit,"QR");
 
       sYlumi->ShiftFillLast(vals[7], std::sqrt(vals[15]), nLumiReset);
       myLinFit->SetParameter(0, sYlumi->getTH1()->GetMean(2));
       myLinFit->SetParameter(1, 0.0);
-      sYlumi->getTH1()->Fit("myLinFit","QR");
+      sYlumi->getTH1()->Fit(myLinFit,"QR");
 
       sZlumi->ShiftFillLast(vals[3], std::sqrt(vals[11]), nLumiReset);
       myLinFit->SetParameter(0, sZlumi->getTH1()->GetMean(2));
       myLinFit->SetParameter(1, 0.0);
-      sZlumi->getTH1()->Fit("myLinFit","QR");
+      sZlumi->getTH1()->Fit(myLinFit,"QR");
 
       dxdzlumi->ShiftFillLast(vals[4], std::sqrt(vals[12]), nLumiReset);
       myLinFit->SetParameter(0, dxdzlumi->getTH1()->GetMean(2));
       myLinFit->SetParameter(1, 0.0);
-      dxdzlumi->getTH1()->Fit("myLinFit","QR");
+      dxdzlumi->getTH1()->Fit(myLinFit,"QR");
 
       dydzlumi->ShiftFillLast(vals[5], std::sqrt(vals[13]), nLumiReset);
       myLinFit->SetParameter(0, dydzlumi->getTH1()->GetMean(2));
       myLinFit->SetParameter(1, 0.0);
-      dydzlumi->getTH1()->Fit("myLinFit","QR");
+      dydzlumi->getTH1()->Fit(myLinFit,"QR");
       
       goodVxCounter->ShiftFillLast((double)counterVx, std::sqrt((double)counterVx), nLumiReset);      
       myLinFit->SetParameter(0, goodVxCounter->getTH1()->GetMean(2));
       myLinFit->SetParameter(1, 0.0);
-      goodVxCounter->getTH1()->Fit("myLinFit","QR");
+      goodVxCounter->getTH1()->Fit(myLinFit,"QR");
 
       if (lastLumiOfFit % prescaleHistory == 0)
 	{
diff --git a/DQM/DTMonitorClient/src/DTLocalTriggerLutTest.cc b/DQM/DTMonitorClient/src/DTLocalTriggerLutTest.cc
index 43db640ae6b4e..47eeca795816c 100644
--- a/DQM/DTMonitorClient/src/DTLocalTriggerLutTest.cc
+++ b/DQM/DTMonitorClient/src/DTLocalTriggerLutTest.cc
@@ -144,13 +144,11 @@ void DTLocalTriggerLutTest::runClientDiagnostic() {
 	    double phiSlope = 0;
 	    double phiCorr  = 0;
 	    try {
-	      PhitkvsPhitrigProf->Fit("pol1","CQO");
-	      TF1 *ffPhi= PhitkvsPhitrigProf->GetFunction("pol1");
-	      if (ffPhi) {
-		phiInt   = ffPhi->GetParameter(0);
-		phiSlope = ffPhi->GetParameter(1);
-		phiCorr  = TrackPhitkvsPhitrig->GetCorrelationFactor();
-	      }
+	      TF1 ffPhi("mypol1","pol1");
+	      PhitkvsPhitrigProf->Fit(&ffPhi,"CQO");
+	      phiInt   = ffPhi.GetParameter(0);
+	      phiSlope = ffPhi.GetParameter(1);
+	      phiCorr  = TrackPhitkvsPhitrig->GetCorrelationFactor();
 	    } catch (cms::Exception& iException) {
 	      edm::LogError(category()) << "[" << testName << "Test]: Error fitting PhitkvsPhitrig for Wheel " << wh 
 					<<" Sector " << sect << " Station " << stat;
@@ -180,13 +178,11 @@ void DTLocalTriggerLutTest::runClientDiagnostic() {
 	    double phibSlope = 0;
 	    double phibCorr  = 0;
 	    try {
-	      PhibtkvsPhibtrigProf->Fit("pol1","CQO");
-	      TF1 *ffPhib= PhibtkvsPhibtrigProf->GetFunction("pol1");
-	      if (ffPhib) {
-		phibInt   = ffPhib->GetParameter(0);
-		phibSlope = ffPhib->GetParameter(1);
-		phibCorr  = TrackPhibtkvsPhibtrig->GetCorrelationFactor();
-	      }
+	      TF1 ffPhib("ffPhib","pol1");
+	      PhibtkvsPhibtrigProf->Fit(&ffPhib,"CQO");
+	      phibInt   = ffPhib.GetParameter(0);
+	      phibSlope = ffPhib.GetParameter(1);
+	      phibCorr  = TrackPhibtkvsPhibtrig->GetCorrelationFactor();
 	    } catch (cms::Exception& iException) {
 	      edm::LogError(category()) << "[" << testName << "Test]: Error fitting PhibtkvsPhibtrig for Wheel " << wh 
 					<<" Sector " << sect << " Station " << stat;
@@ -217,12 +213,10 @@ void DTLocalTriggerLutTest::runClientDiagnostic() {
 	  double phiMean = 0;
 	  double phiRMS  = 0;
 	  try {
-	    PhiResidual->Fit("gaus","CQO","",peak-5,peak+5);
-	    TF1 *ffPhi = PhiResidual->GetFunction("gaus");
-	    if ( ffPhi ) {
-	      phiMean = ffPhi->GetParameter(1);
-	      phiRMS  = ffPhi->GetParameter(2);
-	    }
+	    TF1 ffPhi("ffPhi","gaus");
+	    PhiResidual->Fit(&ffPhi,"CQO","",peak-5,peak+5);
+	    phiMean = ffPhi.GetParameter(1);
+	    phiRMS  = ffPhi.GetParameter(2);
 	  } catch (cms::Exception& iException) {
 	    edm::LogError(category()) << "[" << testName << "Test]: Error fitting PhiResidual for Wheel " << wh 
 				      <<" Sector " << sect << " Station " << stat;
@@ -253,12 +247,10 @@ void DTLocalTriggerLutTest::runClientDiagnostic() {
 	  double phibMean = 0;
 	  double phibRMS  = 0;
 	  try {
-	    PhibResidual->Fit("gaus","CQO","",peak-5,peak+5);
-	    TF1 *ffPhib = PhibResidual->GetFunction("gaus");
-	    if ( ffPhib ) {
-	      phibMean = ffPhib->GetParameter(1);
-	      phibRMS  = ffPhib->GetParameter(2);
-	    }
+	    TF1 ffPhib("ffPhib","gaus");
+	    PhibResidual->Fit(&ffPhib,"CQO","",peak-5,peak+5);
+	    phibMean = ffPhib.GetParameter(1);
+	    phibRMS  = ffPhib.GetParameter(2);
 	  } catch (cms::Exception& iException) {
 	    edm::LogError(category()) << "[" << testName << "Test]: Error fitting PhibResidual for Wheel " << wh 
 				      <<" Sector " << sect << " Station " << stat;
diff --git a/DQM/DTMonitorClient/src/DTLocalTriggerSynchTest.cc b/DQM/DTMonitorClient/src/DTLocalTriggerSynchTest.cc
index cc43dacd97664..5c25c3cf403e8 100644
--- a/DQM/DTMonitorClient/src/DTLocalTriggerSynchTest.cc
+++ b/DQM/DTMonitorClient/src/DTLocalTriggerSynchTest.cc
@@ -125,7 +125,9 @@ void DTLocalTriggerSynchTest::runClientDiagnostic() {
 	  MonitorElement* ratioH = innerME.find(fullName(ratioHistoTag))->second;
 	  makeRatioME(numH,denH,ratioH);
 	  try {
-	    getHisto<TH1F>(ratioH)->Fit("pol8","CQO");
+	    //Need our own copy to avoid threading problems
+	    TF1 mypol8("mypol8","pol8");
+	    getHisto<TH1F>(ratioH)->Fit(&mypol8,"CQO");
 	  } catch (cms::Exception& iException) {
 	    edm::LogPrint(category()) << "[" << testName 
 				     << "Test]: Error fitting " 
@@ -174,7 +176,7 @@ void DTLocalTriggerSynchTest::endJob(){
 
 	TH1F *ratioH     = getHisto<TH1F>(dbe->get(getMEName(ratioHistoTag,"", chId)));    
 	if (ratioH->GetEntries()>minEntries) {	      
-	  TF1 *fitF=ratioH->GetFunction("pol8");
+	  TF1 *fitF=ratioH->GetFunction("mypol8");
 	  if (fitF) { fineDelay=fitF->GetMaximumX(0,bxTime); }
 	} else {
 	  LogInfo(category()) << "[" << testName 
diff --git a/DQM/DTMonitorClient/src/DTResolutionTest.cc b/DQM/DTMonitorClient/src/DTResolutionTest.cc
index 20aa3effb6fae..c3fcc9ce91017 100644
--- a/DQM/DTMonitorClient/src/DTResolutionTest.cc
+++ b/DQM/DTMonitorClient/src/DTResolutionTest.cc
@@ -256,8 +256,10 @@ void DTResolutionTest::endLuminosityBlock(LuminosityBlock const& lumiSeg, EventS
 	  TProfile* prof = res_histo_2D_root->ProfileX();
 	  prof->GetXaxis()->SetRangeUser(0,2);
 	  //prof->Fit("pol1","Q0");
+	  //need our own copy to avoid threading problems
+	  TF1 fitting("mypol1","pol1");
 	  try {
-	    prof->Fit("pol1","Q0");
+	    prof->Fit(&fitting,"Q0");
 	  } catch (cms::Exception& iException) {
 	    edm::LogError ("resolution") << "[DTResolutionTest]: Exception when fitting..."
 					 << "SuperLayer : " << slID << "\n"
@@ -266,8 +268,7 @@ void DTResolutionTest::endLuminosityBlock(LuminosityBlock const& lumiSeg, EventS
 	    SlopeHistos.find(make_pair(slID.wheel(),slID.sector()))->second->setBinContent(BinNumber, -99.);
 	    continue;
 	  }
-	  TF1 *fitting = prof->GetFunction("pol1");
-	  double slope = fitting->GetParameter(1);
+	  double slope = fitting.GetParameter(1);
 	  SlopeHistos.find(make_pair(slID.wheel(),slID.sector()))->second->setBinContent(BinNumber, slope);	
 	}
       }
diff --git a/DQM/EcalPreshowerMonitorClient/src/ESPedestalClient.cc b/DQM/EcalPreshowerMonitorClient/src/ESPedestalClient.cc
index fcc23dc28da26..5e7c3d6084a03 100644
--- a/DQM/EcalPreshowerMonitorClient/src/ESPedestalClient.cc
+++ b/DQM/EcalPreshowerMonitorClient/src/ESPedestalClient.cc
@@ -14,7 +14,8 @@
 using namespace edm;
 using namespace std;
 
-ESPedestalClient::ESPedestalClient(const edm::ParameterSet& ps) {
+ESPedestalClient::ESPedestalClient(const edm::ParameterSet& ps):
+  fg(nullptr) {
   
   verbose_       = ps.getUntrackedParameter<bool>("verbose", true);
   debug_         = ps.getUntrackedParameter<bool>("debug", true);
@@ -34,6 +35,7 @@ ESPedestalClient::ESPedestalClient(const edm::ParameterSet& ps) {
 }
 
 ESPedestalClient::~ESPedestalClient() {
+  delete fg;
 }
 
 void ESPedestalClient::beginJob(DQMStore* dqmStore) {
@@ -79,8 +81,8 @@ void ESPedestalClient::endJob(void) {
 
 	    if (meFit==0) continue;
 	    TH1F *rootHisto = meFit->getTH1F();
-	    rootHisto->Fit("fg", "Q", "", 500, 1800);
-	    rootHisto->Fit("fg", "RQ", "", fg->GetParameter(1)-2.*fg->GetParameter(2),fg->GetParameter(1)+2.*fg->GetParameter(2));
+	    rootHisto->Fit(fg, "Q", "", 500, 1800);
+	    rootHisto->Fit(fg, "RQ", "", fg->GetParameter(1)-2.*fg->GetParameter(2),fg->GetParameter(1)+2.*fg->GetParameter(2));
 	    hPed_[iz][senP_[i]-1][senX_[i]-1][senY_[i]-1]->setBinContent(is+1, (int)(fg->GetParameter(1)+0.5));
 	    hTotN_[iz][senP_[i]-1][senX_[i]-1][senY_[i]-1]->setBinContent(is+1, fg->GetParameter(2));
 
diff --git a/DQM/HcalMonitorTasks/src/HcalDetDiagLaserMonitor.cc b/DQM/HcalMonitorTasks/src/HcalDetDiagLaserMonitor.cc
index d6ba4ead819a1..7478c914c51bf 100644
--- a/DQM/HcalMonitorTasks/src/HcalDetDiagLaserMonitor.cc
+++ b/DQM/HcalMonitorTasks/src/HcalDetDiagLaserMonitor.cc
@@ -1907,7 +1907,7 @@ char str[100];
             S2[i]->Fill(adc2fC[j]-2.5,Raddam_data[i].s2_adc[j]); 
           }
           double parm[3];
-          S1[i]->Fit("fitFunc");
+          S1[i]->Fit(fitFunc);
           S1[i]->GetFunction("fitFunc")->GetParameters(parm);
           Raddam_data[i].S1MEAN=S1[i]->GetMean();
           Raddam_data[i].S1RMS=S1[i]->GetRMS();
@@ -1917,7 +1917,7 @@ char str[100];
           Raddam_data[i].S1CHI2=S1[i]->GetFunction("fitFunc")->GetChisquare();
           Raddam_data[i].S1NDF=S1[i]->GetFunction("fitFunc")->GetNDF();
           Raddam_data[i].S1BINWIDTH=Ws1;
-          S2[i]->Fit("fitFunc");
+          S2[i]->Fit(fitFunc);
           S2[i]->GetFunction("fitFunc")->GetParameters(parm);
           Raddam_data[i].S2MEAN=S2[i]->GetMean();
           Raddam_data[i].S2RMS=S2[i]->GetRMS();
diff --git a/DQM/SiPixelHistoricInfoClient/plugins/SiPixelHistoryDQMService.cc b/DQM/SiPixelHistoricInfoClient/plugins/SiPixelHistoryDQMService.cc
index a091b1bfb54cd..ebbc937b92706 100644
--- a/DQM/SiPixelHistoricInfoClient/plugins/SiPixelHistoryDQMService.cc
+++ b/DQM/SiPixelHistoricInfoClient/plugins/SiPixelHistoryDQMService.cc
@@ -86,10 +86,10 @@ bool SiPixelHistoryDQMService::setDBValuesForUser(std::vector<MonitorElement*>::
     // }
     // else
     // if(  ) {
-    Hist->Fit("pol0");
-    TF1* Fit = Hist->GetFunction("pol0");
-    float FitValue = Fit ? Fit->GetParameter(0) : 0;
-    float FitError = Fit ? Fit->GetParError(0) : 0;
+    TF1 Fit("mypol0","pol0");
+    Hist->Fit(&Fit);
+    float FitValue = Fit.GetParameter(0) ;
+    float FitError = Fit.GetParError(0) ;
     std::cout << "FITERROR: " << FitError << std::endl;
 
     values.push_back( FitValue );
diff --git a/DQM/SiStripCommissioningAnalysis/src/PedsFullNoiseAlgorithm.cc b/DQM/SiStripCommissioningAnalysis/src/PedsFullNoiseAlgorithm.cc
index 35faf469b0356..23e3f122832b5 100644
--- a/DQM/SiStripCommissioningAnalysis/src/PedsFullNoiseAlgorithm.cc
+++ b/DQM/SiStripCommissioningAnalysis/src/PedsFullNoiseAlgorithm.cc
@@ -152,6 +152,9 @@ void PedsFullNoiseAlgorithm::analyse() {
     float n_sum = 0., n_sum2 = 0., n_max = -1.*sistrip::invalid_, n_min = sistrip::invalid_;
     float r_sum = 0., r_sum2 = 0., r_max = -1.*sistrip::invalid_, r_min = sistrip::invalid_;		   
     // Iterate through strips of APV
+
+    //Need our own copy for thread safety
+    TF1 mygaus("mygaus","gaus");
     for ( uint16_t istr = 0; istr < 128; istr++ ) {
     	
       ana->ksProb_[iapv].push_back(0);
@@ -189,8 +192,8 @@ void PedsFullNoiseAlgorithm::analyse() {
         // If the histogram is valid.
         if(noisehist->Integral() > 0){
         	ana->noiseRMS_[iapv][istr] = noisehist->GetRMS();
-        	noisehist->Fit("gaus","Q");                       
-        	ana->noiseGaus_[iapv][istr] 	= noisehist->GetFunction("gaus")->GetParameter(2);
+        	noisehist->Fit(&mygaus,"Q");                       
+        	ana->noiseGaus_[iapv][istr] 	= mygaus.GetParameter(2);
           
           // new Bin84 method
       		std::vector<float> integralFrac;
@@ -222,7 +225,7 @@ void PedsFullNoiseAlgorithm::analyse() {
         	// Compare shape of ADC to a histogram made of Gaus Fit for KSProb, Chi2Prob, Etc...
         	TH1D * FitHisto = new TH1D("FitHisto","FitHisto",noisehist->GetNbinsX(),
                            -noisehist->GetNbinsX()/2,noisehist->GetNbinsX()/2);
-        	FitHisto->Add(noisehist->GetFunction("gaus"));
+        	FitHisto->Add(&mygaus);
         	FitHisto->Sumw2();
         	noisehist->Sumw2();
           
diff --git a/DQMOffline/JetMET/plugins/SusyPostProcessor.cc b/DQMOffline/JetMET/plugins/SusyPostProcessor.cc
index 37022cc274b16..c7c4061d93ab0 100644
--- a/DQMOffline/JetMET/plugins/SusyPostProcessor.cc
+++ b/DQMOffline/JetMET/plugins/SusyPostProcessor.cc
@@ -67,6 +67,9 @@ void SusyPostProcessor::dqmEndJob(DQMStore::IBooker& ibook_, DQMStore::IGetter&
   metFolders.push_back("Uncleaned/");
   metFolders.push_back("Cleaned/");
 
+  //Need our own copy for thread safety
+  TF1 mygaus("mygaus","gaus");
+
   for (int i=0; i<int(Dirs.size()); i++) {
 
     std::string prefix = "dummy";
@@ -83,11 +86,11 @@ void SusyPostProcessor::dqmEndJob(DQMStore::IBooker& ibook_, DQMStore::IGetter&
       MEy = iget_.get(dirName + "/"+"MEy");
 
       if (MEx && MEx->kind() == MonitorElement::DQM_KIND_TH1F) {
-	if (MEx->getTH1F()->GetEntries() > 50) MEx->getTH1F()->Fit("gaus", "q");
+	if (MEx->getTH1F()->GetEntries() > 50) MEx->getTH1F()->Fit(&mygaus, "q");
       }
 
       if (MEy && MEy->kind() == MonitorElement::DQM_KIND_TH1F) {
-	if (MEy->getTH1F()->GetEntries() > 50) MEy->getTH1F()->Fit("gaus", "q");
+	if (MEy->getTH1F()->GetEntries() > 50) MEy->getTH1F()->Fit(&mygaus, "q");
       }
     }
   }
diff --git a/DQMOffline/PFTau/plugins/PFClient.cc b/DQMOffline/PFTau/plugins/PFClient.cc
index d02524768a443..98f6c12ec104f 100644
--- a/DQMOffline/PFTau/plugins/PFClient.cc
+++ b/DQMOffline/PFTau/plugins/PFClient.cc
@@ -288,12 +288,11 @@ void PFClient::getHistogramParameters(MonitorElement* me_slice, double& average,
     rms     = me_slice->getRMS();  
     TH1F* th_slice = me_slice->getTH1F();
     if (th_slice && th_slice->GetEntries() > 0) {
-      th_slice->Fit( "gaus","Q0");
-      TF1* gaus = th_slice->GetFunction( "gaus" );
-      if (gaus) {
-	sigma = gaus->GetParameter(2);
-        mean  = gaus->GetParameter(1);
-      }
+      //need our own copy for thread safety
+      TF1 gaus("mygaus","gaus");
+      th_slice->Fit( &gaus,"Q0");
+      sigma = gaus.GetParameter(2);
+      mean  = gaus.GetParameter(1);
     }
   }
 }
diff --git a/DQMOffline/PFTau/plugins/PFClient_JetRes.cc b/DQMOffline/PFTau/plugins/PFClient_JetRes.cc
index 4e3ecca585abe..ddb521d68fa4b 100644
--- a/DQMOffline/PFTau/plugins/PFClient_JetRes.cc
+++ b/DQMOffline/PFTau/plugins/PFClient_JetRes.cc
@@ -185,12 +185,11 @@ void PFClient_JetRes::getHistogramParameters(MonitorElement* me_slice, double& a
     rms     = me_slice->getRMS();  
     TH1F* th_slice = me_slice->getTH1F();
     if (th_slice && th_slice->GetEntries() > 0) {
-      th_slice->Fit( "gaus","Q0");
-      TF1* gaus = th_slice->GetFunction( "gaus" );
-      if (gaus) {
-	sigma = gaus->GetParameter(2);
-        mean  = gaus->GetParameter(1);
-      }
+      //need our own copy for thread safety
+      TF1 gaus("mygaus","gaus");
+      th_slice->Fit( &gaus,"Q0");
+      sigma = gaus.GetParameter(2);
+      mean  = gaus.GetParameter(1);
     }
   }
 }
diff --git a/DQMOffline/RecoB/src/BTagDifferentialPlot.cc b/DQMOffline/RecoB/src/BTagDifferentialPlot.cc
index d79dc30be23a9..abaa2f8253352 100644
--- a/DQMOffline/RecoB/src/BTagDifferentialPlot.cc
+++ b/DQMOffline/RecoB/src/BTagDifferentialPlot.cc
@@ -334,8 +334,10 @@ BTagDifferentialPlot::getMistag(const double& fixedBEfficiency, TH1F * effPurHis
     // is above or below.
     // Fit a plynomial, and evaluate the mistag at the requested value.
     int fitStatus;
+    //need our own copy for thread safety
+    TF1 myfunc("myfunc","pol4");
     try {
-      fitStatus = effPurHist->Fit("pol4", "q");
+      fitStatus = effPurHist->Fit(&myfunc, "q");
     }catch (cms::Exception& iException){
       return pair<double, double>(effForBEff, effForBEffErr);
     }
@@ -345,9 +347,8 @@ BTagDifferentialPlot::getMistag(const double& fixedBEfficiency, TH1F * effPurHis
       //      edm::LogInfo("BTagDifferentialPlot")<<"Fit OK to hisogram " << effPurHist->GetTitle() << " entries = " << effPurHist->GetEntries();
       return pair<double, double>(effForBEff, effForBEffErr);
     }
-    TF1 *myfunc = effPurHist->GetFunction("pol4");
-    effForBEff = myfunc->Eval(fixedBEfficiency);
-    effPurHist->RecursiveRemove(myfunc);
+    effForBEff = myfunc.Eval(fixedBEfficiency);
+    effPurHist->RecursiveRemove(&myfunc);
     //Error: first non-empty bin on the right and take its error
     for (int i = iBinGet+1; i< effPurHist->GetNbinsX(); ++i) {
       if (effPurHist->GetBinContent(i)!=0) {
diff --git a/DQMServices/Examples/src/DQMClientExample.cc b/DQMServices/Examples/src/DQMClientExample.cc
index cd2d359ac5691..3b32a538ccafb 100644
--- a/DQMServices/Examples/src/DQMClientExample.cc
+++ b/DQMServices/Examples/src/DQMClientExample.cc
@@ -168,7 +168,7 @@ void DQMClientExample::performClient(){
 
       ////---- make the Gaussian fit to this histogram
       TF1 *f1 = new TF1("f1","gaus",1,3);
-      rootHisto->Fit("f1");
+      rootHisto->Fit(f1);
       mean = f1->GetParameter(1);
       rms = f1->GetParameter(2);
     }
diff --git a/RecoHI/HiTracking/plugins/HIPixelMedianVtxProducer.cc b/RecoHI/HiTracking/plugins/HIPixelMedianVtxProducer.cc
index 4da085b30912b..2ed02936bcfe6 100644
--- a/RecoHI/HiTracking/plugins/HIPixelMedianVtxProducer.cc
+++ b/RecoHI/HiTracking/plugins/HIPixelMedianVtxProducer.cc
@@ -156,7 +156,7 @@ void HIPixelMedianVtxProducer::produce
   f1.SetParLimits(2,0.001,0.05);
   f1.SetParLimits(3,0.0,0.005*tracks.size());
     
-  histo.Fit("f1","QN");
+  histo.Fit(&f1,"QN");
     
   LogTrace("MinBiasTracking")
     << "  [vertex position] fitted    = "
diff --git a/RecoLocalCalo/EcalRecAlgos/interface/EcalUncalibRecHitRecAnalFitAlgo.h b/RecoLocalCalo/EcalRecAlgos/interface/EcalUncalibRecHitRecAnalFitAlgo.h
index 2a3f86b3d3cad..c985bbeec418c 100644
--- a/RecoLocalCalo/EcalRecAlgos/interface/EcalUncalibRecHitRecAnalFitAlgo.h
+++ b/RecoLocalCalo/EcalRecAlgos/interface/EcalUncalibRecHitRecAnalFitAlgo.h
@@ -125,14 +125,14 @@ template<class C> class EcalUncalibRecHitRecAnalFitAlgo : public EcalUncalibRecH
 
 
 
-    graph.Fit("pulseShape","QRM");
+    graph.Fit(&pulseShape,"QRM");
     //TF1 *pulseShape2=graph.GetFunction("pulseShape");
 
     if ( std::string(gMinuit->fCstatu.Data()) == std::string("CONVERGED ") ) {
 
       double amplitude_value=pulseShape.GetParameter(0);
 
-      graph.Fit("pedestal","QRL");
+      graph.Fit(&pedestal,"QRL");
       //TF1 *pedestal2=graph.GetFunction("pedestal");
       double pedestal_value=pedestal.GetParameter(0);
 
diff --git a/RecoLocalCalo/EcalRecAlgos/src/ESRecHitFitAlgo.cc b/RecoLocalCalo/EcalRecAlgos/src/ESRecHitFitAlgo.cc
index c01deeb4ce70b..e01b914eb85cf 100644
--- a/RecoLocalCalo/EcalRecAlgos/src/ESRecHitFitAlgo.cc
+++ b/RecoLocalCalo/EcalRecAlgos/src/ESRecHitFitAlgo.cc
@@ -59,7 +59,7 @@ double* ESRecHitFitAlgo::EvalAmplitude(const ESDataFrame& digi, double ped) cons
     double para[10];
     TGraph *gr = new TGraph(3, tx, adc);
     fit_->SetParameters(50, 10);
-    gr->Fit("fit", "MQ");
+    gr->Fit(fit_, "MQ");
     fit_->GetParameters(para); 
     fitresults[0] = para[0]; 
     fitresults[1] = para[1];  
diff --git a/RecoPixelVertexing/PixelLowPtUtilities/plugins/PixelVertexProducerMedian.cc b/RecoPixelVertexing/PixelLowPtUtilities/plugins/PixelVertexProducerMedian.cc
index 3e25272a095d2..9c9bf52395d06 100644
--- a/RecoPixelVertexing/PixelLowPtUtilities/plugins/PixelVertexProducerMedian.cc
+++ b/RecoPixelVertexing/PixelLowPtUtilities/plugins/PixelVertexProducerMedian.cc
@@ -114,7 +114,7 @@ void PixelVertexProducerMedian::produce
     TF1 f1("f1","[0]*exp(-0.5 * ((x-[1])/[2])^2) + [3]");
     f1.SetParameters(10.,0.,0.01, 1.);
   
-    histo.Fit("f1","QN");
+    histo.Fit(&f1,"QN");
   
     LogTrace("MinBiasTracking")
               << "  [vertex position] fitted    = "
diff --git a/RecoVertex/BeamSpotProducer/src/BSFitter.cc b/RecoVertex/BeamSpotProducer/src/BSFitter.cc
index 4fbe90bff668a..244d05e24ae68 100644
--- a/RecoVertex/BeamSpotProducer/src/BSFitter.cc
+++ b/RecoVertex/BeamSpotProducer/src/BSFitter.cc
@@ -360,17 +360,18 @@ reco::BeamSpot BSFitter::Fit_z_chi2(double *inipar) {
 		 //std::cout<<"z0="<<iparam->z0()<<"; sigZ0="<<iparam->sigz0()<<std::endl;
 	}
 
-	h1z->Fit("gaus","QLM0");
+	//Use our own copy for thread safety
+	TF1 fgaus("fgaus","gaus");
+	h1z->Fit(&fgaus,"QLM0");
 	//std::cout << "fitted "<< std::endl;
 	
-	TF1 *fgaus = h1z->GetFunction("gaus");
 	//std::cout << "got function" << std::endl;
-	double fpar[2] = {fgaus->GetParameter(1), fgaus->GetParameter(2) };
+	double fpar[2] = {fgaus.GetParameter(1), fgaus.GetParameter(2) };
 	//std::cout<<"Debug fpar[2] = (" <<fpar[0]<<","<<fpar[1]<<")"<<std::endl;
 	reco::BeamSpot::CovarianceMatrix matrix;
 	// add matrix values.
-	matrix(2,2) = fgaus->GetParError(1) * fgaus->GetParError(1);
-	matrix(3,3) = fgaus->GetParError(2) * fgaus->GetParError(2);
+	matrix(2,2) = fgaus.GetParError(1) * fgaus.GetParError(1);
+	matrix(3,3) = fgaus.GetParError(2) * fgaus.GetParError(2);
 	
 	//delete h1z;
 
@@ -564,16 +565,19 @@ reco::BeamSpot BSFitter::Fit_d0phi() {
 	//LogDebug ("BSFitter") << " d0-phi fit done.";
 	//std::cout<< " d0-phi fit done." << std::endl;
 
-	h1z->Fit("gaus","QLM0","",h1z->GetMean() -2.*h1z->GetRMS(),h1z->GetMean() +2.*h1z->GetRMS());
+	//Use our own copy for thread safety
+	TF1 fgaus("fgaus","gaus");
+	//returns 0 if OK
+	auto status = h1z->Fit(&fgaus,"QLM0","",h1z->GetMean() -2.*h1z->GetRMS(),h1z->GetMean() +2.*h1z->GetRMS());
 
 	//std::cout << "fitted "<< std::endl;
-	TF1 *fgaus = h1z->GetFunction("gaus");
+
 	//std::cout << "got function" << std::endl;
-	if (!fgaus){	
+	if (status){	
 	  edm::LogError("NoBeamSpotFit")<<"gaussian fit failed. no BS d0 fit";		
 	  return reco::BeamSpot();
 	}
-	double fpar[2] = {fgaus->GetParameter(1), fgaus->GetParameter(2) };
+	double fpar[2] = {fgaus.GetParameter(1), fgaus.GetParameter(2) };
     
 	reco::BeamSpot::CovarianceMatrix matrix;
 	// first two parameters
@@ -590,8 +594,8 @@ reco::BeamSpot BSFitter::Fit_d0phi() {
 	}
 
     // Z0 and sigmaZ
-	matrix(2,2) = fgaus->GetParError(1) * fgaus->GetParError(1);
-	matrix(3,3) = fgaus->GetParError(2) * fgaus->GetParError(2);
+	matrix(2,2) = fgaus.GetParError(1) * fgaus.GetParError(1);
+	matrix(3,3) = fgaus.GetParError(2) * fgaus.GetParError(2);
     
 	ftmp = x_result;
 
diff --git a/RecoVertex/BeamSpotProducer/src/PVFitter.cc b/RecoVertex/BeamSpotProducer/src/PVFitter.cc
index d0b107e2ef12a..ed354e9439716 100644
--- a/RecoVertex/BeamSpotProducer/src/PVFitter.cc
+++ b/RecoVertex/BeamSpotProducer/src/PVFitter.cc
@@ -345,13 +345,16 @@ bool PVFitter::runFitter() {
       TH1F *h1PVy = (TH1F*) hPVy->ProjectionX("h1PVy", 0, -1, "e");
       TH1F *h1PVz = (TH1F*) hPVx->ProjectionY("h1PVz", 0, -1, "e");
 
-      h1PVx->Fit("gaus","QLM0");
-      h1PVy->Fit("gaus","QLM0");
-      h1PVz->Fit("gaus","QLM0");
+      //Use our own copy for thread safety
+      TF1 gaus("localGaus","gaus");
 
-      TF1 *gausx = h1PVx->GetFunction("gaus");
-      TF1 *gausy = h1PVy->GetFunction("gaus");
-      TF1 *gausz = h1PVz->GetFunction("gaus");
+      h1PVx->Fit(&gaus,"QLM0");
+      h1PVy->Fit(&gaus,"QLM0");
+      h1PVz->Fit(&gaus,"QLM0");
+
+      TF1 *gausx = h1PVx->GetFunction("localGaus");
+      TF1 *gausy = h1PVy->GetFunction("localGaus");
+      TF1 *gausz = h1PVz->GetFunction("localGaus");
 
       fwidthX = gausx->GetParameter(2);
       fwidthY = gausy->GetParameter(2);
diff --git a/SimG4Core/GFlash/TB/fitMaximumPoint.cc b/SimG4Core/GFlash/TB/fitMaximumPoint.cc
index 8a48303f4a230..1e4e1ec8d921d 100644
--- a/SimG4Core/GFlash/TB/fitMaximumPoint.cc
+++ b/SimG4Core/GFlash/TB/fitMaximumPoint.cc
@@ -149,7 +149,7 @@ Float_t Fit_MaximumPoint( TH2F *h2, Int_t Ireturn = 0){  // xxx
   Double_t fitresults_errors[3]   = {0.};
   TF1 *f1 = new TF1("f1",Pol2,-50.,50.,3);  
   f1->SetParameters( 1.,1.,1.);
-  h1->Fit("f1","Q","",posmin_fit, posmax_fit);
+  h1->Fit(f1,"Q","",posmin_fit, posmax_fit);
   for(int i=0 ; i< 3 ; i++) { 
     fitresults[i]        = f1->GetParameter(i); 
     fitresults_errors[i] = f1->GetParError(i); 
@@ -208,7 +208,7 @@ Float_t Fit_MaximumPoint( TH2F *h2, Int_t Ireturn = 0){  // xxx
   Double_t fitresults0_errors[3]   = {0.};
   TF1 *f0 = new TF1("f0",Pol2_Special,-50.,50.,3);  
   f0->SetParameters( -1.,x_max_guess,5000.);
-  h0->Fit("f0","Q","",posmin_fit, posmax_fit);
+  h0->Fit(f0,"Q","",posmin_fit, posmax_fit);
   for(int i=0 ; i< 3 ; i++) { 
     fitresults0[i]        = f0->GetParameter(i); 
     fitresults0_errors[i] = f0->GetParError(i); 
diff --git a/SimG4Core/GFlash/TB/lateralAllCalib.cpp b/SimG4Core/GFlash/TB/lateralAllCalib.cpp
index 4e0f5dd648008..62429e97b83c5 100644
--- a/SimG4Core/GFlash/TB/lateralAllCalib.cpp
+++ b/SimG4Core/GFlash/TB/lateralAllCalib.cpp
@@ -1183,7 +1183,7 @@ int main ( int argc, char **argv)
       // gaussian fit to initialize
       TF1 *gausa = new TF1 ("gausa","[0]*exp(-1*(x-[1])*(x-[1])/2/[2]/[2])",h_peak-10*h_rms,h_peak+10*h_rms);
       gausa->SetParameters(h_norm,h_peak,h_rms);
-      H_eMatrix[ii][myXib]->Fit("gausa","","",h_peak-3*h_rms,h_peak+3*h_rms);
+      H_eMatrix[ii][myXib]->Fit(&gausa,"","",h_peak-3*h_rms,h_peak+3*h_rms);
       double gausNorm  = gausa->GetParameter(0);
       double gausMean  = gausa->GetParameter(1);
       double gausSigma = fabs(gausa->GetParameter(2));
@@ -1204,7 +1204,7 @@ int main ( int argc, char **argv)
       // cb_p->FixParameter(3, 3.);     // solo x caso 500 no Birk
       cb_p->SetParameter(4, gausNorm);
       cb_p->SetParLimits(2, 0.1, 5.);
-      H_eMatrix[ii][myXib]->Fit("cb_p","lR","",myXmin,myXmax);
+      H_eMatrix[ii][myXib]->Fit(&cb_p,"lR","",myXmin,myXmax);
 
       double matrix_gmean      = cb_p->GetParameter(0);
       double matrix_gsigma     = cb_p->GetParameter(1); 
@@ -1250,7 +1250,7 @@ int main ( int argc, char **argv)
       // gaussian fit to initialize
       TF1 *gausa = new TF1 ("gausa","[0]*exp(-1*(x-[1])*(x-[1])/2/[2]/[2])",h_peak-10*h_rms,h_peak+10*h_rms);
       gausa->SetParameters(h_norm,h_peak,h_rms);
-      H_eRatio[ii][myXib]->Fit("gausa","","",h_peak-3*h_rms,h_peak+3*h_rms);
+      H_eRatio[ii][myXib]->Fit(gausa,"","",h_peak-3*h_rms,h_peak+3*h_rms);
       double gausNorm  = gausa->GetParameter(0);
       double gausMean  = gausa->GetParameter(1);
       double gausSigma = fabs(gausa->GetParameter(2));
@@ -1269,7 +1269,7 @@ int main ( int argc, char **argv)
       cb_p->FixParameter(3, 5.);
       cb_p->SetParameter(4, gausNorm);
       cb_p->SetParLimits(2, 0.1, 5.);
-      H_eRatio[ii][myXib]->Fit("cb_p","lR","",myXmin,myXmax);
+      H_eRatio[ii][myXib]->Fit(cb_p,"lR","",myXmin,myXmax);
 
       double ratio_gmean      = cb_p->GetParameter(0);
       double ratio_gsigma     = cb_p->GetParameter(1); 
diff --git a/Validation/RecoEgamma/plugins/EgammaObjects.cc b/Validation/RecoEgamma/plugins/EgammaObjects.cc
index adecb9251a8cf..d53a70a61ad66 100644
--- a/Validation/RecoEgamma/plugins/EgammaObjects.cc
+++ b/Validation/RecoEgamma/plugins/EgammaObjects.cc
@@ -14,6 +14,8 @@
 
 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
 
+#include "TF1.h"
+
 EgammaObjects::EgammaObjects( const edm::ParameterSet& ps )
 {
   particleID = ps.getParameter<int>("particleID");
@@ -829,16 +831,18 @@ void EgammaObjects::getEfficiencyHistosViaDividing()
 
 void EgammaObjects::fitHistos()
 {
-  hist_EtOverTruth_->Fit("gaus","QEM");
+  //Use our own copy for thread safety
+  TF1 gaus("mygaus","gaus");
+  hist_EtOverTruth_->Fit(&gaus,"QEM");
 //  hist_EtNumRecoOverNumTrue_->Fit("pol1","QEM");
 
-  hist_EOverTruth_->Fit("gaus","QEM");
+  hist_EOverTruth_->Fit(&gaus,"QEM");
 //  hist_ENumRecoOverNumTrue_->Fit("pol1","QEM");
 
-  hist_EtaOverTruth_->Fit("gaus","QEM");
+  hist_EtaOverTruth_->Fit(&gaus,"QEM");
 //  hist_EtaNumRecoOverNumTrue_->Fit("pol1","QEM");
 
-  hist_PhiOverTruth_->Fit("gaus","QEM");
+  hist_PhiOverTruth_->Fit(&gaus,"QEM");
 //  hist_PhiNumRecoOverNumTrue_->Fit("pol1","QEM");
 
   /*
diff --git a/Validation/RecoParticleFlow/src/Comparator.cc b/Validation/RecoParticleFlow/src/Comparator.cc
index 0a5814078ac0a..8bbbdbd14ac0b 100644
--- a/Validation/RecoParticleFlow/src/Comparator.cc
+++ b/Validation/RecoParticleFlow/src/Comparator.cc
@@ -152,7 +152,7 @@ void Comparator::DrawGaussSigmaSlice(const char* key, const int rebinFactor, con
   fitfcndgssrms3->SetLineWidth(3);
   fitfcndgssrms3->SetLineStyle(2);
   fitfcndgssrms3->SetLineColor(4);
-  hrms->Fit("fitfcndgssrms3","0R");
+  hrms->Fit(fitfcndgssrms3,"0R");
 
   TH1D* ha=TH2Ana.SigmaGauss();
   TF1 *fitfcndgsse3 = new TF1("fitfcndgsse3",fitFunction_f,hrms->GetXaxis()->GetBinLowEdge(1),hrms->GetXaxis()->GetBinUpEdge(hrms->GetNbinsX()),4);
@@ -160,7 +160,7 @@ void Comparator::DrawGaussSigmaSlice(const char* key, const int rebinFactor, con
   fitfcndgsse3->SetLineWidth(3);
   fitfcndgsse3->SetLineStyle(1);
   fitfcndgsse3->SetLineColor(4);
-  ha->Fit("fitfcndgsse3","0R");
+  ha->Fit(fitfcndgsse3,"0R");
 
   dir = dir0_;
   dir->cd();
@@ -172,7 +172,7 @@ void Comparator::DrawGaussSigmaSlice(const char* key, const int rebinFactor, con
   fitfcndgssrmsb3->SetLineWidth(3);
   fitfcndgssrmsb3->SetLineStyle(2);
   fitfcndgssrmsb3->SetLineColor(2);
-  hrmsb->Fit("fitfcndgssrmsb3","0R");
+  hrmsb->Fit(fitfcndgssrmsb3,"0R");
 
   TH1D* hb=TH2Anab.SigmaGauss();
   TF1 *fitfcndgsseb3 = new TF1("fitfcndgsseb3",fitFunction_f,hrms->GetXaxis()->GetBinLowEdge(1),hrms->GetXaxis()->GetBinUpEdge(hrms->GetNbinsX()),4);
@@ -180,7 +180,7 @@ void Comparator::DrawGaussSigmaSlice(const char* key, const int rebinFactor, con
   fitfcndgsseb3->SetLineWidth(3);
   fitfcndgsseb3->SetLineStyle(1);
   fitfcndgsseb3->SetLineColor(2);
-  hb->Fit("fitfcndgsseb3","0R");
+  hb->Fit(fitfcndgsseb3,"0R");
 
   Draw(hb,ha,mode);
   //Draw(hrms,ha,mode);
@@ -216,7 +216,7 @@ void Comparator::DrawGaussSigmaOverMeanXSlice(const char* key, const int rebinFa
   fitXfcndgssrms3->SetLineWidth(3);
   fitXfcndgssrms3->SetLineStyle(2);
   fitXfcndgssrms3->SetLineColor(4);
-  hrms->Fit("fitXfcndgssrms3","0R");
+  hrms->Fit(fitXfcndgssrms3,"0R");
 
   TH1D* ha=TH2Ana.SigmaGauss();
   ha->Divide(meanXslice);
@@ -225,7 +225,7 @@ void Comparator::DrawGaussSigmaOverMeanXSlice(const char* key, const int rebinFa
   fitXfcndgsse3->SetLineWidth(3);
   fitXfcndgsse3->SetLineStyle(1);
   fitXfcndgsse3->SetLineColor(4);
-  ha->Fit("fitXfcndgsse3","0R");
+  ha->Fit(fitXfcndgsse3,"0R");
 
   dir = dir0_;
   dir->cd();
@@ -238,7 +238,7 @@ void Comparator::DrawGaussSigmaOverMeanXSlice(const char* key, const int rebinFa
   fitXfcndgssrmsb3->SetLineWidth(3);
   fitXfcndgssrmsb3->SetLineStyle(2);
   fitXfcndgssrmsb3->SetLineColor(2);
-  hrmsb->Fit("fitXfcndgssrmsb3","0R");
+  hrmsb->Fit(fitXfcndgssrmsb3,"0R");
 
   TH1D* hb=TH2Anab.SigmaGauss();
   hb->Divide(meanXslice);
@@ -247,7 +247,7 @@ void Comparator::DrawGaussSigmaOverMeanXSlice(const char* key, const int rebinFa
   fitXfcndgsseb3->SetLineWidth(3);
   fitXfcndgsseb3->SetLineStyle(1);
   fitXfcndgsseb3->SetLineColor(2);
-  hb->Fit("fitXfcndgsseb3","0R");
+  hb->Fit(fitXfcndgsseb3,"0R");
 
   Draw(hb,ha,mode);
   //Draw(hrms,ha,mode);
diff --git a/Validation/RecoParticleFlow/src/TH2Analyzer.cc b/Validation/RecoParticleFlow/src/TH2Analyzer.cc
index 78727c4820485..03a24bac3c637 100644
--- a/Validation/RecoParticleFlow/src/TH2Analyzer.cc
+++ b/Validation/RecoParticleFlow/src/TH2Analyzer.cc
@@ -266,7 +266,7 @@ void TH2Analyzer::ProcessSlice(const int i, TH1D* proj ) const {
       //proj->Draw();
       TF1* f1= new TF1("f1", "gaus", fitmin, fitmax);
       f1->SetParameters(proj->GetRMS(),proj->GetMean(),proj->GetBinContent(proj->GetMaximumBin()));
-      proj->Fit("f1","R", "", proj->GetXaxis()->GetXmin(), proj->GetXaxis()->GetXmax());
+      proj->Fit(f1,"R", "", proj->GetXaxis()->GetXmin(), proj->GetXaxis()->GetXmax());
   
       //std::ostringstream oss;
       //oss << i;