Commit 3f36b18e6a3a8120f37e0dcbef0f449498746cf6

Authored by David Mayerich
1 parent 274c3d2b

Adding planewave object

@@ -5,6 +5,11 @@ project(bimsim) @@ -5,6 +5,11 @@ project(bimsim)
5 5
6 #set the module directory 6 #set the module directory
7 set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} "${CMAKE_SOURCE_DIR}") 7 set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} "${CMAKE_SOURCE_DIR}")
  8 +set(CMAKE_AUTOMOC ON)
  9 +
  10 +# As moc files are generated in the binary dir, tell CMake
  11 + # to always look for includes there:
  12 +set(CMAKE_INCLUDE_CURRENT_DIR ON)
8 13
9 #find BOOST 14 #find BOOST
10 set(Boost_USE_STATIC_LIBS ON) 15 set(Boost_USE_STATIC_LIBS ON)
@@ -12,44 +17,34 @@ set(Boost_USE_MULTITHREADED ON) @@ -12,44 +17,34 @@ set(Boost_USE_MULTITHREADED ON)
12 set(Boost_USE_STATIC_RUNTIME OFF) 17 set(Boost_USE_STATIC_RUNTIME OFF)
13 find_package( Boost 1.46.0 COMPONENTS program_options ) 18 find_package( Boost 1.46.0 COMPONENTS program_options )
14 19
15 -#find the Qt4  
16 -find_package(Qt4 REQUIRED) 20 +#find the Qt5
  21 +find_package(Qt5Widgets REQUIRED)
  22 +find_package(Qt5Core REQUIRED)
  23 +find_package(Qt5Gui REQUIRED)
  24 +#find_package(Qt5OpenGL REQUIRED)
17 include_directories(${QT_INCLUDE_DIRECTORY}) 25 include_directories(${QT_INCLUDE_DIRECTORY})
18 -include(${QT_USE_FILE})  
19 26
20 #set up CUDA 27 #set up CUDA
21 find_package(CUDA) 28 find_package(CUDA)
22 29
23 -#find OpenGL  
24 -#find_package(OpenGL REQUIRED)  
25 -  
26 -#find GLUT  
27 -#set(GLUT_ROOT_PATH $ENV{GLUT_ROOT_PATH})  
28 -#find_package(GLUT REQUIRED)  
29 -  
30 -#find GLEW  
31 -#find_package(GLEW REQUIRED)  
32 -  
33 -#add Qt OpenGL stuff  
34 -#set(QT_USE_QTOPENGL TRUE)  
35 -  
36 #ask the user for the RTS location 30 #ask the user for the RTS location
37 -set(RTS_ROOT_PATH $ENV{RTS_ROOT_PATH})  
38 find_package(RTS REQUIRED) 31 find_package(RTS REQUIRED)
39 32
40 #set the include directories 33 #set the include directories
41 include_directories( 34 include_directories(
42 ${CMAKE_CURRENT_BINARY_DIR} 35 ${CMAKE_CURRENT_BINARY_DIR}
43 - ${QT_INCLUDES}  
44 - ${QT_QTOPENGL_INCLUDE_DIR}  
45 -# ${OPENGL_INCLUDE_DIR}  
46 -# ${GLEW_INCLUDE_PATH}  
47 -# ${GLUT_INCLUDE_DIR} 36 + ${Qt5Widgets_INCLUDE_DIRS}
  37 + ${Qt5Core_INCLUDE_DIRS}
  38 + ${Qt5Gui_INCLUDE_DIRS}
  39 +# ${Qt5OpenGL_INCLUDE_DIRS}
48 ${RTS_INCLUDE_DIR} 40 ${RTS_INCLUDE_DIR}
49 - ${Boost_INCLUDE_DIR} 41 + ${Boost_INCLUDE_DIR}
50 ) 42 )
51 43
52 -#enable warnings 44 +#build position independent code for Qt (-fPIC)
  45 +#set(CMAKE_CXX_FLAGS "-fPIC")
  46 +
  47 +#enable warnings (-Wall)
53 if(CMAKE_COMPILER_IS_GNUCC OR CMAKE_COMPILER_IS_GNUCXX) 48 if(CMAKE_COMPILER_IS_GNUCC OR CMAKE_COMPILER_IS_GNUCXX)
54 add_definitions(-Wall) 49 add_definitions(-Wall)
55 endif() 50 endif()
@@ -57,39 +52,39 @@ endif() @@ -57,39 +52,39 @@ endif()
57 #Assign source files to the appropriate variables 52 #Assign source files to the appropriate variables
58 file(GLOB SRC_CPP "*.cpp") 53 file(GLOB SRC_CPP "*.cpp")
59 file(GLOB SRC_H "*.h") 54 file(GLOB SRC_H "*.h")
60 -#file(GLOB SRC_UI "*.ui")  
61 -#file(GLOB SRC_QRC "*.qrc") 55 +file(GLOB SRC_UI "*.ui")
62 file(GLOB SRC_CU "*.cu") 56 file(GLOB SRC_CU "*.cu")
63 57
64 -#set up copying data files  
65 -#configure_file(rtsSU8_k.txt ${CMAKE_CURRENT_BINARY_DIR}/rtsSU8_k.txt @ONLY)  
66 -#configure_file(rtsSU8_n.txt ${CMAKE_CURRENT_BINARY_DIR}/rtsSU8_n.txt @ONLY)  
67 -#configure_file(SurfaceMagnitude.glsl ${CMAKE_CURRENT_BINARY_DIR}/SurfaceMagnitude.glsl @ONLY)  
68 -#configure_file(SurfaceDeform.glsl ${CMAKE_CURRENT_BINARY_DIR}/SurfaceDeform.glsl @ONLY)  
69 -  
70 #determine which source files have to be moc'd 58 #determine which source files have to be moc'd
71 -#Qt4_wrap_cpp(UI_MOC ${SRC_H})  
72 -#Qt4_wrap_ui(UI_H ${SRC_UI})  
73 -#Qt4_add_resources(ALL_RCC ${ALL_QRC}) 59 +Qt5_wrap_cpp(UI_MOC ${SRC_H})
  60 +Qt5_wrap_ui(UI_H ${SRC_UI})
  61 +Qt5_add_resources(ALL_RCC ${ALL_QRC})
74 62
75 -#moc the necessary files  
76 -#Qt4_automoc(${ALL_CPP})  
77 -  
78 -#source_group(QtMoc FILES ${UI_MOC})  
79 -#source_group(QtUI FILES ${SRC_UI}) 63 +set(CMAKE_CXX_FLAGS "${Qt5Widgets_EXECUTABLE_COMPILE_FLAGS}")
80 64
81 #create an executable 65 #create an executable
82 cuda_add_executable(bimsim 66 cuda_add_executable(bimsim
83 ${SRC_CPP} 67 ${SRC_CPP}
84 ${SRC_H} 68 ${SRC_H}
85 -# ${UI_H}  
86 -# ${UI_MOC}  
87 -# ${ALL_RCC} 69 + ${RTS_SOURCE}
  70 + ${UI_H}
  71 + ${SRC_UI}
88 ${SRC_CU} 72 ${SRC_CU}
89 ) 73 )
90 74
  75 +#specify which qt5 modules to use
  76 +qt5_use_modules(bimsim Core Widgets OpenGL Gui)
  77 +
91 #set the link libraries 78 #set the link libraries
92 -target_link_libraries(bimsim ${QT_LIBRARIES} ${QT_QTOPENGL_LIBRARY} ${OPENGL_gl_LIBRARY} ${OPENGL_glu_LIBRARY} ${GLEW_LIBRARY} ${CUDA_cufft_LIBRARY} ${CUDA_cublas_LIBRARY} ${GLUT_glut_LIBRARY} ${Boost_LIBRARIES}) 79 +target_link_libraries(bimsim
  80 + ${Qt5Widgets_LIBRARIES}
  81 + ${Qt5Core_LIBRARIES}
  82 + ${Qt5Gui_LIBRARIES}
  83 +# ${Qt5OpenGL_LIBRARIES}
  84 + ${CUDA_cufft_LIBRARY}
  85 + ${CUDA_cublas_LIBRARY}
  86 + ${Boost_LIBRARIES}
  87 + )
93 88
94 89
95 90
@@ -10,15 +10,14 @@ @@ -10,15 +10,14 @@
10 IF (WIN32) 10 IF (WIN32)
11 FIND_PATH( GLEW_INCLUDE_PATH GL/glew.h 11 FIND_PATH( GLEW_INCLUDE_PATH GL/glew.h
12 $ENV{PROGRAMFILES}/GLEW/include 12 $ENV{PROGRAMFILES}/GLEW/include
13 - ${PROJECT_SOURCE_DIR}/src/nvgl/glew/include  
14 - DOC "The directory where GL/glew.h resides") 13 + $ENV{GLEW_PATH}/include
  14 + )
15 FIND_LIBRARY( GLEW_LIBRARY 15 FIND_LIBRARY( GLEW_LIBRARY
16 NAMES glew GLEW glew32 glew32s 16 NAMES glew GLEW glew32 glew32s
17 PATHS 17 PATHS
18 $ENV{PROGRAMFILES}/GLEW/lib 18 $ENV{PROGRAMFILES}/GLEW/lib
19 - ${PROJECT_SOURCE_DIR}/src/nvgl/glew/bin  
20 - ${PROJECT_SOURCE_DIR}/src/nvgl/glew/lib  
21 - DOC "The GLEW library") 19 + $ENV{GLEW_PATH}/lib/Release/Win32
  20 + )
