Supplementary Material For: Correcting the Bias of Light Obscuration and Flow Imaging

Supplementary Material for: Correcting the Bias of Light Obscuration and Flow Imaging Particle Counters

Dean C. Ripple1 and Zhishang Hu2

1Biomolecular Measurement Division, National Institute of Standards and Technology, Gaithersburg, MD USA

2Center for Computational and Systems Biology, Institute of Biophysics, Chinese Academy of Sciences, Beijing, China

Overview of the procedure for bias correction of light obscuration measurements

The procedure for correcting the bias of light obscuration measurements has multiple steps, as given in detail in the Methods section of the main paper. These steps are outlined below to provide the reader with a concise summary of the complete procedure.

Procedure:

1.Experimentally determine the effective aperture  of the light obscuration apparatus (this needs to be done only once)

2.Using Mie theory and the value of , calculate the extinction efficiency factor for PSL beads versus diameter, Qpsl(d) (this needs to be done only once)

3.For a given particle suspension type, measure the phase map of a selection of particles

4.Convert phase map to refractive index versus particle diameter, using Eq. 2.

5.From the refractive index versus diameter plot, identify largest diameter with nearly constant refractive index

6.For particle diameters less than the cut-off diameter dc, model the particles as homogeneous spheroids of refractive index n(dc), and use Mie or RGD theory to predict the extinction efficiency factor Q.

7.For particle diameters greater than dc, obtain Q by scaling Q(dc) according to Eq. 7.

8.Obtain a diameter scale correction using Eq. 4.

Of these steps, only steps 3 and 4 are laborious.

Behavior of RGD Scattering Theory for Large Spheres and Spheroids

The extinction efficiency factor of a sphere, for detection by scattering at angles ≥, is given by(20)

