pro m67_color_fit,names,aa,bb,rms,gmrf,dataf,ngood,zdb=zdb,xform=xform ; This routine uses m67_basel_colors to get synthetic colors, RA, Dec for ; all of Latham's stars in M67. It then matches the star positions to ; those in the stars database to obtain magnitudes for all matching ; stars in all available filters. It creates the colors ; u-g (?), g-r, r-i, i-z, z-J, J-H, H-K, g-Gred (?), and g-D51 ; for the observed stars, and for each color it performs a robust fit ; to determine coefficients aa,bb such that, as nearly a possible, ; aa_(g-r) + bb_(g-r)*((g-r_bas)-.5) + (g-r_bas) = (g-r_obs) ; aa_(r-i) + bb_(r-i)*((g-r_bas)-.5) + (r-i_bas) = (r-i_obs) ; etc. ; Thus, if the synthetic colors are an optimal representation of the observed ; colors, then aa = bb = 0 for all colors. ; On return, rms(ncolor) contains the robust rms magnitude of the fit ; residuals. ; Array dataf contains the observed colors for stars that were retained ; in the fit for each color, gmrf contains the basel g-r colors for these stars, ; and ngood(9) contains the number of retained stars. ; The routine makes a plot of the fitted values vs g-r for each color, ; overplotted with the fitted function. The plot is entitled Color_Fit.ps ; If keyword zdb is set, its value overrides zdbase ; If keyword xform is set, the synthetic colors are transformed using ; the coefficients found in file xform before fits are done. ; ;******************************* ; Hacked to provide only constant offsets, not linear fits, 7 July 2008 ; This version used for colblock_cas3 - related reduction, at least. ;******************************* ; constants zdbase=getenv('ZDBASE') if(keyword_set(zdb)) then zdbase=zdb starpath=zdbase+'/survey/stars/' m67db=['t131+010','t131+011','t131+012','t132+010','t132+011','t132+012',$ 't133+010','t133+011','t133+012'] ; tiles to search ntiles=9 rcap=10./60. ; capture radius = 2 arcsec expressed in arcmin astrolib common models,colblock,teffu,logzu,loggu ;restore,'~/d/basel/color_block1.sav' ;restore,'~/d/basel/colblock_castelli_smth.sav' restore,'~/d/basel/colblock_cas3_smth.sav' cblock=colblock ; get the Latham/Montgomery-derived properties of the M67 stars. m67_basel_colors,ra,dec,v,bmv,teff,logg,bas_colors,names,/mix bas_colors=transpose(bas_colors) if(keyword_set(xform)) then xform_colors,bas_colors,xform,bas_colors bas_colors=transpose(bas_colors) nst=n_elements(ra) ; The bas_colors have untraceable provenance, so ignore them and create new ones ; using colors_blk_intrp. This is more consistent, anyway, because this routine ; will be used to provide colors for star classification. nst=n_elements(teff) logzi=fltarr(nst)-.1 ; appropriate for M67 colors_blk_intrp,cblock,teffu,loggu,logzu,teff,logg,logzi,colsout,$ xform=xform red_vector,3.1*.09,dcolor,xform=xform stop colsred=colsout+rebin(dcolor,9,nst) bas_colors1=transpose(colsred) stop ; find matching stars in the 'stars' database mago=fltarr(nst,10)-1. for i=0,ntiles-1 do begin path=starpath+m67db(i) dbopen,path for j=0,nst-1 do begin ll=dbcircle(ra(j)/15.,dec(j),rcap,dis) if(ll(0) gt 0) then begin dbext,ll(0),'u,g,r,i,z,j,h,k,gred,d51',u,g,r,ii,z,j2m,h2m,k2m,gred,d51 mago(j,*)=reform([u,g,r,ii,z,j2m,h2m,k2m,gred,d51],1,10) endif endfor endfor ; strip out unmatched stars s=where(mago(*,1) gt 0.,ns) if(ns gt 0) then begin mago=mago(s,*) ra=ra(s) dec=dec(s) v=v(s) bmv=bmv(s) teff=teff(s) logg=logg(s) bascm=bas_colors1(s,*) endif ; make observed colors corresp to synthetic ones obsc=fltarr(ns,9) obsc(*,0)=mago(*,0)-mago(*,1) obsc(*,1)=mago(*,1)-mago(*,2) obsc(*,2)=mago(*,2)-mago(*,3) obsc(*,3)=mago(*,3)-mago(*,4) obsc(*,4)=mago(*,4)-mago(*,5) obsc(*,5)=mago(*,5)-mago(*,6) obsc(*,6)=mago(*,6)-mago(*,7) obsc(*,7)=mago(*,1)-mago(*,8) obsc(*,8)=mago(*,1)-mago(*,9) s1=where(obsc lt (-2.) or obsc gt 40. or obsc eq 0.,ns1) if(ns1 gt 0) then obsc(s1)=99.9 ; do the fits aa=fltarr(9) bb=fltarr(9) rms=fltarr(9) gmrf=fltarr(300,9) dataf=fltarr(300,9) ngood=lonarr(9) for i=0,8 do begin s=where(obsc(*,i) lt 5. and obsc(*,i) gt (-2.),ns) if(ns gt 5) then begin xx=bascm(s,1)-0.5 yy=obsc(s,i)-bascm(s,i) wt=fltarr(ns)+1. funs=fltarr(ns,2) funs(*,0)=1. funs(*,1)=xx ; cc=lstsqr(yy,funs,wt,2,rmsf,chisq,resid,1) cc=[mean(yy)] resid=yy-cc(0) quartile,resid,med,q,dq s1=where(abs(resid/dq) ge 4./1.349,ns1) if(ns1 gt 0) then wt(s1)=0. ; cc=lstsqr(yy,funs,wt,2,rmsf,chisq,resid,1) sgg=where(wt gt 0.,nsgg) cc=[mean(yy(sgg))] resid=yy-cc(0) rmsf=stddev(resid(sgg)) aa(i)=cc(0) ; bb(i)=cc(1) bb(i)=0. rms(i)=rmsf s2=where(wt ne 0.,ns2) dataf(0:ns2-1,i)=yy(s2) gmrf(0:ns2-1,i)=funs(s2,1)+.5 ngood(i)=ns2 endif endfor ; make the plot psll,name='Color_Fit_cas3.ps' xtit='Basel g-r color' xr=[0.2,1.4] yr=[-0.25,.25] stop !p.multi=[0,3,3] for i=0,8 do begin nn=ngood(i)-1 if(nn gt 10) then begin plot,gmrf(0:nn,i),dataf(0:nn,i),psym=1,xtit=xtit,ytit=names(i),xran=xr,$ yran=yr,/xsty,/ysty,charsiz=1.5,symsiz=0.7 oplot,gmrf(0:nn,i),aa(i)+bb(i)*(gmrf(0:nn,i)-0.5) endif else begin plot,gmrf(0:10,1),dataf(0:10,1),xtit=xtit,ytit=names(i),xran=xr,yran=yr,$ /xsty,/ysty,/nodata,charsiz=1.5 xyouts,0.60,0.,'No Data',charsiz=1.5 endelse endfor psend sun !p.multi=[0,1,1] stop end