; Please cite Schlawin et al. (2010) if you make use of this IDL code ; function ellk,k ; Computes polynomial approximation for the complete elliptic ; integral of the first kind (Hasting's approximation): m1=1.d0-k^2 a0=1.38629436112d0 a1=0.09666344259d0 a2=0.03590092383d0 a3=0.03742563713d0 a4=0.01451196212d0 b0=0.5d0 b1=0.12498593597d0 b2=0.06880248576d0 b3=0.03328355346d0 b4=0.00441787012d0 ek1=a0+m1*(a1+m1*(a2+m1*(a3+m1*a4))) ek2=(b0+m1*(b1+m1*(b2+m1*(b3+m1*b4))))*alog(m1) return,ek1-ek2 end function ellec,k ; Computes polynomial approximation for the complete elliptic ; integral of the second kind (Hasting's approximation): m1=1.d0-k^2 a1=0.44325141463d0 a2=0.06260601220d0 a3=0.04757383546d0 a4=0.01736506451d0 b1=0.24998368310d0 b2=0.09200180037d0 b3=0.04069697526d0 b4=0.00526449639d0 ee1=1.d0+m1*(a1+m1*(a2+m1*(a3+m1*a4))) ee2=m1*(b1+m1*(b2+m1*(b3+m1*b4)))*alog(1.d0/m1) return,ee1+ee2 end function ellpic_bulirsch,n,k ; Computes the complete elliptical integral of the third kind using ; the algorithm of Bulirsch (1965): kc=sqrt(1d0-k*k) & p=n+1d0 if(min(p) lt 0d0) then print,'Negative p' m0=1d0 & c=1d0 & p=sqrt(p) & d=1d0/p & e=kc L1: f = c & c = d/p+f & g = e/p & d = (f*g+temporary(d))*2d0 p = g + temporary(p) & g = m0 & m0 = kc + temporary(m0) if(max(abs(1d0-kc/g)) gt 1.d-13) then begin kc = 2d0*sqrt(e) & e=kc*m0 goto,L1 endif return,.5d0*!dpi*(c*m0+d)/(m0*(m0+p)) end pro chrom_exact,b,p,lc ; Please cite Schlawin et al. (2010) if you make use of this IDL code ; ; IDL code for computing the transit of a chromosphere by a spherical ; planet. ; ; Input: ; b vector of impact parameter values in units of the stellar radius ; p R_p/R_* = ratio of planet radius to stellar radius ; Output: ; lc light curve of a chromospheric transit. ; a=dblarr(n_elements(b)) indx=where(b+p lt 1d0) if(indx[0] ge 0) then begin & k=sqrt(4d0*b[indx]*p/(1d0-(b[indx]-p)^2)) & $ a[indx]=4d0/sqrt(1d0-(b[indx]-p)^2)*(((b[indx]-p)^2-1d0)*ellec(k)$ -(b[indx]^2-p^2)*ellk(k)+(b[indx]+p)/(b[indx]-p)$ *ellpic_bulirsch(4d0*b[indx]*p/(b[indx]-p)^2,k)) & endif indx=where(b+p gt 1d0 and b-p lt 1d0) if(indx[0] ge 0) then begin & k=sqrt((1d0-(b[indx]-p)^2)/4d0/b[indx]/p) & $ a[indx]=2d0/(b[indx]-p)/sqrt(b[indx]*p)*(4d0*b[indx]*p*(p-b[indx])*ellec(k)$ +(-b[indx]+2d0*b[indx]^2*p+p-2d0*p^3)*ellk(k)$ +(b[indx]+p)*ellpic_bulirsch(-1d0+1d0/(b[indx]-p)^2,k)) & endif lc=1d0-(4d0*!dpi*(p gt b)+a)/4d0/!dpi return end