typedef ROOT::Math::IParamMultiFunction IModelFunction;

class KFitResult : public TFitResult {
	public:
	using TFitResult::TFitResult;
	void ResetModelFunction(TF1* func){
		this->SetModelFunction(std::shared_ptr<IModelFunction>(dynamic_cast<IModelFunction*>(ROOT::Math::WrappedMultiTF1(*func).Clone())));
	}
};

void getresult(TF1* gsn, KFitResult& res){
	//print the fit result and error for a point
	double x[1] = {1.0};
	double err[1] = {0.0};
	res.GetConfidenceIntervals(1,1,1,x,err,0.683,false);
	cout << "f(x=" << x[0] << ") = " << gsn->Eval(x[0]) << " +/- " << err[0] << endl;
}

void dofit(){
	TF1* gsn = new TF1("gsn","gaus",-5,5);
	gsn->SetParameters(1,0,3);
	TH1F* hist = new TH1F("hist","",20,-5,5);
	hist->FillRandom("gsn",1000);
	TFitResultPtr res = hist->Fit("gsn","S0Q","",-5,5);
	KFitResult kres(*res.Get());
	getresult(gsn,kres);
	//save fit result
	TFile* fout = TFile::Open("gsnfit2.root","RECREATE");
	gsn->Write();
	res->SetName("err_gsn");
	res->Write();
	fout->Close();
}

void readfit(){
	TFile* fin = TFile::Open("gsnfit2.root");
	TF1* gsn = (TF1*)fin->Get("gsn");
	TFitResult* res = (TFitResult*)fin->Get("err_gsn");
	KFitResult kres(*res);
	//restore base function
	kres.ResetModelFunction(gsn);
	//check result
	getresult(gsn,kres);
}

void testfitresult2(){
	dofit();
	readfit();
}