NB. Polarization across line(s) in a slowly rotating (spherical) NB. star using SYNSPEC model computations. lambda0 is the first NB. wavelength in A. #lambdas is the number of wavelength points. NB. The inclination is fixed at 90 deg. eqV is rotation vel (km/s) NB. 'atmos_id' is file name for stored SYNSPEC I(µ) and Q(µ) values NB. for a particular Temp & gravity which are found at: NB. http://www.astro.umd.edu/~jph/2ND_STELLAR_SPECTRUM/ NB. Usage (lambda0,#lambdas) SPDX 'atmos_id'; eqV NB. E.g., Z=. (1400,200) SPDX 'T20000g3.5';15 NB. 'line lams flx I Q U P'=. Z -->line= 'atmos_id max_vel SPDX=: dyad define 'atmos_id eqV'=: y NB. eqV = equatorial velocity (km/sec) fseq=: eqV%2.9979e5 NB. equatorial velocity/speed of light atmos_id read_atmos_pol x wd 'msgs' NB. clear print buffer th=. cdtr*th=. i.181 NB. 181 theta points, [0,pi] phi=. (1p1%180)*1.5*i: 60 NB. 121 phi points, -pi/2 < phi < pi/2 'mu xi'=. Xi th; phi Dshift=: fseq*(sin phi)*/sin th NB. (D-shift of each point)/lambda 'I Q'=. Atmos_Pol mu Q=. Q* cos +:xi NB. rotate Stokes to z-axis U=. Q* sin +:xi dA=. (sin BIN th)*DEL th dA=. mup*(($mup=.BIN"1 mu)$dA) NB. dA is mu*sin(th)*d(th) II=. 1p_1* +/(DEL phi)*(BIN +/"2 dA*BIN"2 I) QQ=. 1p_1* +/(DEL phi)*(BIN +/"2 dA*BIN"2 Q) UU=. 1p_1* +/(DEL phi)*(BIN +/"2 dA*BIN"2 U) POL=. 100*QQ%II y;lams;flx;II;QQ;UU;POL ) NB. Xi is angle between projected surface normal and NB. Z', the projected Z-axis. NB. Usage: Xi th;ph --> mu;xi Xi=: monad define 'th ph'=. y n=. ph ,&#th NB. dimensions: phi x del mu=. (cos ph)*/sin th NB. ph columns, delta rows xi=. atan (sin ph)*/tan th mu; xi ) NB. Interpolate I and Q for all values of mu. Atmos_Pol=: monad define nn=. #mur=. ,y IT=. QT=. (nn,(#lams))$0 NB. empty 2D array to be filled below mur=. ,mu NB. ravel mu's to single dimension DI=. lams splineA III NB. derivatives for spline interpolation DQ=. lams splineA QQQ i=. _1 while. (i=. >:i)radians NB. Read in atmosphere tables for I,Q and Flux for model NB. Teff & log g, and interpolate at the requested lambda. read_atmos_pol=: dyad define load '/your_path/splines.ijs' dir2=. '/your_path/2ND_STELLAR_SPECTRUM/' Z=. readB dir2,x NB. E.g., x = 'T20000g3.5' mu0=: >0{Z NB. tabulated values of mu lams=: >3{Z NB. tabulated wavelengths (all 161,014 of them!) flx=: >4{Z NB. flux= mu-weighted intensity integrated over mu 'III QQQ POL'=: >(5 6 7){Z NB. Tabulated Stokes I(mu0) & Q(mu0) 'lam0 iin'=. y ii0=. <: lams I. lam0 lams=: lams{~ ii=. ii0+ i. iin print 'Wavelength range:';({.lams);'to';({:lams);'A';' Max Dop shift:';(lam0*fseq) dup=. I. 0= DEL lams NB. find repeated wavelength entrys (every 118 pts) lams=: dup delete lams NB. and remove it/them flx=: dup delete ii{flx III=: |: dup delete ii{III QQQ=: |: dup delete ii{QQQ NB. POL=: |: dup delete ii{POL wd 'msgs' NB. clear print buffer ) NB. verb to delete items from a vector NB. Usage: (deletion indices) dlt (vector) delete=: dyad define k=. _1 wile. (#x)>(k=.>:k) do. y=. (k{x) dlt y x=. <: x end. y ) dlt=: {. , ] }.~ [: >: [ NB. delete one item