Remove calc_distance (last used by the deleted 'calc_lx')

5 files changed, 1 insertions, 252 deletions

diff --git a/mass_profile/Makefile b/mass_profile/Makefile
index 4f037e6..1d7f48b 100644
--- a/mass_profile/Makefile
+++ b/mass_profile/Makefile
@@ -24,7 +24,7 @@ endif
 OPT_UTIL_INC ?= -I../opt_utilities
 
 TARGETS= fit_nfw_sbp fit_dbeta_sbp fit_beta_sbp \
-		fit_direct_beta calc_distance fit_wang2012_model \
+		fit_direct_beta fit_wang2012_model \
 		fit_nfw_mass fit_mt_pl fit_lt_pl fit_mt_bpl \
 		fit_lt_bpl cooling_time calc_lx_dbeta calc_lx_beta
 HEADERS= projector.hpp spline.hpp vchisq.hpp
@@ -57,9 +57,6 @@ calc_lx_dbeta: calc_lx_dbeta.o beta_cfg.o report_error.o
 calc_lx_beta: calc_lx_beta.o beta_cfg.o dump_fit_qdp.o
 	$(CXX) $(CXXFLAGS) $^ -o $@ $(OPT_UTIL_INC)
 
-calc_distance: calc_distance_main.o calc_distance.o
-	$(CXX) $(CXXFLAGS) -o $@ $^
-
 fit_mt_pl: fit_mt_pl.cpp
 	$(CXX) $(CXXFLAGS) $< -o $@ $(OPT_UTIL_INC)
 
@@ -106,9 +103,6 @@ report_error.o: report_error.cpp report_error.hpp
 dump_fit_qdp.o: dump_fit_qdp.cpp dump_fit_qdp.hpp
 	$(CXX) $(CXXFLAGS) -c $< $(OPT_UTIL_INC)
 
-calc_distance.o: calc_distance.cc calc_distance.h adapt_trapezoid.h
-	$(CXX) $(CXXFLAGS) -c $<
-
 %.o: %.cc
 	$(CXX) $(CXXFLAGS) -c $< $(OPT_UTIL_INC)
 
