pro polyfit green = GETCOLOR('green', 100) red = GETCOLOR('red', 100) enkey="B21" fn = '../products/detcnts.'+enkey+'.ALL.resid.fits' table = readfits(fn, header, EXTEN_NO=1, /SILENT) ;;tstart = TBGET( header, table, 'TSTART', /NOSCALE) ;;tstop = TBGET( header, table, 'TSTOP', /NOSCALE) rev = TBGET( header, table, 'REV', /NOSCALE) data = TBGET( header, table, 'RESID', /NOSCALE) time = LINDGEN(N_ELEMENTS(data)) rev_min=min(rev) rev_max=max(rev) plot_x=[] plot_y=[] for r=rev_min, rev_max do begin index = where(rev eq r, count) if(count lt 10) then continue mn = AVG(data(index)) push,plot_x,r push,plot_y,mn print,r,count,mn endfor plot,plot_x,plot_y ;, yrange=[-0.0001,0.0001] return ;;COEFF = ROBUST_LINEFIT(time, rate, YFIT, SIG, COEF_SIG) COEFF = ROBUST_POLY_FIT(time,data,7,YFIT,SIG, /DOUBLE) print,SIG print,COEFF ;mean_polyfit=AVG(YFIT) ;print,'mean fit',mean_polyfit ;print,'resid',SIG/mean_polyfit*100 poly = coeff;/mean_polyfit ;;save,FILENAME = enkey+'_poly.sav', poly, mean_polyfit scale=1.2 ;;yfit0 = (poly[0] + poly[1]*time + poly[2]*time^2 + poly[3]*time^3)*scale plot,time,data;, yrange=[-0.0001,0.0001] oplot,time,YFIT,color=green ;;oplot,time,YFIT0,color=red end