Open segasai opened 3 years ago
This sounds great Sergey. Please provide a report demonstrating that your code does indeed better. You could do this by comparing the calibrated spectra form the same tile from night to night, with both calibration methods.
Thanks Julien, Looking at the code, do I understand correctly that in here https://github.com/desihub/desispec/blob/a4025dc2d7287e094f5dba1f081e317ecf86872d/py/desispec/scripts/stdstars.py#L492 normflux here is the main output that will be used in the next flux calibration steps ? So to the first order I need to substitute it by mine? (I'm trying to come up with smallest possible change at first)
Yes (stellar model spectrum, normalized to r-band magnitude, with radial velocity, and with dust extinction applied).
Here is the quick comparison of stdstars* from the default script, vs what I get with my quick hack https://github.com/segasai/desispec/commit/4fb835ee11ee716a1147e1c6e82afd2bf3e0eaee
There is certainly difference in model spectra. I obviously cant' claim one is better than other without more checks. But what I am quite certain is that I recover the stellar parameters better.
@julienguy Here's my preliminary report. I've taken the tile 80608, which seems to have at least 27 exposures in Blanc. I've put in a very crude implementation of fitting by my code https://github.com/segasai/desispec/tree/rvspecfit I've then separately fitted all the individual exposures of the tile and computed the stdstars-* files for the default pipeline and mine. I then looked at spreads between inferred spectra from different exposures by computing MAD (relative to median spectrum) . The MAD across all the stars/exposures for the default pipeline is 0.0043637753, for my code is 0.003242135, i.e .25% better. I attach the PDF file ( xx.pdf ) showing the comparison star by star. Each page is one star and black lines shows the MAD across many exposures for the current pipeline, red is what I get. It shows IMO that although there are a few stars where I show higher scatter among exposures than the current code (I have dig into that), it does perform better on average. See below some examples of spectral fits that I get (red is my model, blue are pixels that I mask)
TBH I'd hoped to pull more than 25% better, but I also see that the improvement is larger for higher S/N objects. I.e. if I only look at S/N>30 standards , I get it 32% better and for S/N>50 standars it's 41% better. That does make sense as at low S/N the spectrum doesnt' tell you that much anyway, and the default pipeline basically only relies on the color range, and in my case the Teff|color prior kicks in. For high S/N the template interpolation +details matter more. Also in my tests I haven't assesed the absolute calibration, or different fiber vs different fiber calibration, but only the repeatability.
If I had to guess, I think it's possible to improve things further with a bit of tuning. But it probably would only make sense, if there is enough interest in this approach.
This looks promising. Do you have a comparison of the measured colors of the stars from the imaging catalog and the synthetic colors for your best fit stellar models? I've been obtaining rather good prediction of the g-r color with rms typically <0.03, it would be interesting to see what you obtain (conversation continued outside of github ...)
Hi,
While implementing the flux calibration using gaia-only standards in stdstars code, I have finally understood what stdstars is doing underneath.
And after that I was wondering, what would people think about adding to this/replacing/testing with my spectral fitting/radial velocity code that I use for DESI https://github.com/segasai/rvspecfit
A few notes 1) My code takes as input not-flux calibrated spectra in 3 arms and chi^2 fits them to get RV,feh,logg,teff. 2) The code is pure python and is simple to use, and basically returns stellar parameters,RV,and the best-fit stellar model given arrays of fluxes/variances in multiple arms. 3) The code takes 20-40 haswell core-seconds per star and is trivially parallelizable. 4) The code just processed without issues > 10^6 blanc individual exposures/coadds/recent daily spectra.
As far as I can see, the advantages of at least trying this are 1) one can always flux calibrate the field as soon as there are stars observed. 2) More precise flux calibration due to more precise stellar parameters.
Am I missing something ?
Or is the existing code good enough ?