1 / 73

Monte Carlo, Euler and Quaternions

Monte Carlo, Euler and Quaternions. Jessica C. Ramella-Roman. Stokes vector reference. Yesterday, Meridian plane Three steps “Rotate” to scattering plane Rotational matrix Scatter Scattering matrix “ Rotate” to a new plane Rotational matrix. Today Monte Carlo.

booker
Download Presentation

Monte Carlo, Euler and Quaternions

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript


  1. Monte Carlo, Euler and Quaternions Jessica C. Ramella-Roman

  2. Stokes vector reference • Yesterday, Meridian plane • Three steps • “Rotate” to scattering plane • Rotational matrix • Scatter • Scattering matrix • “Rotate” to a new plane • Rotational matrix

  3. Today Monte Carlo • Reference frame is a triplet of unit vectors • Rotation are about an axis and follow a specific order http://www.grc.nasa.gov/WWW/K-12/airplane

  4. Today Monte Carlo • Reference frame is a triplet of unit vectors • Rotation of the frame is done in two steps.

  5. Euler and Quaternion • The only difference between today programs is the way we handle these rotations • Euler angles rotation Quaternionsalgebra

  6. Euler Monte Carlo • The photon polarization reference frame is tracked at any time via a triplet of unit vectors that are rotated by an azimuth and scattering angle according to a predefined order • Introduced by Bartelet al. (*) * S. Bartel, A. Hielsher, “Monte Carlo simulations of the diffuse backscattering Mueller matrix for highly scattering media,” Applied Optics, Vol. 39, No. 10; (2000).

  7. Euler matrices • Euler's rotation theorem: any rotation may be described using three angles. • And three rotation matrices • About z axis • About x’ axis • About z’ axis Drawings, Wolfram, Mathworld

  8. Euler cnt. • We only need TWO vectors and TWO rotations • Two unit vectors vand u • The third unit vector is defined by the cross product of vand u and is calculated only when a photon reaches a boundary.

  9. Euler cnt. • Advantage • Stokes vector is rotated only once for each scattering event instead of twice as in the meridian plane method. • Simple to implement and visualize • Drawback – • Gimbal lock - makes the rotation fail for angles exactly equal to 90˚ (rare event).

  10. Monte Carlo flow chart

  11. Launch • The two unit vector v and u • v and u define the starting Stokes vector reference plane. • The unit vectorurepresents the direction of photon propagation.

  12. Move • The photon is moved a distance Ds, the new coordinates of the photon are

  13. Drop • Reduction of photon weight (W) • According to material albedo • albedo = ms/(ms+ ma) • absorbed = W*(1-albedo)

  14. Scatter • The rejection method establishes the scattering angle q and azimuth anglef. • The Stokes vector must be rotated by an angle f into the scattering plane before scattering can occur.

  15. Referencing to scattering plane • Done as in meridian plane Monte Carlo • Multiply the Stokes vector by matrix R

  16. Scatter cnt • The interaction with a spherical particle is achieved by multiplying the Stokes vector with a scattering matrix M(q)

  17. M(q) elements • s11, s12 , s33 , s34 calculated with Mie theory • Expressed as • S1, S2* – Mie scattering • *http://omlc.ogi.edu/calc/mie_calc.html • C. Bohren and D. R. Huffman, • Absorption and scattering of light by small • particles, (Wiley Science Paperback Series,1998).

  18. Reference frame • After scattering the Stokes vector the reference coordinate system vu must be updated. • The two rotations, for the angles f and q can be obtained using Euler’s rotational matrices.

  19. Rotational matrix ROT • Rodrigues’s rotational matrix ROT • Accomplishes the general case of rotating any vector by an angle y about an axis K • Where K is the rotational axis, • c=cos(y), s=sin(y) and v=1-cos(y).

  20. Vector rotation • First the unit vector v is rotated about the vector u by an angle f • Multiply v by the rotational matrix ROT(u,f); u remains unchanged

  21. Vector rotation • Second u is rotated about the newly generated v by an angle q. • This is done multiplying the unit vector u by the rotational matrix ROT(v,q).

  22. Update of direction cosines • This is done as in standard Monte Carlo • q and f If |uz|≈ 1

  23. Update of direction cosines If |uz|≠ 1

  24. Photon life • The life of a photon ends when the photon passes through a boundary or when its weight W value falls below a threshold. • Roulette is used to terminate the photon packet when W £Wth. • Gives the photon packet one chance of surviving • If the photon packet does not survive the roulette, the photon weight is reduced to zero and the photon is terminated.

  25. Boundaries • Two final rotations of the Stokes vector are necessary to put the photon status of polarization in the detector reference frame. • This will be achieved with matrix multiplication of two rotational matrix by stokes vector

  26. Boundaries –R(e) • w is reconstructed as the cross product of v and u. • An angle e is needed to rotate the Stokes vector into a scattering plane

  27. Boundaries –R(e) • This is the ONLY reason we need the vector u,v, andw

  28. Boundaries –R(e) • This rotation is about the direction of propagation of the photon, i.e. the axis u • Resulting position of v

  29. Boundaries –R(y) • A second rotation of an angle y about the Z axis • Put the photon reference frame in the detector reference frame

  30. Boundaries –R(y) • Transmission • Reflection

  31. Boundaries End • Stokes vector in the detector frame of reference • We need vector u,v, andwto obtain these angles

  32. Quaternion Monte Carlo • Quaternion Monte Carlo simply uses Quaternion Algebra to handle vector rotation. • Advantage • Stokes vector is rotated only once for each scattering event instead of twice as in the meridian plane method. • No issue with Gimbal lock • Optimized for computer simulations

  33. Quaternions • A quaternion is a 4-tuple of real numbers; it is an element of R4. • Quaternion can also be defined as the sum of a scalar part q0 and a vector part Q in R3of the form

  34. Quaternionscnt. • The vector part Q is the rotational axis and the scalar part is the angle of rotation. • Multiplication of a vector t by the quaternion is equivalent to rotating the vector t around the vector Q of an angle q0.

  35. Vector rotation • First the unit vector v is rotated about the vector u by an angle f

  36. Rotation - f • The first rotation is about the vector u by an angle f. • This is done generating the quaternion qf • And then using the quaternion operator • q∗vqon the vector v • q∗ is the complex conjugate of q

  37. Vector rotation • Second u is rotated about the newly generated v by an angle q.

  38. Rotation - q • The second rotation is about the vector v by an angle q. • This is done generating the quaternion qq • And then using the quaternion operator • qq*uqqon the vectoru

  39. Rotations • These steps are repeated for every scattering event. • At the boundaries the last aligning rotations are the same as in Euler Monte Carlo.

  40. Testing • Comparison with Evans Code • In 1991 Evans designed an adding doubling code that can calculate both the radiance and flux of a polarized light beam exiting the atmosphere • plane parallel slab of thickness L = 4/µs, • absorption coefficient µa=0, • unpolarized incident beam • wavelength l = 0.632 nm.

  41. Evans - Reflectance mode Incident I=[1 0 0 0] Reflectance mode, comparison between Evans adding-doubling code and the meridian plane Monte Carlo program. The results do not include the final rotation for a single detector.

  42. Evans – Transmission mode Incident I=[1 0 0 0] Transmission mode, comparison between Evans adding doubling code and the meridian plane Monte Carlo program. The results are not corrected for a single detector

  43. How to run the code • Download the code • http://faculty.cua.edu/ramella/MonteCarlo/index.html

  44. How to run the code • Download the code • http://faculty.cua.edu/ramella/MonteCarlo/index.html

  45. Download the code Programs

  46. Make Command line Compilation/linking step

  47. Run Command line This program will run a full Mueller matrix Monte Carlo launching four vectors [1 1 0 0] H 0 degree polarized [1 -1 0 0] V 90 degrees polarized [1 0 1 0] P 45 degrees polarized [1 0 0 1] R –Right circular

  48. Memory usage • Intel (64-bit) – Macbook pro 2.3Ghz I7 8GB RAM • Real mem 860 KB • Virtual mem 17 MB

  49. Output is a set of Stokes vector images for each launched Stokes vector outHI -> Horizontal incident polarization and I portion of the Stokes vector outHQ -> Horizontal incident polarization and Q portion of the Stokes vector out

  50. Inside the code GNU Licence

More Related