my_nr.h File Reference

collection of useful routines More...

#include "tbci/vector.h"
#include "tbci/matrix.h"
#include "tbci/constants.h"
#include "tbci/solver/gauss_jordan.h"

Include dependency graph for my_nr.h:

Go to the source code of this file.

Macros
#define	DELTA(i, j)   ((i==j)?1:0)
#define	DELTAX   1e-30
#define	CON2   2.0
#define	CON   SQRT2
#define	BIG   1.0e30
#define	NTAB   10
#define	SAFE   2.0

Functions
template<typename T >
NAMESPACE_TBCI void	do_nothing (T dummy)
template<typename T >
int	basis_trafo (long int value, long int basis, Vector< T > &index)
	INST (template< typename T > class Vector friend T partial_del(T(*func)(const Vector< T > &, const Vector< T > &), const Vector< T > &, const Vector< T > &, int, T, T &);) template< typename T > T partial_del(T(*func)(const Vector< T > &
Matrix< T >	a (10, 10)
Vector< T >	buf_plus (xlen)
Vector< T >	buf_minus (xlen)
STD__ cerr<< "\n error: partial_del: derivate variable index out of range \n";hh=h;buf_plus=x;buf_minus=x;buf_plus(mu)+=hh;buf_minus(mu)-=hh;a(0, 0)=((*func)(buf_plus, p)-(*func)(buf_minus, p))/(2.0 *hh);err=1.0e30;for(int i=1;i< 10;i++){hh/=SQRT2;buf_plus=x;buf_minus=x;buf_plus(mu)+=hh;buf_minus(mu)-=hh;a(0, i)=((*func)(buf_plus, p)-(*func)(buf_minus, p))/(2.0 *hh);fac=2.0;for(REGISTER int j=1;j<=i;j++){a(j, i)=(a(j-1, i)*fac-a(j-1, i-1))/(fac-1.0);fac *=2.0;errt=MAX(MATH__ fabs(a(j, i)-a(j-1, i)), MATH__ fabs(a(j, i)-a(j-1, i-1)));if(errt<=err){err=errt;ans=a(j, i);}}if(MATH__ fabs(a(i, i)-a(i-1, i-1)) >=MATH__ fabs(2.0 *err)) break;}return ans;}template< typename T > T	chisq (Vector< T > &data, Vector< T > &func, Vector< T > &sig)
	INST (template< typename T > class Vector friend T chisquare(T(*func)(const Vector< T > &, const Vector< T > &), Vector< T > &, Vector< T > &, Vector< T > &, Vector< T > &);) template< typename T > T chisquare(T(*func)(const Vector< T > &
T chi	ALIGN (MIN_ALIGN)=0
	for (REGISTER int i=0;i< anzpkt;i++)
	INST (template< typename T > class Vector friend TVector< T > grid_min(T(*func)(const Vector< T > &, const Vector< T > &), Vector< T > &, Vector< T > &, Vector< T > &, Vector< T > &, Vector< T > &, Vector< T > &, long int, int, char);) template< typename T > TVector< T > grid_min(T(*func)(const Vector< T > &
Vector< T >	p (pmin)
TVector< T >	min_pos (pmin)
Vector< T >	factor (plen)
Vector< T >	arg (1)
Vector< T >	fval (x.size())
Vector< T >	index (plen)
Vector< T >	b (plen, 0)
Vector< T >	a (plen, 0)
STD__ cerr<< "\n parameter space discretisation in progress ... \n ";factor(j)=(pmax(j)-pmin(j))/((T) res-(T) 1);for(long int i=0;i < len;i++){STD__ cerr<< "\n error: basis trafo failed\n";p(j)=factor(j)*index(j)+pmin(j);for(int l=0;l< max_it;l++){err=0.0;for(int k=0;k< plen;k++){for(int j=0;j < xlen;j++){arg(0)=x(j);b(k)+=(y(j)-(*func)(p, arg))/(sig(j)*sig(j))*partial_del(func, p, arg, k, x(j)*fac, err);}}err=0.0;for(int t=0;t< plen;t++){for(REGISTER int u=0;u< xlen;u++){arg(0)=x(u);a(t)+=sqr(partial_del(func, p, arg, t, x(u)*fac, err))/(sig(u)*sig(u));}}p(q)+=b(q)*lambda(q)/(a(q)+1e-30);}for(int q=0;q< xlen;q++){arg(0)=x(q);fval(q)=(*func)(p, arg);}chi=chisq(y, fval, sig);if(chi< min){min=chi;min_pos=p;}if(i%1000==0) STD__ cerr<< i<< "\r";}STD__ cerr<< " \n done! \n\n\n";return min_pos;}INST(template< typename T > class Vector friend T lev_mar(Vector< T > &, Vector< T > &, Vector< T > &, Vector< T > &, Vector< T > &, Vector< T > &, T(*func)(const Vector< T > &, const Vector< T > &), Vector< T >, Vector< T > &, T, int, char);) template< typename T > T	lev_mar (Vector< T > &x, Vector< T > &y, Vector< T > &sig, Vector< T > &param, Vector< T > &plow, Vector< T > &phigh, T(*func)(const Vector< T > &, const Vector< T > &), Vector< T > lambda, Vector< T > &lscale, T tol, int iter, char ver)
template<typename T >
void	hpsort (Matrix< T > &x, int pos)

Variables
const Vector< T > const Vector < T > &	x
const Vector< T > const Vector < T > const Vector< T > &	p
const Vector< T > const Vector < T > const Vector< T > int	mu
const Vector< T > const Vector < T > const Vector< T > int T	h
const Vector< T > const Vector < T > const Vector< T > int T T &	err
int	xlen = x.size()
const Vector< T > Vector< T > &	param
const Vector< T > Vector< T > Vector< T > Vector< T > &	y
const Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > &	sigma
T chi	wfvaldiff
int	anzpkt = x.size()
return	chi =0
const Vector< T > Vector< T > Vector< T > Vector< T > &	sig
const Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > &	pmin
const Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > &	pmax
const Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > &	lambda
const Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > long int	res
const Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > long int int	max_it
const Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > Vector< T > long int int char	v
long int	len =0
T min	fac =1e-8
	min =chisq(y, fval, sig)

Detailed Description

collection of useful routines

Definition in file my_nr.h.

Macro Definition Documentation

#define BIG   1.0e30

Definition at line 95 of file my_nr.h.

#define CON   SQRT2

Definition at line 94 of file my_nr.h.

#define CON2   2.0

Definition at line 93 of file my_nr.h.

#define DELTA	(		i,
			j
	)		((i==j)?1:0)

Definition at line 44 of file my_nr.h.

Referenced by lev_mar().

#define DELTAX   1e-30

Definition at line 69 of file my_nr.h.

#define NTAB   10

Definition at line 96 of file my_nr.h.

#define SAFE   2.0

Definition at line 97 of file my_nr.h.

Function Documentation

Matrix<T> a	(	10	,
		10
	)

Vector<T> a	(	plen	,
		0
	)

T chi ALIGN ( MIN_ALIGN  )

pure virtual

Vector<T> arg

(

1

)

Vector<T> b	(	plen	,
		0
	)

template<typename T >

int basis_trafo	(	long int	value,
		long int	basis,
		Vector< T > &	index
	)

Definition at line 73 of file my_nr.h.

References basis, i, index, log(), MATH__, TVector< T >::size(), std::sqrt(), and T.

Vector<T> buf_minus

(

xlen

)

Vector<T> buf_plus

(

xlen

)

STD__ cerr<< "\n error: partial_del: derivate variable index out of range \n"; hh = h ; buf_plus = x; buf_minus = x; buf_plus(mu) += hh; buf_minus(mu) -= hh; a(0,0) = ((*func)(buf_plus, p) - (*func)(buf_minus, p)) / (2.0*hh); err = 1.0e30 ; for (int i = 1; i < 10 ; i++) { hh /= SQRT2 ; buf_plus = x; buf_minus = x; buf_plus(mu) += hh; buf_minus(mu) -= hh; a(0,i) = ((*func)(buf_plus, p) - (*func)(buf_minus, p)) / (2.0*hh); fac = 2.0 ; for (REGISTER int j = 1; j <= i; j++) { a(j,i) = (a(j-1,i)*fac - a(j-1,i-1))/(fac-1.0); fac *= 2.0 ; errt = MAX (MATH__ fabs(a(j,i)-a(j-1,i)), MATH__ fabs(a(j,i)-a(j-1,i-1))); if (errt <= err) { err = errt; ans = a(j,i); } } if (MATH__ fabs (a(i,i)-a(i-1,i-1)) >= MATH__ fabs ( 2.0 *err)) break; } return ans;}template <typename T> T chisq	(	Vector< T > &	data,
		Vector< T > &	func,
		Vector< T > &	sig
	)

Definition at line 159 of file my_nr.h.

References ALIGN, chi, i, len, MIN_ALIGN, REGISTER, sig, TVector< T >::size(), sqr(), and T.

template<typename T >

NAMESPACE_TBCI void do_nothing ( T  dummy )

inline

Definition at line 41 of file my_nr.h.

Vector<T> factor

(

plen

)

for

(

)

Definition at line 183 of file my_nr.h.

References arg(), chi, i, param, sigma, sqr(), wfvaldiff, x, and y.

Vector<T> fval

(

x.

size()

)

template<typename T >

void hpsort	(	Matrix< T > &	x,
		int	pos
	)

Definition at line 453 of file my_nr.h.

References TMatrix< T >::columns(), i, TMatrix< T >::rows(), and x.

Vector<T> index

(

plen

)

INST

(

template< typename T > class Vector friend T partial_del(T(*func)(const Vector< T > &, const Vector< T > &), const Vector< T > &, const Vector< T > &, int, T, T &);

)

const

INST

(

template< typename T > class Vector friend T chisquare(T(*func)(const Vector< T > &, const Vector< T > &), Vector< T > &, Vector< T > &, Vector< T > &, Vector< T > &);

)

const

INST

(

template< typename T > class Vector friend TVector< T > grid_min(T(*func)(const Vector< T > &, const Vector< T > &), Vector< T > &, Vector< T > &, Vector< T > &, Vector< T > &, Vector< T > &, Vector< T > &, long int, int, char);

)

const

STD__ cerr<< "\n parameter space discretisation in progress ... \n "; factor(j) = (pmax(j) - pmin(j)) / ((T)res - (T)1); for (long int i=0;i<len;i++) { STD__ cerr << "\n error: basis trafo failed\n"; p(j) = factor(j) * index(j) + pmin(j); for (int l=0;l<max_it;l++) { err=0.0; for (int k=0;k<plen;k++) { for (int j=0;j<xlen;j++) { arg(0)=x(j); b(k)+=(y(j)-(*func)(p, arg))/(sig(j)*sig(j))* partial_del(func, p, arg, k, x(j)*fac, err); } } err = 0.0; for (int t=0; t<plen; t++) { for (REGISTER int u=0; u<xlen; u++) { arg(0) = x(u); a(t) += sqr(partial_del(func, p, arg, t, x(u)*fac, err))/ (sig(u)*sig(u)); } } p(q) += b(q)*lambda(q)/(a(q)+ 1e-30 ); } for (int q=0; q<xlen; q++) { arg(0)=x(q); fval(q)=(*func)(p, arg); } chi=chisq(y, fval, sig); if ( chi<min ) { min=chi; min_pos=p; } if ( i%1000==0 ) STD__ cerr << i << "\r"; } STD__ cerr << " \n done! \n\n\n"; return min_pos;}INST(template <typename T> class Vector friend T lev_mar(Vector<T>&, Vector<T>&, Vector<T>&, Vector<T>&, Vector<T>&, Vector<T>&, T (*func)(const Vector<T>&, const Vector<T>&), Vector<T>, Vector<T>&, T, int, char);)template <typename T> T lev_mar	(	Vector< T > &	x,
		Vector< T > &	y,
		Vector< T > &	sig,
		Vector< T > &	param,
		Vector< T > &	plow,
		Vector< T > &	phigh,
		T(*)(const Vector< T > &, const Vector< T > &)	func,
		Vector< T >	lambda,
		Vector< T > &	lscale,
		T	tol,
		int	iter,
		char	ver
	)

Definition at line 355 of file my_nr.h.

References a, ALIGN, arg(), b, chi, DELTA, err, std::fabs(), fac, fval(), gaussj(), lambda, len, MATH__, MIN_ALIGN, param, REGISTER, sig, TVector< T >::size(), STD__, T, x, and y.

TVector<T> min_pos

(

pmin

)

Vector<T> p

(

pmin

)

Variable Documentation

int anzpkt = x.size()

Definition at line 181 of file my_nr.h.

T min chi =0

Definition at line 189 of file my_nr.h.

T min err

Initial value:

{
   
   T errt ALIGN(MIN_ALIGN), fac, hh, ans(0)

Definition at line 106 of file my_nr.h.

T min fac =1e-8

Definition at line 296 of file my_nr.h.

const Vector<T> const Vector<T> const Vector<T> int T h

Definition at line 101 of file my_nr.h.

const Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> Vector<T>& lambda

Definition at line 277 of file my_nr.h.

len =0

Definition at line 287 of file my_nr.h.

const Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> long int int max_it

Definition at line 277 of file my_nr.h.

min =chisq(y, fval, sig)

Definition at line 300 of file my_nr.h.

const Vector<T> const Vector<T> const Vector<T> int mu

Definition at line 101 of file my_nr.h.

const Vector<T> const Vector<T> const Vector<T>& p

Definition at line 101 of file my_nr.h.

const Vector<T> Vector<T>& param

Definition at line 172 of file my_nr.h.

const Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> Vector<T>& pmax

Definition at line 277 of file my_nr.h.

const Vector<T> Vector<T> Vector<T> Vector<T> Vector<T>& pmin

Definition at line 277 of file my_nr.h.

const Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> long int res

Definition at line 277 of file my_nr.h.

const Vector<T> Vector<T> Vector<T> Vector<T>& sig

Definition at line 277 of file my_nr.h.

const Vector<T> Vector<T> Vector<T> Vector<T> Vector<T>& sigma

Initial value:

{
  Vector<T> arg(1)

Definition at line 177 of file my_nr.h.

const Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> Vector<T> long int int char v

Initial value:

{

   int plen=pmin.size()

Definition at line 283 of file my_nr.h.

T chi wfvaldiff

Definition at line 180 of file my_nr.h.

const Vector<T> Vector<T>& x

Definition at line 101 of file my_nr.h.

int xlen = x.size()

Definition at line 109 of file my_nr.h.

const Vector<T> Vector<T> Vector<T>& y

Definition at line 172 of file my_nr.h.