From 1f716fe1c51904cba8464a4c999110554e98e363 Mon Sep 17 00:00:00 2001
From: astrojhgu <astrojhgu@ed2142bd-67ad-457f-ba7c-d818d4011675>
Date: Tue, 4 Jan 2011 15:20:38 +0000
Subject: git-svn-id: file:///home/svn/opt_utilities@154
 ed2142bd-67ad-457f-ba7c-d818d4011675

---
 error_estimator/error_estimator.hpp | 17 +++++++++++++----
 example/test_fitter.cpp             |  7 +++++++
 2 files changed, 20 insertions(+), 4 deletions(-)

diff --git a/error_estimator/error_estimator.hpp b/error_estimator/error_estimator.hpp
index 4324cd0..a6c6163 100644
--- a/error_estimator/error_estimator.hpp
+++ b/error_estimator/error_estimator.hpp
@@ -1,6 +1,7 @@
 #ifndef ERROR_EST
 #define ERROR_EST
 #include <core/fitter.hpp>
+#include <core/freeze_param.hpp>
 #include <iostream>
 #include <vector>
 #include <string>
@@ -8,7 +9,7 @@
 namespace opt_utilities
 {
   template <typename Ty,typename Tx,typename Tp,typename Ts,typename Tstr>
-  void estimate_error(fitter<Ty,Tx,Tp,Ts>& fit,const Tstr& pname,typename element_type_trait<Tp>::element_type& lower,typename element_type_trait<Tp>::element_type& upper,const Ts& dchi)
+  void estimate_error(fitter<Ty,Tx,Tp,Ts>& fit,const Tstr& pname,typename element_type_trait<Tp>::element_type& lower,typename element_type_trait<Tp>::element_type& upper,const Ts& dchi,const Ts& precision)
   {
     typedef typename element_type_trait<Tp>::element_type Tpe;
     std::vector<Tstr> pnames;
@@ -41,16 +42,24 @@ namespace opt_utilities
     
     Tpe current_value=
       fit.get_param_value(pname);
+    
+    if(lower>=current_value||upper<=current_value)
+      {
+	std::cerr<<"Error, initial lower and upper limits should be smaller and larger than the best fit value, respectively"<<std::endl;
+	return;
+      }
+
     Tp current_params(fit.get_model().get_all_params());
     Ts current_chi=fit.get_statistic_value();
     Tpe upper1=current_value;
     Tpe step=upper-current_value;
+    
     while(1)
       {
 	fit.set_param_value(pname,upper);
 	fit.fit();
 	Ts chi=fit.get_statistic_value();
-	cerr<<upper<<"\t"<<chi<<std::endl;
+	std::cerr<<upper<<"\t"<<chi<<std::endl;
 	if(chi>current_chi+dchi)
 	  {
 	    break;
@@ -59,7 +68,7 @@ namespace opt_utilities
       }
     
     Ts target_chi=current_chi+dchi;
-    while(upper-upper1>1E-10)
+    while(upper-upper1>std::abs(precision))
       {
 	std::cerr<<upper<<"\t"<<upper1<<std::endl;
 	fit.set_param_value(pname,upper1);
@@ -114,7 +123,7 @@ namespace opt_utilities
     fit.set_param_value(pname,current_value);
     fit.fit();
     
-    while(lower1-lower>1E-10)
+    while(lower1-lower>std::abs(precision))
       {
 	fit.set_param_value(pname,lower1);
 	fit.fit();
diff --git a/example/test_fitter.cpp b/example/test_fitter.cpp
index 8da9d86..a658a6c 100644
--- a/example/test_fitter.cpp
+++ b/example/test_fitter.cpp
@@ -1,6 +1,7 @@
 #include <core/optimizer.hpp>
 #include <methods/powell/powell_method.hpp>
 #include <methods/aga/aga.hpp>
+#include <error_estimator/error_estimator.hpp>
 #include <core/fitter.hpp>
 #include <vector>
 #include <iostream>
@@ -56,4 +57,10 @@ int main()
   f.fit();
   cout<<f.get_param_value("k")<<endl;
   cout<<f.get_param_value("b")<<endl;
+  
+  double lower=3;
+  double upper=6;
+  
+  estimate_error(f,std::string("b"),lower,upper,1.,1E-10);
+  std::cout<<lower<<"\t"<<upper<<endl;
 }
-- 
cgit v1.2.2