Help with orientation issues

Summary:

Can't see a difference in the extinction values when I rotate an ellipsoidal shape

Help with orientation issues

Alison Deatsch (guest) 30 Jun 2016 20:03

I have been working well with this program for awhile, but I hit snags every time I try to look at anything more complicated than a single sphere. I have tried both ellipsoids and cylinders, but I am never able to get my program to register the long edge of my shape. We have a software that converts the resulting extinction values as a function of wavelength into a number density as a function of diameter. But even when I rotate the ellipsoid at all angles, I am never able to see the number density peak signifying the long edge.

The best example I can give to help troubleshoot my issue is this ellipsoid (ddscat.par file copied at the bottom). I intended to create an ellipsoid with dimensions 7x7x11 nm and approximately 10,000 dipoles. I want the Qext values over a wavelength range of 200-1000 nm. I have tried this file for theta=0, theta=90, theta=0-90 (with 10 angles), as well as 0-90 for beta and theta, all while holding all the other angles constant at 0.. and every time I get the exact same Qext values. Shouldn't I get different extinction values when I am looking at the ellipsoid head on versus when it is rotated 90 degrees so I am looking at the long edge? Do you notice anything I am missing or doing wrong or do you have any suggestions of what I might try to change??

I've tried making this larger (70x70x110) and I've tried with as many as 100,000 dipoles as well, never with any success at recognizing the longer edge. I am using version 7.3.1, but not the most recent version of it.

Thank you!!

