From b6179de6d7efffaef298858c1a56f9079f21a89a Mon Sep 17 00:00:00 2001 From: dmayerich Date: Fri, 4 Oct 2013 15:40:10 -0500 Subject: [PATCH] added scripts for spectral simulation --- defaults.h | 4 ++-- fileout.cu | 21 ++++++++++++++++++--- fileout.h | 9 ++++++--- options.h | 21 ++++++++++++++++++--- specimage.py | 149 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++------------------------------------------------------------- 5 files changed, 132 insertions(+), 72 deletions(-) diff --git a/defaults.h b/defaults.h index dc5182f..5ed9204 100644 --- a/defaults.h +++ b/defaults.h @@ -75,10 +75,10 @@ #define DEFAULT_INTENSITY_FILE "testappend" #define DEFAULT_TRANSMITTANCE_FILE "" -#define DEFAULT_ABSORBANCE_FILE "" +#define DEFAULT_ABSORBANCE_FILE "out_a" #define DEFAULT_NEAR_FILE "out_n.bmp" #define DEFAULT_FAR_FILE "out_f.bmp" -#define DEFAULT_EXTENDED_SOURCE "" +#define DEFAULT_EXTENDED_SOURCE "einstein_small.jpg" #define DEFAULT_FIELD_TYPE "magnitude" #define DEFAULT_FORMAT fileoutStruct::formatImage #define DEFAULT_COLORMAP "brewer" diff --git a/fileout.cu b/fileout.cu index 8fdd091..31e51c7 100644 --- a/fileout.cu +++ b/fileout.cu @@ -113,7 +113,12 @@ void fileoutStruct::saveDetector(microscopeStruct* scope) scalarslice I = scope->getIntensity(); if(is_binary(intFile)) - I.toEnvi(intFile, scope->nf.lambda, append); + { + if(wavenumber) + I.toEnvi(intFile, 10000.0f/scope->nf.lambda, append); + else + I.toEnvi(intFile, scope->nf.lambda, append); + } else I.toImage(intFile); } @@ -123,7 +128,12 @@ void fileoutStruct::saveDetector(microscopeStruct* scope) scalarslice I = scope->getAbsorbance(); if(is_binary(absFile)) - I.toEnvi(absFile, scope->nf.lambda, append); + { + if(wavenumber) + I.toEnvi(absFile, 10000.0f/scope->nf.lambda, append); + else + I.toEnvi(absFile, scope->nf.lambda, append); + } else I.toImage(absFile); } @@ -133,7 +143,12 @@ void fileoutStruct::saveDetector(microscopeStruct* scope) scalarslice I = scope->getTransmittance(); if(is_binary(transFile)) - I.toEnvi(transFile, scope->nf.lambda, append); + { + if(wavenumber) + I.toEnvi(transFile, 10000.0f/scope->nf.lambda, append); + else + I.toEnvi(transFile, scope->nf.lambda, append); + } else I.toImage(transFile); } diff --git a/fileout.h b/fileout.h index efad6d8..c1d4e58 100644 --- a/fileout.h +++ b/fileout.h @@ -28,13 +28,15 @@ struct fileoutStruct{ field_type field; + //flag for output in wavenumber + bool wavenumber; + //image_source source; //color map info rts::colormap::colormapType colormap; ptype colorMax; - void Save(microscopeStruct* scope); void Simulate(microscopeStruct* scope); @@ -42,9 +44,10 @@ struct fileoutStruct{ bool is_binary(std::string filename); void saveNearField(nearfieldStruct* nf); void saveFarField(microscopeStruct* scope); - void saveDetector(microscopeStruct* scope); - + void saveDetector(microscopeStruct* scope); }; + + #endif diff --git a/options.h b/options.h index d16d520..0c7aee8 100644 --- a/options.h +++ b/options.h @@ -138,14 +138,14 @@ static void loadMaterials(po::variables_map vm) for(int i=0; i newM(vm["lambda"].as(), matVec[i], matVec[i+1]); + rts::material newM(SCOPE->nf.lambda, matVec[i], matVec[i+1]); SCOPE->nf.mVector.push_back(newM); } } else { //add the command line material as the default (material 0) - rts::material newM(vm["lambda"].as(), vm["n"].as(), vm["k"].as()); + rts::material newM(SCOPE->nf.lambda, vm["n"].as(), vm["k"].as()); SCOPE->nf.mVector.push_back(newM); } @@ -178,7 +178,7 @@ static void loadNearfieldParams(po::variables_map vm) SCOPE->nf.planeWave = planeWave; //get the wavelength - SCOPE->nf.lambda = vm["lambda"].as(); + //SCOPE->nf.lambda = vm["lambda"].as(); //get the incident field amplitude SCOPE->nf.A = vm["amplitude"].as(); @@ -328,6 +328,7 @@ static void SetOptions(po::options_description &desc) ("material-file,M", po::value< vector >()->multitoken(), "material file:\n [lambda n k]") ("materials", po::value< vector >()->multitoken(), "materials specified using n, k pairs:\n ex. --materials n1 k1 n2 k2\n (if used --n and --k are ignored)") ("lambda,l", po::value()->default_value(DEFAULT_LAMBDA), "incident wavelength") + ("nu", po::value(), "incident frequency (in cm^-1)\n(if specified, lambda is ignored)") ("theta,t", po::value()->default_value(DEFAULT_K_THETA), "light direction (polar coords)") ("phi,p", po::value()->default_value(DEFAULT_K_PHI)) ("fx", po::value()->default_value(DEFAULT_FOCUS_X), "incident focal point") @@ -374,6 +375,20 @@ static void LoadParameters(int argc, char *argv[]) exit(1); } + //load the wavelength + if(vm.count("nu")) + { + //wavelength is given in wavenumber - transform and flag + SCOPE->nf.lambda = 10000/vm["nu"].as(); + gFileOut.wavenumber = true; + } + //otherwise we are using lambda = wavelength + else + { + SCOPE->nf.lambda = vm["lambda"].as(); + gFileOut.wavenumber = false; + } + //load spheres loadSpheres(vm); diff --git a/specimage.py b/specimage.py index 522ccf3..5450289 100755 --- a/specimage.py +++ b/specimage.py @@ -1,61 +1,88 @@ -#!/usr/bin/python3 -import subprocess - -command = "bimsim" - -#images -intImage = "out_i.bmp" -absImage = "out_a" -#detector specs -dsize = 208 -dsample = 6.5 -res = int(dsize / dsample) -sampling = 4 -#incident field -#nu = 800 -source = "central.bmp" -order = 100 -#sphere -x = 0 -y = 0 -z = 0 -a = 5 -#spectral samples -nuStart = 800 -nuStop = 4000 -nuStep = 10 -iters = int((nuStop - nuStart) / nuStep) - - -#set the position of the image plane -command += " -u " + str(-dsize/2) -command += " -v " + str(-dsize/2) -command += " -w " + str(a) -command += " -U " + str(dsize/2) -command += " -V " + str(dsize/2) -command += " -W " + str(a) -command += " --plane-norm-x " + str(0) -command += " --plane-norm-y " + str(0) -command += " --plane-norm-z " + str(1) -command += " -R " + str(res) -command += " --supersample " + str(sampling) -command += " -X " + source -command += " --field-order " + str(order) -command += " -I " + intImage -command += " -A " + absImage -command += " --append" -command += " -x " + str(x) -command += " -y " + str(y) -command += " -z " + str(z) -command += " -r " + str(a) - -for inu in range(0, iters): - print("Iteration # " + str(inu + 1) + "/" + str(iters)) - nu = nuStart + inu * nuStep - lam = 10000.0/nu - - runcommand = command + " -l " + str(lam) - print(runcommand) - subprocess.call(runcommand, shell=True) - -#print("Hello world!") \ No newline at end of file +#!/usr/bin/python3 + + +import subprocess, sys + +nuStart = 800 +nuStop = 4000 +nuStep = 10 + +if len(sys.argv) > 1: + nuStart = sys.argv[1] +if len(sys.argv) > 2: + nuStop = sys.argv[2] +if len(sys.argv) > 3: + nuStep = sys.argv[3] + +command = "bimsim" + +#images +intImage = "out_i.bmp" +absImage = "out_a" +#detector specs +dsize = 208 +dsample = 1 +res = int(dsize / dsample) +sampling = 1 +padding = 1 +#incident field +source = "" +order = 100 +mc = 400 +#sphere +x = 0 +y = 0 +z = 0 +a = 5 +n = 1.4 +k = 0.0 +#spectral samples +iters = int((nuStop - nuStart) / nuStep) + 1 +#optics +NAin = 0.2 +NAout = 0.5 + + +#set the position of the image plane +command += " -u " + str(-dsize/2) +command += " -v " + str(-dsize/2) +command += " -w " + str(a) +command += " -U " + str(dsize/2) +command += " -V " + str(dsize/2) +command += " -W " + str(a) +command += " --plane-norm-x " + str(0) +command += " --plane-norm-y " + str(0) +command += " --plane-norm-z " + str(1) +command += " -R " + str(res) +command += " --supersample " + str(sampling) + +if source != "": + command += " -X " + source +command += " --field-order " + str(order) +command += " -d " + str(padding) +command += " -s " + str(mc) +command += " -I " + intImage +command += " -A " + absImage +command += " --append" +#sphere +command += " -x " + str(x) +command += " -y " + str(y) +command += " -z " + str(z) +command += " -r " + str(a) +command += " -n " + str(n) +command += " -k " + str(k) +#set the optics +command += " -c " + str(NAin) +command += " -C " + str(NAout) +command += " -o " + str(NAin) +command += " -O " + str(NAout) + +for inu in range(0, iters): + print("Iteration # " + str(inu + 1) + "/" + str(iters)) + nu = nuStart + inu * nuStep + + runcommand = command + " --nu " + str(nu) + print(runcommand) + subprocess.call(runcommand, shell=True) + +#print("Hello world!") -- libgit2 0.21.4