pro vbmv2tg,vbmvfile,vv,bb,lgt,lgg,bcc,err ; This routine attempts to determine the values of logTeff, logg that give ; specified values of V_abs and B-V. ; On input, it accepts the pathname vbmvfile, which is an idl .sav file ; containing vectors logt, logg, and corresponding arrays vabs, bmv, bolometric ; correction bc. ; It also accepts values ; vv = desired value of vabs ; bb = desired value of B-V ; On output, it returns estimates of the corresponding log(T), log(g), and BC. ; If any problem occurred during the calculation, then on return err=-1, else 0. ; constants wtv=0.1 ; range in B-V is about this much times range in Vabs, ; hence this is the relative weighting. dvmax=0.3 ; insist best point matches vabs better than this dbmax=0.06 ; ditto for B-V err=0 restore,vbmvfile ; contains arrays logg(nlg),logt(nlt),vabs(nlt,nlg),bmv(nlt,nlg), ; llogt(nlt,nlg),llogg(nlt,nlg),bc(nlt,nlg) nlg=n_elements(logg) nlt=n_elements(logt) ; find points closest to desired vabs, B-V. Weight B-V distance more strongly dv=vabs-vv db=bmv-bb dist=sqrt((dv*wtv)^2+db^2) md=min(dist,ix) ; find smallest distance ; test for adequate fit dv0=abs(dv(ix)) db0=abs(db(ix)) if(dv0 gt dvmax or db0 gt dbmax) then begin err=-1 goto,fini endif ; fit okay, so where is this point? boundaries of neighbor region? ixy=long(ix/nlt) ixx=ix-ixy*nlt ixbot=(ixx-1) > 0 ixtop=(ixx+1) < (nlt-1) iybot=(ixy-1)> 0 iytop=(ixy+1) < (nlg-1) ; retrieve the deltas within this box, fit planes dv1=dv(ixbot:ixtop,iybot:iytop) db1=db(ixbot:ixtop,iybot:iytop) ltb=llogt(ixbot:ixtop,iybot:iytop) lgb=llogg(ixbot:ixtop,iybot:iytop) bcb=bc(ixbot:ixtop,iybot:iytop) npt=n_elements(dv1) funs=fltarr(npt,3) funs(*,0)=1. funs(*,1)=reform(ltb) funs(*,2)=reform(lgb) wt=fltarr(npt)+1 cc1=lstsqr(reform(dv1,npt),funs,wt,3,rms,chisq,outpv,1) cc2=lstsqr(reform(db1,npt),funs,wt,3,rms,chisq,outpb,1) cc3=lstsqr(reform(bcb,npt),funs,wt,3,rms,chisq,outpc,1) ; now solve the linear equations to force both delta vabs and delta B-V =0 rhs=[-dv(ix),-db(ix)] aa=[[cc1(1:2)],[cc2(1:2)]] aai=invert(aa) sol=reform(rhs#aai) ; desired displacement from central logt, logg in \$ ; units of logt, logg ; check to see that displacements are smaller than one step in each coord. lt0=llogt(ix) lg0=llogg(ix) dlt=abs(logt(1)-llogt(0)) dlg=abs(logg(1)-llogg(0)) if(abs(sol(0)) gt dlt or abs(sol(1)) gt dlg) then begin err=-1 goto,fini endif ; if so, accept them and report result. lgt=lt0+sol(0) lgg=lg0+sol(1) bcc=cc3(0)+lgt*cc3(1)+lgg*cc3(2) fini: end