' ========== Parameter file for v7.3 ==================='
' Preliminaries '
'NOTORQ' = CMDTRQ*6 (DOTORQ, NOTORQ) — either do or skip torque calculations
'PBCGS2' = CMDSOL*6 (PBCGS2, PBCGST, GPBICG, QMRCCG, PETRKP) — CCG method
'GPFAFT' = CMETHD*6 (GPFAFT, FFTMKL) — FFT method
'GKDLDR' = CALPHA*6 (GKDLDR, LATTDR, FLTRCD) — DDA method
'NOTBIN' = CBINFLAG (NOTBIN, ORIBIN, ALLBIN) — binary output?
' Initial Memory Allocation '
100 100 100 = dimensioning allowance for target generation
' Target Geometry and Composition '
'ELLIPSOID' = CSHAPE*9 shape directive
25 25 39.286 = shape parameters 1 - 3
1 = NCOMP = number of dielectric materials
'../diel/polystyeq' = file with refractive index 1
' Additional Nearfield calculation? '
0 = NRFLD (=0 to skip nearfield calc., =1 to calculate nearfield E)
0.0 0.0 0.0 0.0 0.0 0.0 (fract. extens. of calc. vol. in -x,+x,-y,+y,-z,+z)
' Error Tolerance '
1.00e-5 = TOL = MAX ALLOWED (NORM OF |G>=AC|E>-ACA|X>)/(NORM OF AC|E>)
' Maximum number of iterations '
1000 = MXITER
' Integration limiter for PBC calculations '
1.00e-2 = GAMMA (1e-2 is normal, 3e-3 for greater accuracy)
' Angular resolution for calculation of <cos>, etc. '
0.5 = ETASCA (number of angles is proportional to [(3+x)/ETASCA]^2 )
' Wavelengths (micron) '
200 1000 801 'LIN' = wavelengths (1st,last,howmany,how=LIN,INV,LOG,TAB)
' Refractive index of ambient medium '
1.0000 = NAMBIENT
' Effective Radii (micron) '
4.0691 4.0691 1 'LIN' = eff. radii (1st,last,howmany,how=LIN,INV,LOG,TAB)
' Define Incident Polarizations '
(0,0) (1.,0.) (0.,0.) = Polarization state e01 (k along x axis)
2 = IORTH (=1 to do only pol. state e01; =2 to also do orth. pol. state)
' Specify which output files to write '
0 = IWRKSC (=0 to suppress, =1 to write ".sca" file for each target orient.
' Specify Target Rotations '
0. 0. 1 = BETAMI, BETAMX, NBETA (beta=rotation around a1)
0. 90. 10 = THETMI, THETMX, NTHETA (theta=angle between a1 and k)
0. 0. 1 = PHIMIN, PHIMAX, NPHI (phi=rotation angle of a1 around k)
' Specify first IWAV, IRAD, IORI (normally 0 0 0) '
0 0 0 = first IWAV, first IRAD, first IORI (0 0 0 to begin fresh)
' Select Elements of S_ij Matrix to Print '
9 = NSMELTS = number of elements of S_ij to print (not more than 9)
11 12 21 22 31 33 44 34 43 = indices ij of elements to print
' Specify Scattered Directions '
'LFRAME' = CMDFRM (LFRAME, TFRAME for Lab Frame or Target Frame)
1 = NPLANES = number of scattering planes
0. 0. 180. 1 = phi, theta_min, theta_max (deg) for plane A

Re: Help with orientation issues

pcirrus 30 Jun 2016 21:28

If you look at the source code of tarell.f90 (or documentation) vectors a1 =(1,0,0) and a2=(0,1,0). Propagation of the incident k=(1,0,0), i.e. in the direction of x-axis. But your ellipsoid shape is (25,25,39), i.e. it is elongated in the z-direction. By rotating a1 around k (thetami, thatamx, ntheta) you are presenting the same shape andefficiencies will be the same. Try to rotate phi (90 90 1) first and theta as you do now. Now, there is a code TRAGET.f90 which can be used to generate the shape and plot it. It also returns a1, a2 vectors. This can be used to see the geometry of a target. Targets are just sets of numbers (ix,iy,iz) so they are easy to plot.

Alison Deatsch (guest) 01 Jul 2016 17:04

Thank you! That makes a lot of sense. I think maybe I was getting mixed up with some of my directions.

However, I just tried your solution and I got the exact same values for the Qext values in the qtable! What is happening?? Do you think I am just trying to make too small of a shape? I used the vtrconvert and paraview to look at my target, and it looks like an ellipsoid as I expect it should.

Any other ideas?

effective radius and size parameter

pcirrus 01 Jul 2016 17:33

Your effective radius is is about 4 and your wavelengths are 200-1000. The effective size parameter is x=2 pi r/lambda and is small ("Rayleigh regime"). In such a case slightly non-spherical particle will not make much difference for efficiencies (I suppose). But you need to make sure that you indeed rotate them correctly. Addressing another issue - the amount of dipoles by itself is not of such an importance (it to keep mkd small, i.e. < 0.5, but this is computational issue), the size parameter is of importance. If your problem is just to calculate light scattering by ellipsoids you should have a look at analytical codes (based on series) for light scattering by ellipsoids. Even if your problem is complex enough :) for DDSCAT (for example if your particles have inhomogeneous refractive index) it is good to compare DDSCAT to known non-spherical solvers; it will also give you some insight.

saroj bala (guest) 22 Dec 2017 05:54

respected sir, how we calculate effective radii in case of ellipsoid. pls help

draine 25 Dec 2017 14:31

The effective radius is *defined* to be the radius of an equal volume sphere (see eq. 4 of the UserGuide). An ellipsoid with "diameters" A,B,C along the three principal axes has volume V=(pi/6)*A*B*C and a_eff=(3*V/4*pi)^{1/3} = 0.5*(A*B*C)^{1/3}.

The physical value of a_eff is specified in ddscat.par.

The shape of the ellipsoid, and the number of dipoles that will be used, is determined when you set parameters SHPAR_1, SHPAR_2, SHPAR_3 (=A,B,C)

draine 25 Dec 2017 15:55

correction: SHPAR_1, SHPAR_2, SHPAR_3=(A/d, B/d, C/d) where d=dipole spacing.