diff --git a/mass_profile/adapt_trapezoid.h b/mass_profile/adapt_trapezoid.h
deleted file mode 100644
index 35d22b4..0000000
--- a/mass_profile/adapt_trapezoid.h
+++ /dev/null
@@ -1,124 +0,0 @@
-#ifndef ADAPT_TRAPEZOID_H
-#define ADAPT_TRAPEZOID_H
-#include <list>
-#include <utility>
-#include <cassert>
-//#include <iostream>
-
-//using namespace std;
-
-template<typename T1,typename T2,typename T3>
-class triple
-{
- public:
-  T1 first;
-  T2 second;
-  T3 third;
-  triple(T1 x1,T2 x2,T3 x3)
-    :first(x1),second(x2),third(x3)
-    {
-    }
-};
-
-
-template<typename T1,typename T2,typename T3>
-triple<T1,T2,T3> make_triple(T1 x1,T2 x2,T3 x3)
-{
-  return triple<T1,T2,T3>(x1,x2,x3);
-}
-
-template <typename T>
-T trapezoid(T (*fun)(T),T x1,T x2,T err_limit)
-{
-  int n=256;
-  T result;
-  const int max_division=24;
-  T old_value=0;
-  for(int i=1;i<max_division;i++)
-    {
-      result=0.;
-      n*=2;
-      T step=(x2-x1)/n;
-      for(int j=0;j<n;j++)
-	{
-	  result+=(fun(x1+(j+1)*step)+fun(x1+j*step))*step/T(2.);
-	}
-      old_value-=result;
-      old_value=old_value<0?-old_value:old_value;
-      if(old_value<err_limit)
-	{
-	  return result;
-	}
-      old_value=result;
-    }
-}
-
-
-template <typename T>
-T adapt_trapezoid(T (*fun)(T),T x1,T x2,T err_limit)
-{
-  //  const err_limit=.001;
-  typedef triple<T,T,bool> interval;
-  /*T for interval type,
-    bool for state trur for still to be updated,
-    false for do not need to be updated
-   */
-  std::list<interval> interval_list;
-  T current_sum=((fun(x1)+fun(x2))/2.*(x2-x1));
-  interval_list.push_back(make_triple(x1,current_sum,true));
-  interval_list.push_back(make_triple(x2,(T)0.,true));
-  bool int_state=1;
-  int n_intervals=1;
-  while(int_state)
-    {
-      //std::cout<<n_intervals<<std::endl;
-      int_state=0;
-      typename std::list<interval>::iterator i1=interval_list.begin();
-      typename std::list<interval>::iterator i2=interval_list.begin();
-      i2++;
-      for(;i2!=interval_list.end();i1++,i2=i1,i2++)
-	{
-	  //cout<<i1->first<<"\t"<<i2->first<<endl;
-	  //assert(i2->first>i1->first);
-	  if(i1->third)
-	    {
-	      interval new_interval((i1->first+i2->first)/2,0,true);
-	      
-	      T sum1,sum2;
-	      sum1=(fun(new_interval.first)+fun(i1->first))/2*(new_interval.first-i1->first);
-	      sum2=(fun(new_interval.first)+fun(i2->first))/2*(i2->first-new_interval.first);
-	      new_interval.second=sum2;
-	      T err;
-	      err=i1->second-sum1-sum2;
-	      err=err<0?-err:err;
-
-	      if(err>err_limit/n_intervals)
-		{
-		  i1->second=sum1;
-		  interval_list.insert(i2,new_interval);
-		  n_intervals++;
-		  if(n_intervals>10e6)
-		    {
-		      
-		      break;
-		    }
-		}
-	      else
-		{
-		  i1->third=false;
-		}
-	      int_state=1;
-	    }
-	}
-      
-    }
-  T result=0;
-  for(typename std::list<interval>::iterator i=interval_list.begin();i!=interval_list.end();i++)
-    {
-      result+=i->second;
-    }
-  return result;
-}
-
-#endif
-//end of the file
diff --git a/mass_profile/calc_distance.cc b/mass_profile/calc_distance.cc
deleted file mode 100644
index 91f2ea3..0000000
--- a/mass_profile/calc_distance.cc
+++ /dev/null
@@ -1,61 +0,0 @@
-#include <iostream>
-#include <cmath>
-#include <cstdlib>
-#include <cstddef>
-#include <cassert>
-#include "adapt_trapezoid.h"
-
-//calc_distance
-//usage:
-//calc_distance z
-
-using namespace std;
-
-static double cm=1;
-static double s=1;
-static double km=1000*100;
-static double Mpc=3.08568e+24*cm;
-//static double kpc=3.08568e+21*cm;
-//static double yr=365.*24.*3600.;
-//static double Gyr=1e9*yr;
-static double H=71.*km/s/Mpc;
-static const double c=299792458.*100.*cm;
-//const double c=3e8*100*cm;
-static const double omega_m=0.27;
-static const double omega_l=0.73;
-static const double arcsec2arc_ratio=1./60/60/180*3.1415926;
-
-
-double E(double z)
-{
-  double omega_k=1-omega_m-omega_l;
-  return sqrt(omega_m*(1+z)*(1+z)*(1+z)+omega_k*(1+z)*(1+z)+omega_l);
-}
-
-double f_dist(double z)
-{
-  return 1/E(z);
-}
-
-double f_age(double z)
-{
-  return f_dist(1/z)/(z*z);
-}
-
-
-
-double calc_angular_distance(double z)
-{
-  //return c/H*integer(f_dist,0,z)/(1+z);
-  //return c/H*adapt_trapezoid(f_dist,0.,z,1e-4)/(1+z);
-  return adapt_trapezoid(f_dist,0.,z,1e-4)/(1+z);
-}
-
-double calc_luminosity_distance(double z)
-{
-  //return c/H*integer(f_dist,0,z)/(1+z);
-  return c/H*adapt_trapezoid(f_dist,0.,z,1e-4)*(1+z);
-}
-
-
-//EOF
diff --git a/mass_profile/calc_distance.h b/mass_profile/calc_distance.h
deleted file mode 100644
index ba01316..0000000
--- a/mass_profile/calc_distance.h
+++ /dev/null
@@ -1,8 +0,0 @@
-#ifndef CALC_DISTANCE_H
-#define CALC_DISTANCE_H
-extern double calc_angular_distance(double z);
-extern double calc_luminosity_distance(double z);
-extern double E(double z);
-#endif
-
-//EOF
diff --git a/mass_profile/calc_distance_main.cc b/mass_profile/calc_distance_main.cc
deleted file mode 100644
index 28bdb70..0000000
--- a/mass_profile/calc_distance_main.cc
+++ /dev/null
@@ -1,52 +0,0 @@
-#include "calc_distance.h"
-#include <cstdlib>
-#include <cmath>
-#include <iostream>
-using namespace std;
-
-static double cm=1;
-static double s=1;
-static double km=1000*100;
-static double Mpc=3.08568e+24*cm;
-static double kpc=3.08568e+21*cm;
-//static double yr=365.*24.*3600.;
-//static double Gyr=1e9*yr;
-static double H=71.*km/s/Mpc;
-static const double c=299792458.*100.*cm;
-//const double c=3e8*100*cm;
-static const double pi=4*atan(1);
-static const double omega_m=0.27;
-static const double omega_l=0.73;
-static const double arcsec2arc_ratio=1./60/60/180*pi;
-
-
-
-
-int main(int argc,char* argv[])
-{
-  if(argc<2)
-    {
-      cerr<<"Usage:"<<argv[0]<<" z"<<endl;
-      exit(-1);
-    }
-  if(argc==3)
-    {
-      H=atof(argv[2])*km/s/Mpc;
-    }
-  double z=atof(argv[1]);
-  double d=c/H*calc_angular_distance(z);
-  //double age=calc_age(z);
-  //cout<<d<<endl;
-  cout<<"d_a_cm= "<<d<<" #angular distance in cm"<<endl;
-  cout<<"d_a_mpc= "<<d/Mpc<<"#angular distance in Mpc"<<endl;
-  cout<<"d_l_cm= "<<(1+z)*(1+z)*d<<" #luminosity distance in cm"<<endl;
-  cout<<"d_l_mpc= "<<(1+z)*(1+z)*d/Mpc<<" #luminosity in Mpc"<<endl;
-  cout<<"kpc_per_sec= "<<d/kpc*arcsec2arc_ratio<<" #kpc per arcsec"<<endl;
-  cout<<"For Chandra:"<<endl;
-  cout<<"kpc_per_pixel= "<<d/kpc*0.492*arcsec2arc_ratio<<" #kpc per chandra pixel"<<endl;
-  cout<<"cm_per_pixel= "<<d*.492*arcsec2arc_ratio<<" #cm per chandra pixel"<<endl;
-  cout<<"norm= "<<1E-14/(4*pi*pow(d*(1+z),2))<<" #norm used to calc cooling function"<<endl;
-  //cout<<ddivid(calc_distance,d,0,1,.0001)<<endl;
-  cout<<"E(z)= "<<E(z)<<endl;
-}
-