22 ELSE (WIN32) 21 ELSE (WIN32)
23 FIND_PATH( GLEW_INCLUDE_PATH GL/glew.h 22 FIND_PATH( GLEW_INCLUDE_PATH GL/glew.h
24 /usr/include 23 /usr/include
@@ -42,10 +41,4 @@ IF (GLEW_INCLUDE_PATH) @@ -42,10 +41,4 @@ IF (GLEW_INCLUDE_PATH)
42 SET( GLEW_FOUND 1 CACHE STRING "Set to 1 if GLEW is found, 0 otherwise") 41 SET( GLEW_FOUND 1 CACHE STRING "Set to 1 if GLEW is found, 0 otherwise")
43 ELSE (GLEW_INCLUDE_PATH) 42 ELSE (GLEW_INCLUDE_PATH)
44 SET( GLEW_FOUND 0 CACHE STRING "Set to 1 if GLEW is found, 0 otherwise") 43 SET( GLEW_FOUND 0 CACHE STRING "Set to 1 if GLEW is found, 0 otherwise")
45 -ENDIF (GLEW_INCLUDE_PATH)  
46 -  
47 -MARK_AS_ADVANCED(  
48 - GLEW_FOUND  
49 - GLEW_INCLUDE_PATH  
50 - GLEW_LIBRARY  
51 -)  
52 \ No newline at end of file 44 \ No newline at end of file
  45 +ENDIF (GLEW_INCLUDE_PATH)
53 \ No newline at end of file 46 \ No newline at end of file
1 # Tries to find the RTS include directory 1 # Tries to find the RTS include directory
2 2
3 - FIND_PATH( RTS_INCLUDE_DIR NAMES rts_glShaderProgram.h  
4 - PATHS  
5 - ${CMAKE_CURRENT_SOURCE_DIR}/rts  
6 - ${RTS_ROOT_PATH}  
7 -) 3 +find_path( RTS_INCLUDE_DIR rts_glShaderProgram.h HINTS $ENV{RTS_PATH})
8 4
9 -IF (RTS_FOUND)  
10 - #The following deprecated settings are for backwards compatibility with CMake1.4  
11 - SET (RTS_INCLUDE_PATH ${RTS_INCLUDE_DIR})  
12 -ENDIF(RTS_FOUND) 5 +set (RTS_INCLUDE_DIRS ${RTS_INCLUDE_DIR})
  6 +file (GLOB RTS_SOURCE ${RTS_INCLUDE_DIR}/rts/source/*.cpp)
13 7
14 FIND_PACKAGE_HANDLE_STANDARD_ARGS(RTS REQUIRED_VARS TRUE RTS_INCLUDE_DIR) 8 FIND_PACKAGE_HANDLE_STANDARD_ARGS(RTS REQUIRED_VARS TRUE RTS_INCLUDE_DIR)
@@ -32,10 +32,10 @@ extern bool verbose; @@ -32,10 +32,10 @@ extern bool verbose;
32 #include "rts/math/point.h" 32 #include "rts/math/point.h"
33 #include "rts/math/quad.h" 33 #include "rts/math/quad.h"
34 34
35 -typedef rts::rtsComplex<ptype> bsComplex;  
36 -typedef rts::rtsVector<ptype, 3> bsVector;  
37 -typedef rts::rtsPoint<ptype, 3> bsPoint;  
38 -typedef rts::rtsQuad<ptype, 3> bsRect; 35 +typedef rts::complex<ptype> bsComplex;
  36 +typedef rts::vector<ptype, 3> bsVector;
  37 +typedef rts::point<ptype, 3> bsPoint;
  38 +typedef rts::quad<ptype, 3> bsRect;
39 39
40 40
41 #endif 41 #endif
1 #include "fileout.h" 1 #include "fileout.h"
2 -//#include "scalarfield.h"  
3 -  
4 -  
5 -/*void fileoutStruct::saveMag(fieldslice* U, std::string filename, rts::colormap::colormapType cmap)  
6 -{  
7 - int Rx = U->R[0];  
8 - int Ry = U->R[1];  
9 -  
10 - //allocate space for a scalar field on the GPU  
11 - ptype* gpuScalar;  
12 - int memsize = sizeof(ptype) * Rx * Ry;  
13 - HANDLE_ERROR(cudaMalloc((void**) &gpuScalar, memsize));  
14 - HANDLE_ERROR(cudaMemset(gpuScalar, 0, memsize));  
15 - U->Mag(gpuScalar);  
16 -  
17 -  
18 - rts::colormap::gpu2image<ptype>(gpuScalar, filename, Rx, Ry, 0, colorMax, cmap);  
19 -  
20 - HANDLE_ERROR(cudaFree(gpuScalar));  
21 -}  
22 -  
23 -void fileoutStruct::saveReal_scalar(fieldslice* U, std::string filename, rts::colormap::colormapType cmap)  
24 -{  
25 - //returns the real component  
26 - scalarslice sf = U->Real();  
27 - sf.toImage(filename, false, cmap);  
28 -  
29 -}  
30 -  
31 -void fileoutStruct::saveImag_scalar(fieldslice* U, std::string filename, rts::colormap::colormapType cmap)  
32 -{  
33 - //returns the imaginary component of a field (assumed scalar)  
34 - scalarslice sf = U->Imag();  
35 - sf.toImage(filename, false, cmap);  
36 -}  
37 -  
38 -void fileoutStruct::saveIntensity(fieldslice* U, std::string filename, rts::colormap::colormapType cmap)  
39 -{  
40 - //get the intensity of the field  
41 - scalarslice sf = U->Intensity();  
42 - sf.toImage(filename, true, cmap);  
43 -}  
44 -  
45 -void fileoutStruct::saveAngularSpectrum(fieldslice* U, std::string filename, rts::colormap::colormapType cmap)  
46 -{  
47 - ptype* gpuScalar;  
48 - int memsize = sizeof(ptype) * U->R[0] * U->R[1];  
49 - HANDLE_ERROR(cudaMalloc((void**) &gpuScalar, memsize));  
50 - HANDLE_ERROR(cudaMemset(gpuScalar, 0, memsize));  
51 -  
52 - //convert the field slice to its angular spectrum  
53 - U->toAngularSpectrum();  
54 -  
55 - //convert the angular spectrum to a scalar field  
56 - U->Mag(gpuScalar);  
57 -  
58 - //save the color image  
59 - rts::colormap::gpu2image<ptype>(gpuScalar, filename, U->R[0], U->R[1], 0, colorMax, cmap);  
60 -  
61 - HANDLE_ERROR(cudaFree(gpuScalar));  
62 -  
63 -}*/  
64 2
65 void fileoutStruct::saveNearField(nearfieldStruct* nf) 3 void fileoutStruct::saveNearField(nearfieldStruct* nf)
66 { 4 {
@@ -221,7 +159,7 @@ void fileoutStruct::Save(microscopeStruct* scope) @@ -221,7 +159,7 @@ void fileoutStruct::Save(microscopeStruct* scope)
221 159
222 } 160 }
223 161
224 - 162 +
225 163
226 164
227 } 165 }
@@ -5,7 +5,7 @@ @@ -5,7 +5,7 @@
5 //#include "defaults.h" 5 //#include "defaults.h"
6 #include "dataTypes.h" 6 #include "dataTypes.h"
7 7
8 -#include "rts/graphics/colormap.h" 8 +#include "rts/visualization/colormap.h"
9 #include "fieldslice.h" 9 #include "fieldslice.h"
10 #include "nearfield.h" 10 #include "nearfield.h"
11 #include "microscope.h" 11 #include "microscope.h"
@@ -23,6 +23,12 @@ microscopeStruct* SCOPE; @@ -23,6 +23,12 @@ microscopeStruct* SCOPE;
23 23
24 #include "warnings.h" 24 #include "warnings.h"
25 25
  26 +#include "planewave.h"
  27 +
  28 +//user interface
  29 +#include "qtMainDialog.h"
  30 +bool gui = false;
  31 +
26 fileoutStruct gFileOut; 32 fileoutStruct gFileOut;
27 bool verbose = false; 33 bool verbose = false;
28 using namespace std; 34 using namespace std;
@@ -33,6 +39,15 @@ int cbessjyva(double v,complex&lt;double&gt; z,double &amp;vm,complex&lt;double&gt;*cjv, @@ -33,6 +39,15 @@ int cbessjyva(double v,complex&lt;double&gt; z,double &amp;vm,complex&lt;double&gt;*cjv,
33 int main(int argc, char *argv[]) 39 int main(int argc, char *argv[])
34 { 40 {
35 41
  42 + //benchtest planewave class
  43 + rts::vector<ptype, 3> k(1, 0, 0);
  44 + rts::vector<ptype, 3> E(2, 2, 0);
  45 + planewave<ptype> P(k, E);
  46 +
  47 + std::cout<<P<<std::endl;
  48 +
  49 + exit(1);
  50 +
36 SCOPE = new microscopeStruct(); 51 SCOPE = new microscopeStruct();
37 52
38 LoadParameters(argc, argv); 53 LoadParameters(argc, argv);
@@ -40,12 +55,26 @@ int main(int argc, char *argv[]) @@ -40,12 +55,26 @@ int main(int argc, char *argv[])
40 //initialize GPU memory for fields 55 //initialize GPU memory for fields
41 SCOPE->init(); 56 SCOPE->init();
42 57
43 - gFileOut.Save(SCOPE); 58 + if(gui)
  59 + {
  60 + QApplication app(argc, argv);
  61 + qtMainDialog bsDialog;
  62 +
  63 + //populate the user interface with the default and command-line values
  64 + bsDialog.populateUi();
  65 + bsDialog.show();
  66 + return app.exec();
  67 + }
  68 + else
  69 + {
  70 +
  71 + gFileOut.Save(SCOPE);
44 72
45 - if(verbose)  
46 - OutputOptions(); 73 + if(verbose)
  74 + OutputOptions();
47 75
48 - SCOPE->destroy(); 76 + SCOPE->destroy();
  77 + }
49 78
50 79
51 80
@@ -4,7 +4,7 @@ @@ -4,7 +4,7 @@
4 #include "rts/tools/progressbar.h" 4 #include "rts/tools/progressbar.h"
5 #include "rts/cuda/timer.h" 5 #include "rts/cuda/timer.h"
6 #include "dataTypes.h" 6 #include "dataTypes.h"
7 -#include "rts/graphics/colormap.h" 7 +#include "rts/visualization/colormap.h"
8 8
9 #include <QImage> 9 #include <QImage>
10 10
@@ -49,6 +49,8 @@ microscopeStruct::microscopeStruct() @@ -49,6 +49,8 @@ microscopeStruct::microscopeStruct()
49 49
50 void microscopeStruct::init() 50 void microscopeStruct::init()
51 { 51 {
  52 +
  53 + nf.scalarSim = scalarSim;
52 //Ud.scalarField = scalarSim; 54 //Ud.scalarField = scalarSim;
53 //Ufd.scalarField = scalarSim; 55 //Ufd.scalarField = scalarSim;
54 56
@@ -137,7 +139,7 @@ void microscopeStruct::clearDetector() @@ -137,7 +139,7 @@ void microscopeStruct::clearDetector()
137 //flag for a vector simulation 139 //flag for a vector simulation
138 void microscopeStruct::setPos(bsPoint pMin, bsPoint pMax, bsVector normal) 140 void microscopeStruct::setPos(bsPoint pMin, bsPoint pMax, bsVector normal)
139 { 141 {
140 - pos = rts::rtsQuad<ptype, 3>(pMin, pMax, normal); 142 + pos = rts::quad<ptype, 3>(pMin, pMax, normal);
141 } 143 }
142 144
143 void microscopeStruct::setRes(int x_res, int y_res, int pad, int supersampling) 145 void microscopeStruct::setRes(int x_res, int y_res, int pad, int supersampling)
@@ -31,7 +31,7 @@ struct microscopeStruct @@ -31,7 +31,7 @@ struct microscopeStruct
31 ptype objective[2]; 31 ptype objective[2];
32 32
33 //image position and orientation in world space 33 //image position and orientation in world space
34 - rts::rtsQuad<ptype, 3> pos; 34 + rts::quad<ptype, 3> pos;
35 35
36 vector<sourcePoint> focalPoints; 36 vector<sourcePoint> focalPoints;
37 37
@@ -20,15 +20,15 @@ void mcSampleNA(bsVector* samples, int N, bsVector k, ptype NAin, ptype NAout) @@ -20,15 +20,15 @@ void mcSampleNA(bsVector* samples, int N, bsVector k, ptype NAin, ptype NAout)
20 //get the axis of rotation for transforming (0, 0, 1) to k 20 //get the axis of rotation for transforming (0, 0, 1) to k
21 //k = -k; 21 //k = -k;
22 ptype cos_angle = k.dot(bsVector(0, 0, 1)); 22 ptype cos_angle = k.dot(bsVector(0, 0, 1));
23 - rts::rtsMatrix<ptype, 3> rotation; 23 + rts::matrix<ptype, 3> rotation;
24 if(cos_angle != 1.0) 24 if(cos_angle != 1.0)
25 { 25 {
26 bsVector axis = bsVector(0, 0, 1).cross(k).norm(); 26 bsVector axis = bsVector(0, 0, 1).cross(k).norm();
27 27
28 ptype angle = acos(cos_angle); 28 ptype angle = acos(cos_angle);
29 - rts::rtsQuaternion<ptype> quat; 29 + rts::quaternion<ptype> quat;
30 quat.CreateRotation(angle, axis); 30 quat.CreateRotation(angle, axis);
31 - rotation = quat.toMatrix(); 31 + rotation = quat.toMatrix3();
32 } 32 }
33 33
34 //find the phi values associated with the cassegrain ring 34 //find the phi values associated with the cassegrain ring
@@ -39,7 +39,7 @@ void nearfieldStruct::destroy() @@ -39,7 +39,7 @@ void nearfieldStruct::destroy()
39 39
40 void nearfieldStruct::setPos(bsPoint pMin, bsPoint pMax, bsVector normal) 40 void nearfieldStruct::setPos(bsPoint pMin, bsPoint pMax, bsVector normal)
41 { 41 {
42 - pos = rts::rtsQuad<ptype, 3>(pMin, pMax, normal); 42 + pos = rts::quad<ptype, 3>(pMin, pMax, normal);
43 } 43 }
44 44
45 void nearfieldStruct::setRes(int x_res, int y_res) 45 void nearfieldStruct::setRes(int x_res, int y_res)
@@ -115,7 +115,7 @@ void nearfieldStruct::calcSpheres() @@ -115,7 +115,7 @@ void nearfieldStruct::calcSpheres()
115 115
116 //set the refractive index for the sphere 116 //set the refractive index for the sphere
117 int imat = sVector[i].iMaterial; 117 int imat = sVector[i].iMaterial;
118 - rts::rtsComplex<ptype> n = mVector[imat](lambda); 118 + rts::complex<ptype> n = mVector[imat](lambda);
119 119
120 //calculate the scattering coefficients 120 //calculate the scattering coefficients
121 sVector[i].calcCoeff(lambda, n); 121 sVector[i].calcCoeff(lambda, n);
@@ -137,19 +137,47 @@ void nearfieldStruct::calcSpheres() @@ -137,19 +137,47 @@ void nearfieldStruct::calcSpheres()
137 void nearfieldStruct::calcUs() 137 void nearfieldStruct::calcUs()
138 { 138 {
139 139
140 -  
141 - if(lut_us)  
142 - scalarUpLut(); 140 + if(scalarSim)
  141 + {
  142 + if(lut_us)
  143 + scalarUpLut();
  144 + else
  145 + scalarUs();
  146 +
  147 + if(lut_us)
  148 + std::cout<<"Using LUT for Us"<<std::endl;
  149 + std::cout<<"Calculate Us Scalar Sim."<<std::endl;
  150 + }
143 else 151 else
144 - scalarUs(); 152 + {
  153 + std::cout<<"Calculate Us Vector Sim."<<std::endl;
  154 + }
145 } 155 }
146 156
147 void nearfieldStruct::calcUf() 157 void nearfieldStruct::calcUf()
148 { 158 {
149 - if(lut_uf)  
150 - scalarUfLut();  
151 - else  
152 - scalarUf(); 159 + if(scalarSim)
  160 + {
  161 + if(lut_uf)
  162 + scalarUfLut();
  163 + else
  164 + scalarUf();
  165 +
  166 + //if(lut_uf)
  167 + // std::cout<<"Using LUT for Uf"<<std::endl;
  168 + //std::cout<<"Calculate Us Scalar Sim."<<std::endl;
  169 + }
  170 + else
  171 + {
  172 + std::cout<<"Calculating Uf Vector Sim..."<<std::endl;
  173 +
  174 + if(lut_uf)
  175 + vectorUfLut();
  176 + else
  177 + vectorUf();
  178 +
  179 +
  180 + }
153 } 181 }
154 182
155 void nearfieldStruct::Simulate() 183 void nearfieldStruct::Simulate()
@@ -22,12 +22,15 @@ struct nearfieldStruct @@ -22,12 +22,15 @@ struct nearfieldStruct
22 //amplitude of the incident field 22 //amplitude of the incident field
23 ptype A; 23 ptype A;
24 24
  25 + //incident field polarization;
  26 + bsVector E;
  27 +
25 //position of the focus in 3D space 28 //position of the focus in 3D space
26 bsVector k; //cartesian coordinates, normalized 29 bsVector k; //cartesian coordinates, normalized
27 bsPoint focus; 30 bsPoint focus;
28 31
29 //slice position and orientation in world space 32 //slice position and orientation in world space
30 - rts::rtsQuad<ptype, 3> pos; 33 + rts::quad<ptype, 3> pos;
31 34
32 //slices for the focused field 35 //slices for the focused field
33 fieldslice Uf; 36 fieldslice Uf;
@@ -89,9 +92,15 @@ struct nearfieldStruct @@ -89,9 +92,15 @@ struct nearfieldStruct
89 92
90 //this function re-computes the focused field 93 //this function re-computes the focused field
91 void calcUf(); 94 void calcUf();
  95 +
  96 + //scalar functions
92 void scalarUf(); 97 void scalarUf();
93 void scalarUfLut(); 98 void scalarUfLut();
94 99
  100 + //vector functions
  101 + void vectorUf();
  102 + void vectorUfLut();
  103 +
95 void calcBesselLut(ptype* j, ptype d_min, ptype d_max, int dR); 104 void calcBesselLut(ptype* j, ptype d_min, ptype d_max, int dR);
96 105
97 //compute the field scattered by all of the materials 106 //compute the field scattered by all of the materials
@@ -178,6 +178,8 @@ void nearfieldStruct::scalarUf() @@ -178,6 +178,8 @@ void nearfieldStruct::scalarUf()
178 else 178 else
179 { 179 {
180 //pre-compute the cosine of the obscuration and objective angles 180 //pre-compute the cosine of the obscuration and objective angles
  181 + //cout<<"Condenser angle in: "<<asin(condenser[0])<<std::endl;
  182 + //cout<<"Condenser angle out: "<<asin(condenser[1])<<std::endl;
181 ptype cosAlpha = cos(asin(condenser[0])); 183 ptype cosAlpha = cos(asin(condenser[0]));
182 ptype cosBeta = cos(asin(condenser[1])); 184 ptype cosBeta = cos(asin(condenser[1]));
183 //compute the scalar Uf field (this will be in the x_hat channel of Uf) 185 //compute the scalar Uf field (this will be in the x_hat channel of Uf)
nfVectorUf.cu 0 โ†’ 100644
  1 +#include "nearfield.h"
  2 +#include "rts/math/spherical_bessel.h"
  3 +#include "rts/math/legendre.h"
  4 +#include <stdlib.h>
  5 +#include "rts/cuda/error.h"
  6 +#include "rts/cuda/timer.h"
  7 +
  8 +//Incident field for a single plane wave
  9 +__global__ void gpuVectorUfp(bsComplex* Uf, bsVector k, ptype kmag, bsPoint f, ptype A, bsRect ABCD, int uR, int vR)
  10 +{
  11 + /*Compute the scalar focused field using Debye focusing
  12 + k = direction of focused light, where |k| = 2*pi/lambda
  13 + P = rect struct describing the field slice
  14 + rX, rY = resolution of the field slice
  15 + cNAin = inner NA of the condenser
  16 + cNAout = outer NA of the condenser
  17 + */
  18 +
  19 + //get the current coordinate in the plane slice
  20 + int iu = blockIdx.x * blockDim.x + threadIdx.x;
  21 + int iv = blockIdx.y * blockDim.y + threadIdx.y;
  22 +
  23 + //make sure that the thread indices are in-bounds
  24 + if(iu >= uR || iv >= vR) return;
  25 +
  26 + //compute the index (easier access to the scalar field array)
  27 + int i = iv*uR + iu;
  28 +
  29 + //compute the parameters for u and v
  30 + ptype u = (ptype)iu / uR;
  31 + ptype v = (ptype)iv / vR;
  32 +
  33 + //get the rtsPoint in world space and then the r vector
  34 + bsPoint p = ABCD(u, v);
  35 + bsVector r = p - f;
  36 + //ptype d = r.len();
  37 +
  38 + ptype k_dot_r = kmag * k.dot(r);
  39 + bsComplex d(0, k_dot_r);
  40 +
  41 + Uf[i] = exp(d) * A;
  42 +
  43 +}
  44 +
  45 +//Incident field for a focused point source
  46 +__global__ void gpuVectorUf(bsComplex* Uf, bsVector k, ptype kmag, bsPoint f, ptype A, bsRect ABCD, int uR, int vR, ptype cosAlpha, ptype cosBeta, int nl, ptype j_conv = 1.4)
  47 +{
  48 + //Compute the scalar focused field using Debye focusing
  49 + // k = direction of focused light, where |k| = 2*pi/lambda
  50 + // P = rect struct describing the field slice
  51 + // rX, rY = resolution of the field slice
  52 + // cNAin = inner NA of the condenser
  53 + // cNAout = outer NA of the condenser
  54 +
  55 +
  56 + //get the current coordinate in the plane slice
  57 + int iu = blockIdx.x * blockDim.x + threadIdx.x;
  58 + int iv = blockIdx.y * blockDim.y + threadIdx.y;
  59 +
  60 + //make sure that the thread indices are in-bounds
  61 + if(iu >= uR || iv >= vR) return;
  62 +
  63 + //compute the index (easier access to the scalar field array)
  64 + int i = iv*uR + iu;
  65 +
  66 + //compute the parameters for u and v
  67 + ptype u = (ptype)iu / (uR);
  68 + ptype v = (ptype)iv / (vR);
  69 +
  70 + //get the rtsPoint in world space and then the r vector
  71 + bsPoint p = ABCD(u, v);
  72 + bsVector r = p - f;
  73 + ptype d = r.len();
  74 + if(d < EPSILON_FLOAT)
  75 + {
  76 + Uf[i] = A * 2 * PI * (cosAlpha - cosBeta);
  77 + return;
  78 + }
  79 +
  80 + //get info for the light direction and frequency
  81 + //k = k.norm();
  82 + r = r.norm();
  83 +
  84 + //compute the imaginary factor i^l
  85 + bsComplex im = bsComplex(0, 1);
  86 + bsComplex il = bsComplex(1, 0);
  87 +
  88 + //Bessel and Legendre functions are computed dynamically to save memory
  89 + //initialize the Bessel and Legendre functions
  90 + ptype j[2];
  91 + ptype kd = kmag * d;
  92 + rts::init_sbesselj<ptype>(kd, j);
  93 +
  94 + ptype P[2];
  95 + //get the angle between k and r (light direction and position vector)
  96 + ptype cosTheta;
  97 + cosTheta = k.dot(r);
  98 +
  99 + //deal with the degenerate case where r == 0
  100 + //if(isnan(cosTheta))
  101 + // cosTheta = 0;
  102 + rts::init_legendre<ptype>(cosTheta, P[0], P[1]);
  103 +
  104 + //initialize legendre functions for the cassegrain angles
  105 + ptype Palpha[3];
  106 + //ptype cosAlpha = cos(asin(cNAin));
  107 + rts::init_legendre<ptype>(cosAlpha, Palpha[0], Palpha[1]);
  108 + Palpha[2] = 1;
  109 +
  110 + ptype Pbeta[3];
  111 + //ptype cosBeta = cos(asin(cNAout));
  112 + rts::init_legendre<ptype>(cosBeta, Pbeta[0], Pbeta[1]);
  113 + Pbeta[2] = 1;
  114 +
  115 + //for each order l
  116 + bsComplex sumUf(0.0, 0.0);
  117 + ptype jl = 0.0;
  118 + ptype Pl;
  119 + for(int l = 0; l<=nl; l++)
  120 + {
  121 +
  122 + if(l==0)
  123 + {
  124 +
  125 + jl = j[0];
  126 + Pl = P[0];
  127 + }
  128 + else if(l==1)
  129 + {
  130 + jl = j[1];
  131 + Pl = P[1];
  132 +
  133 + //adjust the cassegrain Legendre function
  134 + Palpha[2] = Palpha[0];
  135 + rts::shift_legendre<ptype>(l+1, cosAlpha, Palpha[0], Palpha[1]);
  136 + Pbeta[2] = Pbeta[0];
  137 + rts::shift_legendre<ptype>(l+1, cosBeta, Pbeta[0], Pbeta[1]);
  138 + }
  139 + else
  140 + {
  141 + rts::shift_sbesselj<ptype>(l, kd, j);//, j_conv);
  142 + rts::shift_legendre<ptype>(l, cosTheta, P[0], P[1]);
  143 +
  144 + jl = j[1];
  145 + Pl = P[1];
  146 +
  147 + //adjust the cassegrain outer Legendre function
  148 + Palpha[2] = Palpha[0];
  149 + rts::shift_legendre<ptype>(l+1, cosAlpha, Palpha[0], Palpha[1]);
  150 + Pbeta[2] = Pbeta[0];
  151 + rts::shift_legendre<ptype>(l+1, cosBeta, Pbeta[0], Pbeta[1]);
  152 + }
  153 +
  154 + sumUf += il * jl * Pl * (Palpha[1] - Palpha[2] - Pbeta[1] + Pbeta[2]);
  155 +
  156 + il *= im;
  157 + }
  158 +
  159 + Uf[i] = sumUf * 2 * PI * A;
  160 +
  161 +}
  162 +
  163 +
  164 +void nearfieldStruct::vectorUf()
  165 +{
  166 +
  167 +
  168 + gpuStartTimer();
  169 +
  170 + //create one thread for each pixel of the field slice
  171 + dim3 dimBlock(SQRT_BLOCK, SQRT_BLOCK);
  172 + dim3 dimGrid((Uf.R[0] + SQRT_BLOCK -1)/SQRT_BLOCK, (Uf.R[1] + SQRT_BLOCK - 1)/SQRT_BLOCK);
  173 +
  174 + //if we are computing a plane wave, call the gpuScalarUfp function
  175 + if(planeWave)
  176 + {
  177 + std::cout<<"Calculating vector plane wave..."<<std::endl;
  178 + gpuVectorUfp<<<dimGrid, dimBlock>>>(Uf.x_hat, k, 2 * PI / lambda, focus, A, pos, Uf.R[0], Uf.R[1]);
  179 + }
  180 + //otherwise compute the condenser info and create a focused field
  181 + else
  182 + {
  183 + //pre-compute the cosine of the obscuration and objective angles
  184 + ptype cosAlpha = cos(asin(condenser[0]));
  185 + ptype cosBeta = cos(asin(condenser[1]));
  186 + //compute the scalar Uf field (this will be in the x_hat channel of Uf)
  187 + gpuVectorUf<<<dimGrid, dimBlock>>>(Uf.x_hat, k, 2 * PI / lambda, focus, A, pos, Uf.R[0], Uf.R[1], cosAlpha, cosBeta, m);
  188 + }
  189 +
  190 + t_Uf = gpuStopTimer();
  191 +}
nfVectorUfLut.cu 0 โ†’ 100644
  1 +#include "nearfield.h"
  2 +
  3 +#include "rts/math/legendre.h"
  4 +#include "rts/cuda/error.h"
  5 +#include "rts/cuda/timer.h"
  6 +
  7 +void nearfieldStruct::vectorUfLut()
  8 +{
  9 +
  10 +}
@@ -5,7 +5,7 @@ @@ -5,7 +5,7 @@
5 5
6 #include "nearfield.h" 6 #include "nearfield.h"
7 #include "microscope.h" 7 #include "microscope.h"
8 -#include "rts/graphics/colormap.h" 8 +#include "rts/visualization/colormap.h"
9 #include "fileout.h" 9 #include "fileout.h"
10 //extern nearfieldStruct* NF; 10 //extern nearfieldStruct* NF;
11 extern microscopeStruct* SCOPE; 11 extern microscopeStruct* SCOPE;
@@ -24,11 +24,18 @@ using namespace std; @@ -24,11 +24,18 @@ using namespace std;
24 namespace po = boost::program_options; 24 namespace po = boost::program_options;
25 25
26 extern bool verbose; 26 extern bool verbose;
  27 +extern bool gui;
27 28
28 29
29 30
30 static void lNearfield(po::variables_map vm) 31 static void lNearfield(po::variables_map vm)
31 { 32 {
  33 + //test to see if we are running a vector field simulation
  34 + bool vectorField = false;
  35 + if(vm.count("vector"))
  36 + vectorField = true;
  37 + SCOPE->scalarSim = !vectorField;
  38 +
32 //test to see if we are simulating a plane wave 39 //test to see if we are simulating a plane wave
33 bool planeWave = DEFAULT_PLANEWAVE; 40 bool planeWave = DEFAULT_PLANEWAVE;
34 if(vm.count("plane-wave")) 41 if(vm.count("plane-wave"))
@@ -176,6 +183,10 @@ void lFlags(po::variables_map vm, po::options_description desc) @@ -176,6 +183,10 @@ void lFlags(po::variables_map vm, po::options_description desc)
176 { 183 {
177 SCOPE->nf.lut_uf = true; 184 SCOPE->nf.lut_uf = true;
178 } 185 }
  186 +
  187 + //gui
  188 + if(vm.count("gui"))
  189 + gui = true;
179 } 190 }
180 191
181 void lWavelength(po::variables_map vm) 192 void lWavelength(po::variables_map vm)
@@ -315,6 +326,8 @@ void lSpheres(po::variables_map vm) @@ -315,6 +326,8 @@ void lSpheres(po::variables_map vm)
315 { 326 {
316 //otherwise output an error 327 //otherwise output an error
317 cout<<"BIMSIM Error - A material is not loaded for sphere "<<s+1<<"."<<endl; 328 cout<<"BIMSIM Error - A material is not loaded for sphere "<<s+1<<"."<<endl;
  329 + cout<<"Material requested: "<<SCOPE->nf.sVector[s].iMaterial + 1<<endl;
  330 + cout<<"Number of materials: "<<SCOPE->nf.mVector.size()<<endl;
318 exit(1); 331 exit(1);
319 } 332 }
320 } 333 }
@@ -353,13 +366,13 @@ static void lMaterials(po::variables_map vm) @@ -353,13 +366,13 @@ static void lMaterials(po::variables_map vm)
353 for(unsigned int i=0; i<filenames.size(); i++) 366 for(unsigned int i=0; i<filenames.size(); i++)
354 { 367 {
355 //load the file into a string 368 //load the file into a string
356 - std::ifstream ifs(filenames[i].c_str()); 369 + //std::ifstream ifs(filenames[i].c_str());
357 370
358 - std::string instr((std::istreambuf_iterator<char>(ifs)), std::istreambuf_iterator<char>()); 371 + //std::string instr((std::istreambuf_iterator<char>(ifs)), std::istreambuf_iterator<char>());
359 372
360 //load the list of spheres from a string 373 //load the list of spheres from a string
361 - rts::material<ptype> newM;  
362 - newM.fromStr(instr, ""); 374 + rts::material<ptype> newM(filenames[i].c_str());
  375 + //newM.fromStr(instr, "");
363 SCOPE->nf.mVector.push_back(newM); 376 SCOPE->nf.mVector.push_back(newM);
364 } 377 }
365 } 378 }
@@ -497,19 +510,24 @@ static void SetOptions(po::options_description &amp;desc) @@ -497,19 +510,24 @@ static void SetOptions(po::options_description &amp;desc)
497 { 510 {
498 desc.add_options() 511 desc.add_options()
499 ("help", "prints this help") 512 ("help", "prints this help")
500 - ("verbose", "verbose output\n") 513 + ("gui", "run using the Qt GUI")
  514 + ("verbose", "verbose output\n\nOutput Parameters\n--------------------------")
501 515
  516 + ("vector", "run a vector field simulation")
502 ("intensity", po::value<string>()->default_value(DEFAULT_INTENSITY_FILE), "output measured intensity (filename)") 517 ("intensity", po::value<string>()->default_value(DEFAULT_INTENSITY_FILE), "output measured intensity (filename)")
503 ("absorbance", po::value<string>()->default_value(DEFAULT_ABSORBANCE_FILE), "output measured absorbance (filename)") 518 ("absorbance", po::value<string>()->default_value(DEFAULT_ABSORBANCE_FILE), "output measured absorbance (filename)")
504 ("transmittance", po::value<string>()->default_value(DEFAULT_TRANSMITTANCE_FILE), "output measured transmittance (filename)") 519 ("transmittance", po::value<string>()->default_value(DEFAULT_TRANSMITTANCE_FILE), "output measured transmittance (filename)")
505 ("far-field", po::value<string>()->default_value(DEFAULT_FAR_FILE), "output far-field at detector (filename)") 520 ("far-field", po::value<string>()->default_value(DEFAULT_FAR_FILE), "output far-field at detector (filename)")
506 ("near-field", po::value<string>()->default_value(DEFAULT_NEAR_FILE), "output field at focal plane (filename)") 521 ("near-field", po::value<string>()->default_value(DEFAULT_NEAR_FILE), "output field at focal plane (filename)")
507 - ("extended-source", po::value<string>()->default_value(DEFAULT_EXTENDED_SOURCE), "image of source at focus (filename)\n") 522 + ("extended-source", po::value<string>()->default_value(DEFAULT_EXTENDED_SOURCE), "image of source at focus (filename)")
  523 + ("output-type", po::value<string>()->default_value(DEFAULT_FIELD_TYPE), "output field value:\n magnitude, polarization, real, imaginary, angular-spectrum")
  524 + ("colormap", po::value<string>()->default_value(DEFAULT_COLORMAP), "colormap: gray, brewer")
  525 + ("append", "append result to an existing file\n (binary files only)\n\nSphere Parameters\n--------------------------")
508 526
509 ("spheres", po::value< vector<ptype> >()->multitoken(), "sphere position: x y z a m") 527 ("spheres", po::value< vector<ptype> >()->multitoken(), "sphere position: x y z a m")
510 ("sphere-file", po::value< vector<string> >()->multitoken(), "sphere file:\n [x y z radius material]") 528 ("sphere-file", po::value< vector<string> >()->multitoken(), "sphere file:\n [x y z radius material]")
511 ("materials", po::value< vector<ptype> >()->multitoken(), "refractive indices as n, k pairs:\n ex. -m n0 k0 n1 k1 n2 k2") 529 ("materials", po::value< vector<ptype> >()->multitoken(), "refractive indices as n, k pairs:\n ex. -m n0 k0 n1 k1 n2 k2")
512 - ("material-file", po::value< vector<string> >()->multitoken(), "material file:\n [lambda n k]\n") 530 + ("material-file", po::value< vector<string> >()->multitoken(), "material file:\n [lambda n k]\n\nOptics\n--------------------------")
513 531
514 ("lambda", po::value<ptype>()->default_value(DEFAULT_LAMBDA), "incident wavelength") 532 ("lambda", po::value<ptype>()->default_value(DEFAULT_LAMBDA), "incident wavelength")
515 ("nu", po::value<ptype>(), "incident frequency (in cm^-1)\n(if specified, lambda is ignored)") 533 ("nu", po::value<ptype>(), "incident frequency (in cm^-1)\n(if specified, lambda is ignored)")
@@ -518,7 +536,7 @@ static void SetOptions(po::options_description &amp;desc) @@ -518,7 +536,7 @@ static void SetOptions(po::options_description &amp;desc)
518 ("condenser", po::value< vector<ptype> >()->multitoken(), "condenser numerical aperature\nA pair of values can be used to specify an inner obscuration: -c NAin NAout") 536 ("condenser", po::value< vector<ptype> >()->multitoken(), "condenser numerical aperature\nA pair of values can be used to specify an inner obscuration: -c NAin NAout")
519 ("objective", po::value< vector<ptype> >()->multitoken(), "objective numerical aperature\nA pair of values can be used to specify an inner obscuration: -c NAin NAout") 537 ("objective", po::value< vector<ptype> >()->multitoken(), "objective numerical aperature\nA pair of values can be used to specify an inner obscuration: -c NAin NAout")
520 ("focus", po::value< vector<ptype> >()->multitoken(), "focal position for the incident point source\n (default = --focus 0 0 0)") 538 ("focus", po::value< vector<ptype> >()->multitoken(), "focal position for the incident point source\n (default = --focus 0 0 0)")
521 - ("plane-wave", "simulates an incident plane wave\n") 539 + ("plane-wave", "simulates an incident plane wave\n\n\nImaging Parameters\n--------------------------")
522 540
523 ("resolution", po::value<unsigned int>()->default_value(DEFAULT_SLICE_RES), "resolution of the detector") 541 ("resolution", po::value<unsigned int>()->default_value(DEFAULT_SLICE_RES), "resolution of the detector")
524 ("plane-lower-left", po::value< vector<ptype> >()->multitoken(), "lower-left position of the image plane") 542 ("plane-lower-left", po::value< vector<ptype> >()->multitoken(), "lower-left position of the image plane")
@@ -526,7 +544,7 @@ static void SetOptions(po::options_description &amp;desc) @@ -526,7 +544,7 @@ static void SetOptions(po::options_description &amp;desc)
526 ("plane-normal", po::value< vector<ptype> >()->multitoken(), "normal for the image plane") 544 ("plane-normal", po::value< vector<ptype> >()->multitoken(), "normal for the image plane")
527 ("xy", po::value< vector<ptype> >()->multitoken(), "specify an x-y image plane\n (standard microscope)") 545 ("xy", po::value< vector<ptype> >()->multitoken(), "specify an x-y image plane\n (standard microscope)")
528 ("xz", po::value< vector<ptype> >()->multitoken(), "specify a x-z image plane\n (cross-section of the focal volume)") 546 ("xz", po::value< vector<ptype> >()->multitoken(), "specify a x-z image plane\n (cross-section of the focal volume)")
529 - ("yz", po::value< vector<ptype> >()->multitoken(), "specify a y-z image plane\n (cross-section of the focal volume)\n") 547 + ("yz", po::value< vector<ptype> >()->multitoken(), "specify a y-z image plane\n (cross-section of the focal volume)\n\nSampling Parameters\n--------------------------")
530 548
531 ("samples", po::value<int>()->default_value(DEFAULT_SAMPLES), "Monte-Carlo samples used to compute Us") 549 ("samples", po::value<int>()->default_value(DEFAULT_SAMPLES), "Monte-Carlo samples used to compute Us")
532 ("padding", po::value<unsigned int>()->default_value(DEFAULT_PADDING), "FFT padding for the objective bandpass") 550 ("padding", po::value<unsigned int>()->default_value(DEFAULT_PADDING), "FFT padding for the objective bandpass")
@@ -536,10 +554,6 @@ static void SetOptions(po::options_description &amp;desc) @@ -536,10 +554,6 @@ static void SetOptions(po::options_description &amp;desc)
536 ("recursive", "evaluate all Bessel functions recursively\n") 554 ("recursive", "evaluate all Bessel functions recursively\n")
537 ("recursive-us", "evaluate scattered-field Bessel functions recursively\n") 555 ("recursive-us", "evaluate scattered-field Bessel functions recursively\n")
538 ("lut-uf", "evaluate the focused-field using a look-up table\n") 556 ("lut-uf", "evaluate the focused-field using a look-up table\n")
539 -  
540 - ("output-type", po::value<string>()->default_value(DEFAULT_FIELD_TYPE), "output field value:\n magnitude, polarization, real, imaginary, angular-spectrum")  
541 - ("colormap", po::value<string>()->default_value(DEFAULT_COLORMAP), "colormap: gray, brewer")  
542 - ("append", "append result to an existing file\n (binary files only)")  
543 ; 557 ;
544 } 558 }
545 559
@@ -563,7 +577,6 @@ static void LoadParameters(int argc, char *argv[]) @@ -563,7 +577,6 @@ static void LoadParameters(int argc, char *argv[])
563 lWavelength(vm); 577 lWavelength(vm);
564 578
565 //load materials 579 //load materials
566 - //loadMaterials(vm);  
567 lMaterials(vm); 580 lMaterials(vm);
568 581
569 //load the sphere data 582 //load the sphere data
@@ -575,19 +588,11 @@ static void LoadParameters(int argc, char *argv[]) @@ -575,19 +588,11 @@ static void LoadParameters(int argc, char *argv[])
575 //load the position and orientation of the image plane 588 //load the position and orientation of the image plane
576 lImagePlane(vm); 589 lImagePlane(vm);
577 590
578 - //load spheres  
579 - //loadSpheres(vm);  
580 -  
581 -  
582 591
583 lNearfield(vm); 592 lNearfield(vm);
584 593
585 loadOutputParams(vm); 594 loadOutputParams(vm);
586 595
587 - //loadMicroscopeParams(vm);  
588 -  
589 - //loadSliceParams(vm);  
590 -  
591 //if an extended source will be used 596 //if an extended source will be used
592 if(vm["extended-source"].as<string>() != "") 597 if(vm["extended-source"].as<string>() != "")
593 { 598 {
planewave.h 0 โ†’ 100644
  1 +#ifndef __PLANEWAVE__
  2 +#define __PLANEWAVE__
  3 +
  4 +#include <iostream>
  5 +#include <sstream>
  6 +
  7 +#include "rts/math/vector.h"
  8 +
  9 +template<class T>
  10 +class planewave
  11 +{
  12 + rts::vector<T, 3> k;
  13 + rts::vector<T, 3> E;
  14 +
  15 + public:
  16 +
  17 + //constructor, initialize to an x-polarized wave propagating along z
  18 + planewave()
  19 + {
  20 + k = rts::vector<T, 3>(0, 0, 1);
  21 + E = rts::vector<T, 3>(1, 0, 0);
  22 + }
  23 +
  24 + planewave(rts::vector<T, 3> k_vec, rts::vector<T, 3> E0)
  25 + {
  26 + k = k_vec;
  27 +
  28 + //enforce k \dot E = 0
  29 + rts::vector<T, 3> s = E0.cross(k);
  30 + rts::vector<T, 3> E_hat;
  31 +
  32 + if(s.len() == 0)
  33 + E_hat = rts::vector<T, 3>(0, 0, 0);
  34 + else
  35 + E_hat = (s.cross(k)).norm();
  36 +
  37 + E = E_hat * (E_hat.dot(E0));
  38 + }
  39 +
  40 + // refract will bend the wave vector k to correspond to the normalized vector v
  41 + refract(rts::vector<T, 3> v)
  42 + {
  43 + //make sure that v is normalized
  44 + v = v.norm();
  45 +
  46 +
  47 + }
  48 +
  49 + std::string toStr()
  50 + {
  51 + std::stringstream ss;
  52 +
  53 + ss<<"k = "<<k<<std::endl;
  54 + ss<<"E = "<<E<<std::endl;
  55 +
  56 + return ss.str();
  57 + }
  58 +
  59 +
  60 +
  61 +};
  62 +
  63 +template <typename T>
  64 +std::ostream& operator<<(std::ostream& os, planewave<T> p)
  65 +{
  66 + os<<p.toStr();
  67 + return os;
  68 +}
  69 +
  70 +
  71 +
  72 +#endif
qtMainDialog.cpp 0 โ†’ 100755
  1 +#include "qtMainDialog.h"
  2 +#include <iostream>
  3 +
  4 +qtMainDialog::qtMainDialog(QWidget *parent, Qt::WindowFlags flags)
  5 + : QMainWindow(parent, flags)
  6 +{
  7 + ui.setupUi(this);
  8 +
  9 + outfile = "ui-out.bmp";
  10 +}
  11 +
  12 +qtMainDialog::~qtMainDialog()
  13 +{
  14 + updating = false;
  15 +}
  16 +
  17 +void qtMainDialog::closeEvent(QCloseEvent *event)
  18 +{
  19 + std::cout<<"Exiting"<<std::endl;
  20 + exit(0);
  21 +
  22 +}
  23 +
  24 +
  25 +/********************************
  26 +Populate the user interface
  27 +********************************/
  28 +void qtMainDialog::populateUi()
  29 +{
  30 + updating = true;
  31 +
  32 + //get the normal for the image plane
  33 + bsVector n = SCOPE->nf.pos.n();
  34 +
  35 + ui.spinNx->setValue(n[0]);
  36 + ui.spinNy->setValue(n[1]);
  37 + ui.spinNz->setValue(n[2]);
  38 +
  39 + //get the center point for the image plane
  40 + bsPoint c = SCOPE->nf.pos.p(0.5, 0.5);
  41 +
  42 + ui.spinCx->setValue(c[0]);
  43 + ui.spinCy->setValue(c[1]);
  44 + ui.spinCz->setValue(c[2]);
  45 +
  46 + //get the plane size (in microns)
  47 + ptype S = SCOPE->nf.pos.X.len();
  48 + ptype pad_div = SCOPE->padding * 2 + 1;
  49 + S /= pad_div;
  50 +
  51 + ui.spinS->setValue(S);
  52 +
  53 + //get the detector resolution
  54 + ui.spinR->setValue(SCOPE->Ud.R[0]);
  55 +
  56 + updating = false;
  57 +}
  58 +
  59 +/*********************************
  60 +Change the image plane position
  61 +*********************************/
  62 +void qtMainDialog::positionImage()
  63 +{
  64 + if(updating) return;
  65 +
  66 + //get the plane normal
  67 + bsVector n(ui.spinNx->value(), ui.spinNy->value(), ui.spinNz->value());
  68 +
  69 + //get the plane center point
  70 + bsPoint c(ui.spinCx->value(), ui.spinCy->value(), ui.spinCz->value());
  71 +
  72 + //get the plane orientation
  73 + ptype theta = ui.spinTheta->value();
  74 +
  75 + //get the plane size
  76 + ptype S = ui.spinS->value() * (2 * SCOPE->padding + 1);
  77 +
  78 + //create a new image plane
  79 + SCOPE->nf.pos = rts::quad<ptype, 3>(c, n, S, S, theta);
  80 +
  81 +
  82 +}
  83 +
  84 +/*********************************
  85 +Render an image
  86 +*********************************/
  87 +void qtMainDialog::renderImage()
  88 +{
  89 + //run the near-field simulation
  90 + SCOPE->SimulateScattering();
  91 +
  92 +
  93 + //determine the colormap type
  94 + rts::colormapType cmap = rts::cmGrayscale;
  95 +
  96 + if( ui.cmbColormap->currentText() == tr("brewer") )
  97 + cmap = rts::cmBrewer;
  98 + else
  99 + cmap = rts::cmGrayscale;
  100 +
  101 + //near field rendering
  102 + if( ui.radDisplayNearfield->isChecked() )
  103 + {
  104 + scalarslice S;
  105 + bool positive_vals = false;
  106 +
  107 + if( ui.cmbDisplayNearfield->currentText() == tr("magnitude") )
  108 + {
  109 + std::cout<<"magnitude"<<std::endl;
  110 + S = SCOPE->nf.U.Mag();
  111 + positive_vals = true;
  112 + //S.toImage(outfile.toStdString(), positive_vals, cmap);
  113 + }
  114 + else if( ui.cmbDisplayNearfield->currentText() == tr("real") )
  115 + S = SCOPE->nf.U.Real();
  116 + else if( ui.cmbDisplayNearfield->currentText() == tr("imaginary") )
  117 + S = SCOPE->nf.U.Imag();
  118 +
  119 + S.toImage(outfile.toStdString(), positive_vals, cmap);
  120 + }
  121 +
  122 + //run the far-field simulation
  123 + SCOPE->SimulateImaging();
  124 +
  125 + //far field rendering
  126 + if( ui.radDisplayFarfield->isChecked() )
  127 + {
  128 + scalarslice S;
  129 + bool positive_vals = false;
  130 +
  131 + if( ui.cmbDisplayFarfield->currentText() == tr("magnitude") )
  132 + {
  133 + S = SCOPE->Ud.Mag();
  134 + positive_vals = true;
  135 + }
  136 + else if( ui.cmbDisplayFarfield->currentText() == tr("real") )
  137 + S = SCOPE->Ud.Real();
  138 + else if( ui.cmbDisplayFarfield->currentText() == tr("imaginary") )
  139 + S = SCOPE->Ud.Imag();
  140 +
  141 + S.toImage(outfile.toStdString(), positive_vals, cmap);
  142 + }
  143 +
  144 + //detector rendering
  145 + if( ui.radDisplayDetector->isChecked() )
  146 + {
  147 + scalarslice I;
  148 + bool positive_vals = true;
  149 +
  150 + if( ui.cmbDisplayDetector->currentText() == tr("intensity") )
  151 + I = SCOPE->getIntensity();
  152 +
  153 + else if( ui.cmbDisplayDetector->currentText() == tr("absorbance") )
  154 + I = SCOPE->getAbsorbance();
  155 +
  156 + I.toImage(outfile.toStdString(), positive_vals, cmap);
  157 + }
  158 +
  159 +}
qtMainDialog.h 0 โ†’ 100755
  1 +#ifndef VolumeSpiderDialog_H
  2 +#define VolumeSpiderDialog_H
  3 +
  4 +#include <QtWidgets/QMainWindow>
  5 +#include <QDragEnterEvent>
  6 +#include <qmimedata.h>
  7 +#include <qfiledialog.h>
  8 +#include <qinputdialog.h>
  9 +#include "ui_qtMainDialog.h"
  10 +
  11 +//simulation parameters
  12 +#include "microscope.h"
  13 +extern microscopeStruct* SCOPE;
  14 +
  15 +//#include "fileout.h"
  16 +//extern fileoutStruct gFileOut;
  17 +#include "rts/visualization/colormap.h"
  18 +
  19 +
  20 +class qtMainDialog : public QMainWindow
  21 +{
  22 +Q_OBJECT
  23 +
  24 +public:
  25 +qtMainDialog(QWidget *parent = 0, Qt::WindowFlags flags = 0);
  26 +~qtMainDialog();
  27 +
  28 +void closeEvent(QCloseEvent *event);
  29 +bool updating;
  30 +
  31 +QString outfile;
  32 +
  33 +void refreshUI()
  34 +{
  35 + updating = false;
  36 +}
  37 +
  38 +void populateUi();
  39 +void renderImage();
  40 +void positionImage();
  41 +
  42 +private:
  43 +Ui::qtMainDialogUI ui;
  44 +
  45 +public slots:
  46 +
  47 +/*************
  48 +Buttons
  49 +*************/
  50 +void on_btnRender_pressed()
  51 +{
  52 + renderImage();
  53 +}
  54 +void on_btnClearDetector_pressed()
  55 +{
  56 + SCOPE->clearDetector();
  57 +}
  58 +
  59 +/***************
  60 +Spinners
  61 +***************/
  62 +void on_spinCx_valueChanged(double d)
  63 +{
  64 + positionImage();
  65 +}
  66 +void on_spinCy_valueChanged(double d)
  67 +{
  68 + positionImage();
  69 +}
  70 +void on_spinCz_valueChanged(double d)
  71 +{
  72 + positionImage();
  73 +}
  74 +void on_spinNx_valueChanged(double d)
  75 +{
  76 + positionImage();
  77 +}
  78 +void on_spinNy_valueChanged(double d)
  79 +{
  80 + positionImage();
  81 +}
  82 +void on_spinNz_valueChanged(double d)
  83 +{
  84 + positionImage();
  85 +}
  86 +void on_spinS_valueChanged(double d)
  87 +{
  88 + positionImage();
  89 +}
  90 +void on_spinTheta_valueChanged(double d)
  91 +{
  92 + positionImage();
  93 +}
  94 +
  95 +
  96 +};
  97 +
  98 +#endif // INTERACTIVEMIE_H
qtMainDialog.ui 0 โ†’ 100644
  1 +<?xml version="1.0" encoding="UTF-8"?>
  2 +<ui version="4.0">
  3 + <class>qtMainDialogUI</class>
  4 + <widget class="QMainWindow" name="qtMainDialogUI">
  5 + <property name="geometry">
  6 + <rect>
  7 + <x>0</x>
  8 + <y>0</y>
  9 + <width>800</width>
  10 + <height>600</height>
  11 + </rect>
  12 + </property>
  13 + <property name="windowTitle">
  14 + <string>MainWindow</string>
  15 + </property>
  16 + <widget class="QWidget" name="centralwidget">
  17 + <widget class="QDoubleSpinBox" name="spinCx">
  18 + <property name="geometry">
  19 + <rect>
  20 + <x>80</x>
  21 + <y>40</y>
  22 + <width>81</width>
  23 + <height>27</height>
  24 + </rect>
  25 + </property>
  26 + <property name="minimum">
  27 + <double>-99.989999999999995</double>
  28 + </property>
  29 + <property name="singleStep">
  30 + <double>0.010000000000000</double>
  31 + </property>
  32 + </widget>
  33 + <widget class="QDoubleSpinBox" name="spinCy">
  34 + <property name="geometry">
  35 + <rect>
  36 + <x>170</x>
  37 + <y>40</y>
  38 + <width>81</width>
  39 + <height>27</height>
  40 + </rect>
  41 + </property>
  42 + <property name="minimum">
  43 + <double>-99.989999999999995</double>
  44 + </property>
  45 + <property name="singleStep">
  46 + <double>0.010000000000000</double>
  47 + </property>
  48 + </widget>
  49 + <widget class="QDoubleSpinBox" name="spinCz">
  50 + <property name="geometry">
  51 + <rect>
  52 + <x>260</x>
  53 + <y>40</y>
  54 + <width>81</width>
  55 + <height>27</height>
  56 + </rect>
  57 + </property>
  58 + <property name="minimum">
  59 + <double>-99.989999999999995</double>
  60 + </property>
  61 + <property name="singleStep">
  62 + <double>0.010000000000000</double>
  63 + </property>
  64 + </widget>
  65 + <widget class="QDoubleSpinBox" name="spinNy">
  66 + <property name="geometry">
  67 + <rect>
  68 + <x>170</x>
  69 + <y>70</y>
  70 + <width>81</width>
  71 + <height>27</height>
  72 + </rect>
  73 + </property>
  74 + <property name="minimum">
  75 + <double>-99.989999999999995</double>
  76 + </property>
  77 + <property name="singleStep">
  78 + <double>0.010000000000000</double>
  79 + </property>
  80 + </widget>
  81 + <widget class="QDoubleSpinBox" name="spinNz">
  82 + <property name="geometry">
  83 + <rect>
  84 + <x>260</x>
  85 + <y>70</y>
  86 + <width>81</width>
  87 + <height>27</height>
  88 + </rect>
  89 + </property>
  90 + <property name="minimum">
  91 + <double>-99.989999999999995</double>
  92 + </property>
  93 + <property name="singleStep">
  94 + <double>0.010000000000000</double>
  95 + </property>
  96 + </widget>
  97 + <widget class="QDoubleSpinBox" name="spinNx">
  98 + <property name="geometry">
  99 + <rect>
  100 + <x>80</x>
  101 + <y>70</y>
  102 + <width>81</width>
  103 + <height>27</height>
  104 + </rect>
  105 + </property>
  106 + <property name="minimum">
  107 + <double>-99.989999999999995</double>
  108 + </property>
  109 + <property name="singleStep">
  110 + <double>0.010000000000000</double>
  111 + </property>
  112 + </widget>
  113 + <widget class="QLabel" name="label">
  114 + <property name="geometry">
  115 + <rect>
  116 + <x>10</x>
  117 + <y>40</y>
  118 + <width>66</width>
  119 + <height>17</height>
  120 + </rect>
  121 + </property>
  122 + <property name="text">
  123 + <string>pc = </string>
  124 + </property>
  125 + <property name="alignment">
  126 + <set>Qt::AlignRight|Qt::AlignTrailing|Qt::AlignVCenter</set>
  127 + </property>
  128 + </widget>
  129 + <widget class="QLabel" name="label_2">
  130 + <property name="geometry">
  131 + <rect>
  132 + <x>10</x>
  133 + <y>70</y>
  134 + <width>66</width>
  135 + <height>17</height>
  136 + </rect>
  137 + </property>
  138 + <property name="text">
  139 + <string>n = </string>
  140 + </property>
  141 + <property name="alignment">
  142 + <set>Qt::AlignRight|Qt::AlignTrailing|Qt::AlignVCenter</set>
  143 + </property>
  144 + </widget>
  145 + <widget class="QDoubleSpinBox" name="spinS">
  146 + <property name="geometry">
  147 + <rect>
  148 + <x>80</x>
  149 + <y>130</y>
  150 + <width>81</width>
  151 + <height>27</height>
  152 + </rect>
  153 + </property>
  154 + <property name="singleStep">
  155 + <double>0.010000000000000</double>
  156 + </property>
  157 + </widget>
  158 + <widget class="QLabel" name="label_3">
  159 + <property name="geometry">
  160 + <rect>
  161 + <x>10</x>
  162 + <y>130</y>
  163 + <width>66</width>
  164 + <height>17</height>
  165 + </rect>
  166 + </property>
  167 + <property name="text">
  168 + <string>S = </string>
  169 + </property>
  170 + <property name="alignment">
  171 + <set>Qt::AlignRight|Qt::AlignTrailing|Qt::AlignVCenter</set>
  172 + </property>
  173 + </widget>
  174 + <widget class="QGroupBox" name="groupBox">
  175 + <property name="geometry">
  176 + <rect>
  177 + <x>70</x>
  178 + <y>210</y>
  179 + <width>331</width>
  180 + <height>161</height>
  181 + </rect>
  182 + </property>
  183 + <property name="title">
  184 + <string>Display</string>
  185 + </property>
  186 + <widget class="QRadioButton" name="radDisplayNearfield">
  187 + <property name="geometry">
  188 + <rect>
  189 + <x>113</x>
  190 + <y>30</y>
  191 + <width>116</width>
  192 + <height>22</height>
  193 + </rect>
  194 + </property>
  195 + <property name="text">
  196 + <string>near field</string>
  197 + </property>
  198 + <property name="checked">
  199 + <bool>true</bool>
  200 + </property>
  201 + </widget>
  202 + <widget class="QRadioButton" name="radDisplayFarfield">
  203 + <property name="geometry">
  204 + <rect>
  205 + <x>113</x>
  206 + <y>60</y>
  207 + <width>116</width>
  208 + <height>22</height>
  209 + </rect>
  210 + </property>
  211 + <property name="text">
  212 + <string>far field</string>
  213 + </property>
  214 + </widget>
  215 + <widget class="QComboBox" name="cmbDisplayNearfield">
  216 + <property name="geometry">
  217 + <rect>
  218 + <x>10</x>
  219 + <y>30</y>
  220 + <width>101</width>
  221 + <height>27</height>
  222 + </rect>
  223 + </property>
  224 + <item>
  225 + <property name="text">
  226 + <string>magnitude</string>
  227 + </property>
  228 + </item>
  229 + <item>
  230 + <property name="text">
  231 + <string>real</string>
  232 + </property>
  233 + </item>
  234 + <item>
  235 + <property name="text">
  236 + <string>imaginary</string>
  237 + </property>
  238 + </item>
  239 + </widget>
  240 + <widget class="QComboBox" name="cmbDisplayDetector">
  241 + <property name="geometry">
  242 + <rect>
  243 + <x>10</x>
  244 + <y>90</y>
  245 + <width>101</width>
  246 + <height>27</height>
  247 + </rect>
  248 + </property>
  249 + <item>
  250 + <property name="text">
  251 + <string>intensity</string>
  252 + </property>
  253 + </item>
  254 + <item>
  255 + <property name="text">
  256 + <string>absorbance</string>
  257 + </property>
  258 + </item>
  259 + </widget>
  260 + <widget class="QRadioButton" name="radDisplayDetector">
  261 + <property name="geometry">
  262 + <rect>
  263 + <x>113</x>
  264 + <y>90</y>
  265 + <width>116</width>
  266 + <height>22</height>
  267 + </rect>
  268 + </property>
  269 + <property name="text">
  270 + <string>detector</string>
  271 + </property>
  272 + </widget>
  273 + <widget class="QComboBox" name="cmbDisplayFarfield">
  274 + <property name="geometry">
  275 + <rect>
  276 + <x>10</x>
  277 + <y>60</y>
  278 + <width>101</width>
  279 + <height>27</height>
  280 + </rect>
  281 + </property>
  282 + <item>
  283 + <property name="text">
  284 + <string>magnitude</string>
  285 + </property>
  286 + </item>
  287 + <item>
  288 + <property name="text">
  289 + <string>real</string>
  290 + </property>
  291 + </item>
  292 + <item>
  293 + <property name="text">
  294 + <string>imaginary</string>
  295 + </property>
  296 + </item>
  297 + </widget>
  298 + <widget class="QComboBox" name="cmbColormap">
  299 + <property name="geometry">
  300 + <rect>
  301 + <x>220</x>
  302 + <y>30</y>
  303 + <width>91</width>
  304 + <height>27</height>
  305 + </rect>
  306 + </property>
  307 + <item>
  308 + <property name="text">
  309 + <string>brewer</string>
  310 + </property>
  311 + </item>
  312 + <item>
  313 + <property name="text">
  314 + <string>grayscale</string>
  315 + </property>
  316 + </item>
  317 + </widget>
  318 + <widget class="QPushButton" name="btnRender">
  319 + <property name="geometry">
  320 + <rect>
  321 + <x>10</x>
  322 + <y>130</y>
  323 + <width>98</width>
  324 + <height>27</height>
  325 + </rect>
  326 + </property>
  327 + <property name="text">
  328 + <string>Render</string>
  329 + </property>
  330 + </widget>
  331 + <widget class="QPushButton" name="btnClearDetector">
  332 + <property name="geometry">
  333 + <rect>
  334 + <x>220</x>
  335 + <y>90</y>
  336 + <width>51</width>
  337 + <height>27</height>
  338 + </rect>
  339 + </property>
  340 + <property name="text">
  341 + <string>Clear</string>
  342 + </property>
  343 + </widget>
  344 + </widget>
  345 + <widget class="QSpinBox" name="spinR">
  346 + <property name="geometry">
  347 + <rect>
  348 + <x>260</x>
  349 + <y>130</y>
  350 + <width>81</width>
  351 + <height>27</height>
  352 + </rect>
  353 + </property>
  354 + <property name="maximum">
  355 + <number>999999</number>
  356 + </property>
  357 + <property name="value">
  358 + <number>0</number>
  359 + </property>
  360 + </widget>
  361 + <widget class="QLabel" name="label_4">
  362 + <property name="geometry">
  363 + <rect>
  364 + <x>190</x>
  365 + <y>130</y>
  366 + <width>66</width>
  367 + <height>17</height>
  368 + </rect>
  369 + </property>
  370 + <property name="text">
  371 + <string>R = </string>
  372 + </property>
  373 + <property name="alignment">
  374 + <set>Qt::AlignRight|Qt::AlignTrailing|Qt::AlignVCenter</set>
  375 + </property>
  376 + </widget>
  377 + <widget class="QDoubleSpinBox" name="spinTheta">
  378 + <property name="geometry">
  379 + <rect>
  380 + <x>80</x>
  381 + <y>100</y>
  382 + <width>81</width>
  383 + <height>27</height>
  384 + </rect>
  385 + </property>
  386 + <property name="maximum">
  387 + <double>6.280000000000000</double>
  388 + </property>
  389 + <property name="singleStep">
  390 + <double>0.010000000000000</double>
  391 + </property>
  392 + </widget>
  393 + <widget class="QLabel" name="label_5">
  394 + <property name="geometry">
  395 + <rect>
  396 + <x>10</x>
  397 + <y>100</y>
  398 + <width>66</width>
  399 + <height>17</height>
  400 + </rect>
  401 + </property>
  402 + <property name="text">
  403 + <string>theta = </string>
  404 + </property>
  405 + <property name="alignment">
  406 + <set>Qt::AlignRight|Qt::AlignTrailing|Qt::AlignVCenter</set>
  407 + </property>
  408 + </widget>
  409 + </widget>
  410 + <widget class="QMenuBar" name="menubar">
  411 + <property name="geometry">
  412 + <rect>
  413 + <x>0</x>
  414 + <y>0</y>
  415 + <width>800</width>
  416 + <height>25</height>
  417 + </rect>
  418 + </property>
  419 + </widget>
  420 + <widget class="QStatusBar" name="statusbar"/>
  421 + </widget>
  422 + <resources/>
  423 + <connections/>
  424 +</ui>
@@ -11,20 +11,50 @@ scalarslice::scalarslice(int x, int y) @@ -11,20 +11,50 @@ scalarslice::scalarslice(int x, int y)
11 R[1] = y; 11 R[1] = y;
12 12
13 //allocate memory on the GPU 13 //allocate memory on the GPU
14 - HANDLE_ERROR(cudaMalloc( (void**)&S, sizeof(ptype) * x * y )); 14 + HANDLE_ERROR(cudaMalloc( (void**)&S, sizeof(ptype) * x * y ));
  15 + //std::cout<<"Scalerslice created."<<std::endl;
15 } 16 }
16 17
17 scalarslice::scalarslice() 18 scalarslice::scalarslice()
18 { 19 {
19 R[0] = R[1] = 0; 20 R[0] = R[1] = 0;
20 S = NULL; 21 S = NULL;
  22 +
  23 + //std::cout<<"Scalerslice created (default)."<<std::endl;
21 } 24 }
22 25
23 scalarslice::~scalarslice() 26 scalarslice::~scalarslice()
24 { 27 {
25 if(S != NULL) 28 if(S != NULL)
26 HANDLE_ERROR(cudaFree(S)); 29 HANDLE_ERROR(cudaFree(S));
27 - S = NULL; 30 + S = NULL;
  31 + R[0] = R[1] = 0;
  32 +
  33 + //std::cout<<"Scalerslice destroyed."<<std::endl;
  34 +}
  35 +
  36 +//assignment operator
  37 +scalarslice & scalarslice::operator= (const scalarslice & rhs)
  38 +{
  39 + //de-allocate any existing GPU memory
  40 + if(S != NULL)
  41 + HANDLE_ERROR(cudaFree(S));
  42 +
  43 + //copy the slice resolution
  44 + R[0] = rhs.R[0];
  45 + R[1] = rhs.R[1];
  46 +
  47 + //allocate the necessary memory
  48 + HANDLE_ERROR(cudaMalloc(&S, sizeof(ptype) * R[0] * R[1]));
  49 +
  50 + //copy the slice
  51 + HANDLE_ERROR(cudaMemcpy(S, rhs.S, sizeof(ptype) * R[0] * R[1], cudaMemcpyDeviceToDevice));
  52 +
  53 +
  54 + std::cout<<"Assignment operator."<<std::endl;
  55 +
  56 + return *this;
  57 +
28 } 58 }
29 59
30 void scalarslice::toImage(std::string filename, ptype vmin, ptype vmax, rts::colormapType cmap) 60 void scalarslice::toImage(std::string filename, ptype vmin, ptype vmax, rts::colormapType cmap)
@@ -2,7 +2,7 @@ @@ -2,7 +2,7 @@
2 #define RTS_SCALAR_SLICE 2 #define RTS_SCALAR_SLICE
3 3
4 #include "dataTypes.h" 4 #include "dataTypes.h"
5 -#include "rts/graphics/colormap.h" 5 +#include "rts/visualization/colormap.h"
6 6
7 struct scalarslice 7 struct scalarslice
8 { 8 {
@@ -19,7 +19,10 @@ struct scalarslice @@ -19,7 +19,10 @@ struct scalarslice
19 19
20 void toImage(std::string filename, ptype vmin, ptype vmax, rts::colormapType cmap = rts::cmBrewer); 20 void toImage(std::string filename, ptype vmin, ptype vmax, rts::colormapType cmap = rts::cmBrewer);
21 void toImage(std::string filename, bool positive = true, rts::colormapType cmap = rts::cmBrewer); 21 void toImage(std::string filename, bool positive = true, rts::colormapType cmap = rts::cmBrewer);
22 - void toEnvi(std::string filename, ptype wavelength = 0, bool append = false); 22 + void toEnvi(std::string filename, ptype wavelength = 0, bool append = false);
  23 +
  24 + //assignment operator
  25 + scalarslice & operator= (const scalarslice & rhs);
23 26
24 }; 27 };
25 28
@@ -6,7 +6,7 @@ @@ -6,7 +6,7 @@
6 #include <stdlib.h> 6 #include <stdlib.h>
7 #include <fstream> 7 #include <fstream>
8 8
9 -using namespace rts; 9 +//using namespace rts;
10 using namespace std; 10 using namespace std;
11 11
12 int cbessjyva(double v,complex<double> z,double &vm,complex<double>*cjv, 12 int cbessjyva(double v,complex<double> z,double &vm,complex<double>*cjv,
@@ -18,7 +18,7 @@ int cbessjyva_sph(int v,complex&lt;double&gt; z,double &amp;vm,complex&lt;double&gt;*cjv, @@ -18,7 +18,7 @@ int cbessjyva_sph(int v,complex&lt;double&gt; z,double &amp;vm,complex&lt;double&gt;*cjv,
18 int bessjyv_sph(int v, double z, double &vm, double* cjv, 18 int bessjyv_sph(int v, double z, double &vm, double* cjv,
19 double* cyv, double* cjvp, double* cyvp); 19 double* cyv, double* cjvp, double* cyvp);
20 20
21 -void sphere::calcCoeff(ptype lambda, rtsComplex<ptype> ri) 21 +void sphere::calcCoeff(ptype lambda, bsComplex ri)
22 { 22 {
23 /* These calculations are done at high-precision on the CPU 23 /* These calculations are done at high-precision on the CPU
24 since they are only required once for each wavelength. 24 since they are only required once for each wavelength.
@@ -219,7 +219,7 @@ void sphere::calcLut(bsComplex* j, bsComplex* h, ptype lambda, bsComplex n, int @@ -219,7 +219,7 @@ void sphere::calcLut(bsComplex* j, bsComplex* h, ptype lambda, bsComplex n, int
219 calcHankelLut(h, k, rR); 219 calcHankelLut(h, k, rR);
220 } 220 }
221 221
222 -void sphere::calcUp(ptype lambda, bsComplex n, rts::rtsQuad<ptype, 3> nfPlane, unsigned int R) 222 +void sphere::calcUp(ptype lambda, bsComplex n, rts::quad<ptype, 3> nfPlane, unsigned int R)
223 { 223 {
224 //calculate the parameters of the lookup table 224 //calculate the parameters of the lookup table
225 225
@@ -52,7 +52,7 @@ struct sphere @@ -52,7 +52,7 @@ struct sphere
52 52
53 //assignment operator 53 //assignment operator
54 sphere & operator=(const sphere &rhs); 54 sphere & operator=(const sphere &rhs);
55 - 55 +
56 //copy constructor 56 //copy constructor
57 sphere(const sphere &rhs); 57 sphere(const sphere &rhs);
58 58
@@ -81,7 +81,7 @@ struct sphere @@ -81,7 +81,7 @@ struct sphere
81 void calcHankelLut(bsComplex* h, ptype k, int rR); 81 void calcHankelLut(bsComplex* h, ptype k, int rR);
82 82
83 //calculate the scattering domain Us(theta, r) 83 //calculate the scattering domain Us(theta, r)
84 - void calcUp(ptype lambda, bsComplex n, rts::rtsQuad<ptype, 3> nfPlane, unsigned int R); 84 + void calcUp(ptype lambda, bsComplex n, rts::quad<ptype, 3> nfPlane, unsigned int R);
85 85
86 void scalarUsp(bsComplex* h, int rR, int thetaR); 86 void scalarUsp(bsComplex* h, int rR, int thetaR);
87 void scalarUip(bsComplex* j, int aR, int thetaR); 87 void scalarUip(bsComplex* j, int aR, int thetaR);