oblique_mc

Monte Carlo in Single Layer Biological Tissue
git clone anongit@rnpnr.xyz:oblique_mc.git
Log | Files | Refs | Feed | README | LICENSE

Commit: 804980d1a2f6f18a0b01e710d393f2501bbcfb23
Parent: e98db29ddb5c4dc67104c0e4c5adfb3aff342110
Author: Randy Palamar
Date:   Tue, 26 Mar 2024 11:17:54 -0600

add build.sh and fix build

Diffstat:
MPhoton.cpp | 2+-
DUtilities.cpp | 94-------------------------------------------------------------------------------
DUtilities.h | 140-------------------------------------------------------------------------------
Abuild.sh | 5+++++
Mmain.cpp | 4++--
Autils.cpp | 93+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Autils.h | 141+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
7 files changed, 242 insertions(+), 237 deletions(-)

diff --git a/Photon.cpp b/Photon.cpp @@ -3,7 +3,7 @@ ////////////////////////////////////////////////////////////////////// #include "Photon.h" -#include "Utilities.h" +#include "utils.h" #include <math.h> ///////////////////////////////////////// diff --git a/Utilities.cpp b/Utilities.cpp @@ -1,94 +0,0 @@ - -#include "Utilities.h" - -/*********************************************************** -* A random number generator from Numerical Recipes in C. ****/ -#define MBIG 1000000000 -#define MSEED 161803398 -#define MZ 0 -#define FAC 1.0E-9 -float ran3(int &idum) -{ - static int inext,inextp; - static long ma[56]; - static int iff=0; - long mj,mk; - int i,ii,k; - if (idum < 0 || iff == 0) - { - iff=1; - mj=MSEED-(idum < 0 ? -idum : idum); - mj %= MBIG; - ma[55]=mj; - mk=1; - for (i=1;i<=54;i++) - { - ii=(21*i) % 55; - ma[ii]=mk; - mk=mj-mk; - if (mk < MZ) mk += MBIG; - mj=ma[ii]; - } - for (k=1;k<=4;k++) - for (i=1;i<=55;i++) - { - ma[i] -= ma[1+(i+30) % 55]; - if (ma[i] < MZ) ma[i] += MBIG; - } - inext=0; - inextp=31; - idum=1; - } - if (++inext == 56) - inext=1; - if (++inextp == 56) - inextp=1; - mj=ma[inext]-ma[inextp]; - if (mj < MZ) mj += MBIG; - ma[inext]=mj; - return mj*FAC; -} -#undef MBIG -#undef MSEED -#undef MZ -#undef FAC - - -/*********************************************************** - * Generate a random number between 0 and 1. Take a - * number as seed the first time entering the function. - * The seed is limited to 1<<15. - * Wang et al found that when idum is too large, ran3 may return - * numbers beyond 0 and 1. - ****/ -//#define STANDARDTEST 0 /* testing program using fixed rnd seed. */ -double RandomNum(void) -{ - static bool first_time=1; - static int idum; /* seed for ran3. */ - if(first_time) - { - #if STANDARDTEST /* Use fixed seed to test the program. */ - idum = - 1; - #else - idum = -(int)time(NULL)%(1<<15); /* use 16-bit integer as the seed. */ - #endif - ran3(idum); - first_time = 0; - idum = 1; - } - return( (double)ran3(idum) ); -} - -///////////////////////////////////////////// -int min(int d1,int d2) -{ - if (d1<d2) - { - return d1; - } - else - { - return d2; - } -} diff --git a/Utilities.h b/Utilities.h @@ -1,140 +0,0 @@ -/* Providing global variables and functions*/ -#include <math.h> -#include <new> -#include <time.h> - -//global constants -#ifndef PI - #define PI 3.14159265 -#endif - -#ifndef ZERO - #define ZERO 1.0e-6 //for small deflection angle -#endif - -#ifndef n0 - #define n0 1.0 //refractive index of ambient medium -#endif - -#ifndef BIGVAL - #define BIGVAL 1e6//big value -#endif - -//global variables -#ifdef GLOBALORIGIN - #define GLOBAL -#else - #define GLOBAL extern -#endif -GLOBAL double dx; -GLOBAL double dy; -GLOBAL int Nx; -GLOBAL int Ny; -GLOBAL double zc; -GLOBAL double n; -GLOBAL double mua; -GLOBAL double mus; -GLOBAL double mut; -GLOBAL double g; -GLOBAL double d; - - -GLOBAL double th_in; // incident angle in degree - -GLOBAL double x_offset; -GLOBAL double y_offset; - -//generate random data -float ran3(int *idum); -double RandomNum(void); - - -//other auxilliary functions -int min(int d1, int d2); -#define SIGN(x) ((x)>=0?1:-1) - - -////////allocte dynamic array///////////////////////// -using namespace std; -template <class T> -T* alloc1D(int dim1, T t) -{ - T* p; - try - { - p=new T[dim1]; - } - catch(bad_alloc) - { - cerr<<"Fail to allocate 1D array!"<<endl; - exit(1); - } - for (int tmp1=0;tmp1<dim1;tmp1++) - { - p[tmp1]=t; - } - return p; -} -template <class T> -void dealloc1D(T* p) -{ - delete [] p; -} -/////////// -template <class T> -T** alloc2D(int dim1,int dim2,T t) -{ - T** p; - try - { - p=new T*[dim1]; - for (int tmp1=0;tmp1<dim1;tmp1++) - { - p[tmp1]=alloc1D(dim2,t); - } - } - catch (bad_alloc) - { - cerr<<"Fail to allocate 2D array!"; - } - return p; -} -template <class T> -void dealloc2D(int dim1,T** p) -{ - for (int tmp1=0;tmp1<dim1;tmp1++) - { - dealloc1D(p[tmp1]); - } - delete [] p; -} -///////// -template <class T> -T*** alloc3D(int dim1,int dim2,int dim3,T t) -{ - T*** p; - try - { - p=new T**[dim1]; - for (int tmp1=0;tmp1<dim1;tmp1++) - { - p[tmp1]=alloc2D(dim2,dim3,t); - } - } - catch (bad_alloc) - { - cerr<<"Fail to allocate 3D array!"; - } - return p; -} -template <class T> -void dealloc3D(int dim1,int dim2,T*** p) -{ - for (int tmp1=0;tmp1<dim1;tmp1++) - { - dealloc2D(dim2,p[tmp1]); - } - delete [] p; -} -//////////////////////////////////////////////// - diff --git a/build.sh b/build.sh @@ -0,0 +1,5 @@ +#!/bin/sh + +srcs="main.cpp Photon.cpp utils.cpp" + +c++ -O3 -Wall -march=native $srcs -o mcml diff --git a/main.cpp b/main.cpp @@ -13,14 +13,14 @@ refrence: Wang, "Biomedical optics", chpater3&6 */ #include <string> #define GLOBALORIGIN -#include "Utilities.h" +#include "utils.h" #undef GLOABLORIGIN #include "Photon.h" //////////////////////////////////////////////// -void main() +int main(void) { ////////Input Interface//////////////////// cout<<"*****************************************"<<endl; diff --git a/utils.cpp b/utils.cpp @@ -0,0 +1,93 @@ +#include "utils.h" + +/*********************************************************** +* A random number generator from Numerical Recipes in C. ****/ +#define MBIG 1000000000 +#define MSEED 161803398 +#define MZ 0 +#define FAC 1.0E-9 +float ran3(int &idum) +{ + static int inext,inextp; + static long ma[56]; + static int iff=0; + long mj,mk; + int i,ii,k; + if (idum < 0 || iff == 0) + { + iff=1; + mj=MSEED-(idum < 0 ? -idum : idum); + mj %= MBIG; + ma[55]=mj; + mk=1; + for (i=1;i<=54;i++) + { + ii=(21*i) % 55; + ma[ii]=mk; + mk=mj-mk; + if (mk < MZ) mk += MBIG; + mj=ma[ii]; + } + for (k=1;k<=4;k++) + for (i=1;i<=55;i++) + { + ma[i] -= ma[1+(i+30) % 55]; + if (ma[i] < MZ) ma[i] += MBIG; + } + inext=0; + inextp=31; + idum=1; + } + if (++inext == 56) + inext=1; + if (++inextp == 56) + inextp=1; + mj=ma[inext]-ma[inextp]; + if (mj < MZ) mj += MBIG; + ma[inext]=mj; + return mj*FAC; +} +#undef MBIG +#undef MSEED +#undef MZ +#undef FAC + + +/*********************************************************** + * Generate a random number between 0 and 1. Take a + * number as seed the first time entering the function. + * The seed is limited to 1<<15. + * Wang et al found that when idum is too large, ran3 may return + * numbers beyond 0 and 1. + ****/ +//#define STANDARDTEST 0 /* testing program using fixed rnd seed. */ +double RandomNum(void) +{ + static bool first_time=1; + static int idum; /* seed for ran3. */ + if(first_time) + { + #if STANDARDTEST /* Use fixed seed to test the program. */ + idum = - 1; + #else + idum = -(int)time(NULL)%(1<<15); /* use 16-bit integer as the seed. */ + #endif + ran3(idum); + first_time = 0; + idum = 1; + } + return( (double)ran3(idum) ); +} + +///////////////////////////////////////////// +int min(int d1,int d2) +{ + if (d1<d2) + { + return d1; + } + else + { + return d2; + } +} diff --git a/utils.h b/utils.h @@ -0,0 +1,141 @@ +/* Providing global variables and functions*/ +#include <math.h> +#include <new> +#include <iostream> +#include <time.h> + +//global constants +#ifndef PI + #define PI 3.14159265 +#endif + +#ifndef ZERO + #define ZERO 1.0e-6 //for small deflection angle +#endif + +#ifndef n0 + #define n0 1.0 //refractive index of ambient medium +#endif + +#ifndef BIGVAL + #define BIGVAL 1e6//big value +#endif + +//global variables +#ifdef GLOBALORIGIN + #define GLOBAL +#else + #define GLOBAL extern +#endif +GLOBAL double dx; +GLOBAL double dy; +GLOBAL int Nx; +GLOBAL int Ny; +GLOBAL double zc; +GLOBAL double n; +GLOBAL double mua; +GLOBAL double mus; +GLOBAL double mut; +GLOBAL double g; +GLOBAL double d; + + +GLOBAL double th_in; // incident angle in degree + +GLOBAL double x_offset; +GLOBAL double y_offset; + +//generate random data +float ran3(int *idum); +double RandomNum(void); + + +//other auxilliary functions +int min(int d1, int d2); +#define SIGN(x) ((x)>=0?1:-1) + + +////////allocte dynamic array///////////////////////// +using namespace std; +template <class T> +T* alloc1D(int dim1, T t) +{ + T* p; + try + { + p=new T[dim1]; + } + catch(bad_alloc) + { + cerr<<"Fail to allocate 1D array!"<<endl; + exit(1); + } + for (int tmp1=0;tmp1<dim1;tmp1++) + { + p[tmp1]=t; + } + return p; +} +template <class T> +void dealloc1D(T* p) +{ + delete [] p; +} +/////////// +template <class T> +T** alloc2D(int dim1,int dim2,T t) +{ + T** p; + try + { + p=new T*[dim1]; + for (int tmp1=0;tmp1<dim1;tmp1++) + { + p[tmp1]=alloc1D(dim2,t); + } + } + catch (bad_alloc) + { + cerr<<"Fail to allocate 2D array!"; + } + return p; +} +template <class T> +void dealloc2D(int dim1,T** p) +{ + for (int tmp1=0;tmp1<dim1;tmp1++) + { + dealloc1D(p[tmp1]); + } + delete [] p; +} +///////// +template <class T> +T*** alloc3D(int dim1,int dim2,int dim3,T t) +{ + T*** p; + try + { + p=new T**[dim1]; + for (int tmp1=0;tmp1<dim1;tmp1++) + { + p[tmp1]=alloc2D(dim2,dim3,t); + } + } + catch (bad_alloc) + { + cerr<<"Fail to allocate 3D array!"; + } + return p; +} +template <class T> +void dealloc3D(int dim1,int dim2,T*** p) +{ + for (int tmp1=0;tmp1<dim1;tmp1++) + { + dealloc2D(dim2,p[tmp1]); + } + delete [] p; +} +//////////////////////////////////////////////// +