(S1

withand. As explained in the main text, the finite detector acceptance angle assures that for sufficiently large values of d, u1 for all angles from to . For , is an oscillatory function with a magnitude that scales(40) as u4. Substituting Cu4, with C a constant, for in Eq. S1, one obtains:

,(S2

which is independent of d. Thus, in the limit of large d, Q approaches a valuethat is dependent on  but independent of diameter, andproportional to , or equivalently,(n)2.

For a spheroid, Eq. S1 still applies, provided that the parameter u is multiplied by a factor dependent on the aspect ratio of the spheroid and the angular orientation of the spheroid axis(38). Because this factor does not depend on the average size of the spheroid,Qapproaches a value independent of diameterin the limit of large spheroids by a similar analysis.

Figure S1. Theoretical variation with aperture angle of the ratio of indicated diameter for PMMA beads immersed in oil to the diameter of the same beads immersed in water. Beads were nominally 70 µm diameter; calculations and measurements were performed at a vacuum wavelength of 670 nm. Comparison with the experimentally measured ratio gives an empirical value for the effective aperture angle.

Figure S2. Particle counts for agitated HSA particles, as a function of diameter, with diameters transformed for each data set using silica bead calibration for FI and Mie and RGD theories with constant nfor LO. Error bars indicate standard error of the meanSEM of the counts of three different sample vials.

Figure S3. Workflow for determining particle refractive index difference Dn and effective diameter.

Sample code for RGD theory calculations of the extinction efficiency factor Q, for a light obscuration system with finite detector aperture. The code is written in Fortran 90. The function QTRAP(FUNC,A,B,TOLER,S,ERR) performs numerical integration of the function FUNC from limits A to B, to a tolerance TOLER, returning the value in S and an estimate of actual error in ERR. QTRAP, apart from slight modifications, was obtained from W. H. Press, B. P. Flannery, S. A. Teukolsky, and W. T. Vetterling, Numerical Recipes, Cambridge University Press, New York (1986), but is not included here due to copyright restrictions.

! LOModelSimple.f90

!****************************************************************************

!

! PROGRAM: LOModelSimple

program LOModelSimple

! Program finds the extinction efficiency factor, Q, for scattering from a

! homogeneous prolate spheroid, using Rayleigh-Gans-Debye theory and with a

! finite detector acceptance angle.

! Q is obtained by integrating the scattering intensity over the scattering

! angle.

!

! Reference:

! Barber PW and Wang D-S, "Rayleigh-Gans Debye applicability to scattering by

! nonspherical particles" (1978) Applied Optics 17:797-803

! Notation differs from Barber and Wang:

! Barber and Wang variable Program variable

! alpha theta0

! phi -phi0

! theta theta

! eta rmu

!

! User specifies orientation, aspect ratio, and refractive index of the

! spheroids. Program steps through a set of specified diameters and calculates

! the extinction efficiency factor.

!

! Spheroid orientation may be random, aligned with the fluid flow, or at a fixed

! angle. For fixed angle, the output value of Q is based not

! on the true projected area but instead on the projected area of a randomly

! oriented particle.

!

! Integrand within the scattering intensity integral varies over many orders of

! magnitude and has many near discontinuities. For this reason, higher

! order integral formulas are not suitable, and program uses a simple

! trapezoidal formula combined with fairly

! stringent convergence criteria (convergence needed for two calls in a row,

! and minimum number of calls).

!

! Program used in conjunction with Mie scattering results for spherical PSL

! beads to generate calibration curves for light obscuration

! apparatus.

!

! Equation validity is not built in--valid for generally small refractive index

! (<1.5 should be fine)

!

! Standard scattering theory applies, except integrals of scattering angle (to

! find Q) are from acceptance angle to pi, and not from zero.

!

! Variables:

! Although all variables are explicitly declared, some real variables begin

! with a lower case r as a reminder that

! they are real and not integers.

! Real variables

! a(42): semimajor axis of prolate spheroid (micrometers)

! b(42): semiminor axis of prolate spheroid (micrometers)

! rlambda: wavelength, entered as the vacuum wavelength, but then scaled

! to be wavelength in the fluid medium

! deltan: difference between the refractive index of the particle and the

! fluid medium

! rk: wave number (reciprocal micrometers)

! pi: constant pi

! pi2: pi/2

! rn: refractive index of the particle, scaled by dividing by the refractive

! index of the fluid medium

! theta: scattering angle (radians)

! aT,thetaT,theta0T,radiusT: temporary values of a, theta, theta0, and

! radius

! ThetaMax: maximum angle of integration, beyond which integrand is negligible

! (radians)

! Q: extinction efficiency factor

! phi0: azimuthal angle of particle axis

! Q1,Q2,Q3: components of Q for subranges of the total integration over

! theta

! ThetaMaxFr,FractionInteg: Parameters for breaking up the total range of

! integration into three pieces

! riMed: refractive index of the medium (coded as water in the visible)

! err: error in Q

! err1,err2,err3: error in Q1, Q2, and Q3

! FracError: fractional error in Q

! eps: small number used in check if particle is spherical or not

! rproj2: square of radius of sphere of same projected area as randomly

! oriented particle

! prefactor: grouped prefactor outside of integral for scattering intensity

! InitToler: absolute tolerance for error of Q integration for first call

! Qold: stored value of Q from previous calculation

! AAngle: acceptance angle of detector = minimum polar angle for light

! scattering to not impinge on obscuration detector

! radius(42): array of radius values at which calculations are performed

! rmu: reciprocal of aspect ratio (greek letter mu) of particles, rmu

! defined as a/b > 1

! toler: tolerance for integration of scattering intensity

! FracToler: fractional tolerance for integration of scattering intensity

! theta0: polar angle of major axis of spheroid

! Integer variables

! i,j: index variables in loops

! dummy: dummy index, used in data reads to skip over a value

! npts: number of points to calculate

! nradius: default number of radius values

! nradiusM: optional maximum number of radius values, to reduce calculation

! time at large radius values

!

! Character variables

! corient: parameter indicating whether particles have fixed (specified

! orientation), random (average over all angles), or aligned (axis along axis

! of fluid flow)

! orientation (f, r, or a):

! cdummy: dummy variable, used for character inputs that are not used in

!calculation

!

! Real functions

! FuncSIRGD: Scattering intensity for a prolate spheroid at fixed orientation

! FuncSIRanRGD: Scattering intensity for a prolate spheroid at random

! orientation

!****************************************************************************

implicit none

real(8) a(42),b(42),rlambda,deltan,rk,pi,rn,theta

real(8) aT,thetaT,theta0T,radiusT,ThetaMax,FuncSI,FuncSIRan,Q,phi0,Q1, Q2,Q3,ThetaMaxFr,FractionInteg

real(8) riMed,err,err1,err2,err3,FracError,eps,rproj2

real(8) FuncSIRGD,FuncSIRanRGD

real(8) prefactor,InitToler,Qold

real(8) AAngle,radius(42),rmu,toler,FracToler,theta0,pi2

integer i,j,dummy,npts,nradius,nradiusM

character*40 filenm

character corient,cdummy

external FuncSIRGD,FuncSIRanRGD

!Set up common block to have transfer values of theta, theta0, radius, all suffix T

common aT,thetaT,theta0T,radiusT,rmu,rn,rk,rlambda,phi0

data radius/0.20,0.24,0.29,0.35,0.42,0.50,0.60,0.75,0.90,1.10,1.30,1.60,1.90,2.30,2.70,3.30,3.90,4.70,5.70,6.90,8.3,10.,12.,14.,16.,18.,20.,22.,24.,26.,28.,30.,32.,34.,36.,38.,40.,42.,44.,46.,48.,50./

pi=3.14159265358979d0

pi2=pi/2.d0

eps=1.d-4

FractionInteg=10.d0

! Absolute tolerance for integral of Q, smallest radius = InitToler

InitToler=1.d-4

! Subsequent fractional tolerance on Q = FracToler; dropped down below by a

! factor of 10 for fixed orientation.

FracToler=1.d-3

! Wavelength of light in surrounding medium = rlambdamed

riMed=1.3309

write(*,*) 'Enter reciprocal of aspect ratio (>1): '

read(*,*) rmu

write(*,*) 'Enter in vacuo wavelength & absolute refractive index difference: '

read(*,*) rlambda,deltan

! Check on reasonable input value

if(deltan .gt. 1.d0) then

write(*,*) 'Enter refr. index DIFFERENCE from unity. (enter any character to continue)'

read(*,*) cdummy

stop

end if

write(*,*) 'Fixed, random, or aligned orientation (f, r, or a): '

read(*,*) corient

if (corient .eq. 'f' .or. corient .eq. 'a') FracToler=FracToler/10.d0

write(*,*) 'Enter instrument acceptance angle (rad): '

read(*,*) AAngle

write(*,*) 'Enter file name for output: '

read(*,*) filenm

open(unit=9,file=filenm)

write(9,'(A)') 'RGD theory'

nradius=42

write(9,*) 'Aspect ratio, in vacuo wavelength, abs. refr. index delta: '

write(9,'(f5.1,'', '',f6.2,'', '',f6.3)') rmu,rlambda,deltan

rlambda=rlambda/riMed

deltan=deltan/riMed

rn=1.d0+deltan

rk=2*pi/rlambda

write(9,'(A,f8.4)') 'Acceptance angle (rad): ',AAngle

if (corient .eq. 'f') then

write(9,*)'Major axis angle of spheroid: fixed'

write(9,*) 'S integrand based on prefactor with random proj. area. '

! choose theta0 = polar angle of major axis of spheroid

write(*,*) 'Enter angle of major axis of spheroid (deg.): '

read(*,*) theta0

write(9,'(A,f6.1)') 'Major axis angle of spheroid (deg.): ',theta0

! convert to rad

theta0=theta0*(pi/180.d0)

if (abs(theta0) .lt. 1.d-6) then

theta0=1.d-6

end if

elseif (corient .eq. 'a') then

write(9,*)'Major axis angle of spheroid: aligned'

write(9,*) 'S integrand based on prefactor with true proj. area. '

! spheroid aligned with flow field

theta0=90.d0

write(9,'(A,f6.1)') 'Major axis angle of spheroid (deg.): ',theta0

theta0=theta0*(pi/180.d0)

else

write(9,*)'Major axis angle of spheroid: random'

write(9,*) 'S integrand based on prefactor with true random proj. area. '

end if

write(*,*) 'Enter max. number of radii (42 max): '

read(*,*) nradiusM

nradius=min(nradiusM,nradius)

write(*,*) 'Index, eff. radius, squared proj. radius, Q, frac. error: '

write(9,*) 'eff. radius, squared proj. radius, Q, frac. error: '

! Loop over radii

do i=1,nradius

! calculate major & minor radii:

a(i)=rmu**(2.d0/3.d0)*radius(i)

b(i)=a(i)/rmu

radiusT=radius(i)

aT=a(i)

thetamax=(pi/180.d0)*min(90.,1000*rlambda/radiusT)+AAngle

! Integral done between acceptance angle and Thetamax. Integral is split into

! forward scattering subrange, & remainder out to ThetaMax,

! breaking at ThetamaxFr

ThetamaxFr=AAngle+(ThetaMax-AAngle)/FractionInteg

if ((corient .eq. 'r' .or. corient .eq. 'f') .and. rmu .gt. (1.d0+eps)) then

! rproj = radius of sphere of same proj. area as randomly oriented particle; see

! Jennings Proc. Roy Soc. A 419:137 (1988)

rproj2=(a(i)**2/(2*rmu))*(1/rmu + (rmu/sqrt(rmu**2-1))*atan(sqrt(rmu**2-1)))

else !corient .eq. 'a' or particle is a sphere

rproj2=a(i)*b(i)

endif

Prefactor=rk**4*a(i)**2*b(i)**4*((rn**2-1)**2/(rn**2+2)**2)/(pi*rproj2)

! Use last value of Q to set tolerance, unless i=1

if (i .eq. 1) then

toler=InitToler/Prefactor

else

toler=FracToler*Qold/Prefactor

end if

!Fixed particle orientation

if (corient .eq. 'f' .or. corient .eq. 'a') Then

!******************

! Integrate S integrand (FuncSIRGD)

theta0T=theta0

call QTRAP(FuncSIRGD,AAngle,ThetaMaxFr,toler,Q1,err1)

call QTRAP(FuncSIRGD,ThetaMaxFr,thetaMax,toler,Q2,err2)

call QTRAP(FuncSIRGD,ThetaMax,pi,toler,Q3,err3)

else

! Random orientation--average over all theta0 values.

! Integrate S integrand (FuncSIRanRGD) from acceptance angle to pi in three

! steps.

call QTRAP(FuncSIRanRGD,AAngle,ThetaMaxFr,toler,Q1,err1)

call QTRAP(FuncSIRanRGD,ThetaMaxFr,thetaMax,toler,Q2,err2)

call QTRAP(FuncSIRanRGD,ThetaMax,pi,toler,Q3,err3)

end if

FracError=(err1+err2+err3)/(Q1+Q2+Q3)

Q=Prefactor*(Q1+Q2+Q3)

Qold=Q

write(*,'(i3,'', '',F8.4,'', '',g11.3,'', '',g11.3,'', '',f6.3)') i,radius(i),rproj2,Q,FracError

! write(*,'(3(g11.3,'', ''))') Q1, Q2, Q3

write(9,'(f8.3,'', '',g12.4,'', '',g13.5,'', '',f6.3)') radius(i),rproj2,Q,FracError

end do

stop

end program LOModelSimple

real(8) function FuncSIRGD(theta)

! FuncSIRGD = Integral of phi dependent terms in scattering intensity

! Scattering angle theta and spheroid angle theta0 are fixed.

implicit none

real(8) aT,theta,thetaT,theta0T,radiusT,rmu,rn,rk,rlambda,toler

real(8) phiIntegral

!Set up common block to have transfer values of theta, theta0, radius, all suffix T

common aT,thetaT,theta0T,radiusT,rmu,rn,rk,rlambda

thetaT=theta

! In this function, spheroid angle theta0 is fixed--only need to call RGD integrand

call RGDphi(phiIntegral)

FuncSIRGD=phiIntegral

return

end

real(8) function FuncSIRanRGD(theta)

! Like FuncSIRGD, but integrate over all orientations of theta0

! FuncSIRanRGD = Integral of all theta0 and phi dependent terms in scattering

! intensity

! Scattering angle theta is fixed.

implicit none

real(8) sinRGD,sinRGDint,aT,theta,thetaT,theta0T,radiusT,rmu,rn,rk,rlambda

real(8) phiIntegral,pi2,pi4,toler,err

!Set up common block to have transfer values of theta, theta0, radius, all suffix T

common aT,thetaT,theta0T,radiusT,rmu,rn,rk,rlambda

external sinRGD

! pi2 = pi/2

! pi4 = pi/4

data pi2/1.570796327/

data pi4/0.785398163/

thetaT=theta

! Find typical magnitude of sinRGDint, and use to set tolerance

toler=1.d-2*sinRGD(pi4)

! Integrate over spheroid angle theta0 from 0 to pi/2

! err = error of integration

call qtrap(sinRGD,1.d-6,pi2,toler,sinRGDint,err)

! Integral of sin(theta0) from 0 to pi/2 = 1, so normalized sinRGDint is same as

! sinRGDint

! No factor of sin(theta), since that's carried in phiIntegral for RGD formulas.

FuncSIRanRGD=sinRGDint

return

end

real(8) function sinRGD(theta0)

! Integrand for averaging (integral of RGD scattering over phi) over theta0

! Parameters:

! theta0: scattering angle

! phiIntegral: integral of form factor squared and theta dependent terms over phi

! sinRGD: product of sin(theta0) and phiIntegral

! RGDphi: integral of scattering intensity over phi

implicit none

real(8) theta0,phiIntegral

real(8) aT,thetaT,theta0T,radiusT,rmu,rn,rk,rlambda

common aT,thetaT,theta0T,radiusT,rmu,rn,rk,rlambda

theta0t=theta0

call RGDphi(phiIntegral)

sinRGD=sin(theta0)*phiIntegral

return

end

subroutine RGDphi(phiIntegral)

! Average scattering integral over phi, but keep theta and theta0 fixed.

! Second-lowest integrand.

implicit none

real(8) toler, rk,aT,thetaT,theta0T,fIntegrand,phiIntegral,err,pi

real(8) radiusT,rmu,rn,rlambda

common aT,thetaT,theta0T,radiusT,rmu,rn,rk,rlambda

external fintegrand

data pi/3.1415925d0/

! estimate typical magnitude of integrand; assume angle-dependent terms have

! magnitude approximately 0.5.

toler=1.d-2*rk*aT*sin(thetaT/2.)

! perform phi integral

call QTRAP(fintegrand,0.d0,pi,toler,phiIntegral,err)

return

end

real(8) function fintegrand(phi)

! Returns theta and phi dependent integrand equal to form factor squared, times

! theta-dependent terms.

! Lowest level of integrand

! beta: angle defined by Barber and Wang paper

! theta0T: scattering angle

! thetaT: angle of spheroid major axis from z axis

! phi: azimuthal angle of spheroid, relative to plane of incident and scattered light

! f: form factor (3/u^3)(sin u - u cos u), same as Barber and Wang, Eq. 3

! fintegrand: f^2 (1 + cos^2(theta) sin(theta)

implicit none

real(8) rk,aT,thetaT,beta,theta0T,phi,u,f,radiusT,rmu,rn,rlambda

common aT,thetaT,theta0T,radiusT,rmu,rn,rk,rlambda

! Eqs. 5 and 4 of Barber and Wang:

beta=acos(-cos(theta0T)*sin(thetaT/2.) + sin(theta0T)*cos(thetaT/2.)*cos(phi))

u=2*rk*aT*sin(thetaT/2)*sqrt(cos(beta)**2 + (sin(beta)/rmu)**2)

if (u .lt. 1.d-2) then

f=1.d0

else

f=(3/u**3)*(sin(u) - u*cos(u))

end if

fIntegrand=(1 + cos(thetaT)**2)*sin(thetaT)*f**2

return

end

This paper is a contribution of the United States government and is not subject to copyright