#ifndef RTS_BEAM #define RTS_BEAM #include "../math/vector.h" #include "../math/function.h" #include "../optics/planewave.h" #include namespace stim{ template class beam : public planewave

{ public: enum beam_type {Uniform, Bartlett, Hamming, Hanning}; private: P _na[2]; //numerical aperature of the focusing optics vec

f; //focal point function apod; //apodization function unsigned int apod_res; //resolution of apodization filter functions void apod_uniform() { apod = (P)1; } void apod_bartlett() { apod = (P)1; apod.insert((P)1, (P)0); } void apod_hanning() { apod = (P)0; P x, y; for(unsigned int n=0; n k = rts::vec

(0, 0, rtsTAU), vec

_E0 = rts::vec

(1, 0, 0), beam_type _apod = Uniform) : planewave

(k, _E0) { _na[0] = (P)0.0; _na[1] = (P)1.0; f = vec

( (P)0, (P)0, (P)0 ); apod_res = 256; //set the default resolution for apodization filters set_apod(_apod); //set the apodization function type } beam

refract(rts::vec

kn) const{ beam

new_beam; new_beam._na[0] = _na[0]; new_beam._na[1] = _na[1]; rts::planewave

pw = planewave

::bend(kn); //std::cout< > mc(unsigned int N = 100000, unsigned int seed = 0) const { /*Create Monte-Carlo samples of a cassegrain objective by performing uniform sampling of a sphere and projecting these samples onto an inscribed sphere. seed = seed for the random number generator */ srand(seed); //seed the random number generator vec

k_hat = beam::k.norm(); ///compute the rotation operator to transform (0, 0, 1) to k P cos_angle = k_hat.dot(rts::vec

(0, 0, 1)); rts::matrix rotation; //if the cosine of the angle is -1, the rotation is just a flip across the z axis if(cos_angle == -1){ rotation(2, 2) = -1; } else if(cos_angle != 1.0) { rts::vec

r_axis = rts::vec

(0, 0, 1).cross(k_hat).norm(); //compute the axis of rotation P angle = acos(cos_angle); //compute the angle of rotation rts::quaternion

quat; //create a quaternion describing the rotation quat.CreateRotation(angle, r_axis); rotation = quat.toMatrix3(); //compute the rotation matrix } //find the phi values associated with the cassegrain ring P PHI[2]; PHI[0] = (P)asin(_na[0]); PHI[1] = (P)asin(_na[1]); //calculate the z-axis cylinder coordinates associated with these angles P Z[2]; Z[0] = cos(PHI[0]); Z[1] = cos(PHI[1]); P range = Z[0] - Z[1]; std::vector< planewave

> samples; //create a vector of plane waves //draw a distribution of random phi, z values P z, phi, theta; for(int i=0; i spherical(1, theta, phi); //convert from spherical to cartesian coordinates rts::vec

cart = spherical.sph2cart(); vec

k_prime = rotation * cart; //create a sample vector //store a wave refracted along the given direction //std::cout<<"k prime: "<::refract(k_prime) * apod(phi/PHI[1])); } return samples; } std::string str() { std::stringstream ss; ss<<"Beam:"<