Doing the gaussian fitting in gnuplot, which turns out to be very easy indeed. Commands such as
gauss1(x) = A1*exp( -((x-x1)/sig1)**2 ) gauss2(x) = A2*exp( -((x-x2)/sig2)**2 ) double(x) = gauss1(x) + gauss2(x) FIT_LIMIT = 1e-6 fit double(x) 'oi-sp1d-11.001' via A1, x1, sig1, A2, x2, sig2
... gnuplot does the fitting ...
Final set of parameters Asymptotic Standard Error ======================= ========================== A1 = 1.73948 +/- 0.009477 (0.5448%) x1 = 29.5711 +/- 0.05439 (0.1839%) sig1 = 11.2308 +/- 0.08459 (0.7532%) A2 = 0.977997 +/- 0.01126 (1.152%) x2 = 5.95878 +/- 0.08081 (1.356%) sig2 = 7.55852 +/- 0.117 (1.548%) plot 'oi-sp1d-11.001', gauss1(x), gauss2(x), double(x)
Also, the fit parameters all get written to fit.log
Now, need to write a script to automatize the fitting to all the lines.
OK, done that - it works for most of the slits but has problems for a few. Wrote parse_fit_log.f90, which makes a table of the results (oi-fit-sky.log.table). Edited this to give a new table (oi-fit-sky-newstart.dat) with sane starting values for all slits. Idea is to use this to generate initial parameters for the fits, and to keep sig1 fixed at 8.0 km/s (or 6.4 for the 070 slit).
This is the command line I originally used (last week) to generate the 1D spectra that I have been using. It seems to predate the option of whether or not to use ufiddle - I think that means that it did:
for i in $(seq -w 1 60) ; do printf '%s\n' SPMoi jw4.csv $i '0 0' 2 1 \ 0 3 '-50 -20' '10 40' 0.05 3 '-150 150' oi-sp1d-$i n n \ oi-sp1d-$i.png/png | src/wspecplot1d ; done
Now have done the further-refined fits, using oi-fit-sky-again.sh, which uses parameter files (oi-fit-sky-??.start) generated by oi_fit_setup.f90 from oi-fit-sky-newstart.dat. The final fits are written to oi-fit-sky-??.save. I ended up using 7 km/s for the sky line widths.
On second thoughts, it would be better to rerun wspecplot1d without using ufiddle (also, say no to smoothing):
for i in $(seq -w 1 60) ; do printf '%s\n' SPMoi jw4.csv $i n '0 0' 2 1 \ 0 n 3 '-50 -20' '10 40' 0.05 3 '-150 150' oi-sp1d-nf-$i n n \ oi-sp1d-nf-$i.png/png | src/wspecplot1d ; done
This is now done. In the sii_shift_scrap.gnumeric spreadsheet I have calculated the implied ufiddle = helio_corr - x1. This is pasted into the dataset file jw4-sky.{gnumeric,csv}
for i in $(seq -w 1 60) ; do printf '%s\n' SPMoi jw4-sky.csv $i y '0 0' 2 1 \ 0 n 3 '-50 -20' '10 40' 0.05 3 '-150 150' oi-sp1d-$i n n \ oi-sp1d-$i.png/png | src/wspecplot1d ; done for i in $(seq -w 1 60) ; do printf '%s\n' SPMoi jw4-sky.csv $i 3 \ 2 1 0 n '-20 60' 20 oi-sky-$i /null | src/wmomplot ; done for i in $(seq -w 1 60) ; do printf '%s\n' SPMoi jw4.csv $i 3 \ 2 1 0 n '-20 60' 20 oi-tere-$i /null | src/wmomplot ; done
That made the individual moment files for each slit, with running average of 20 pixels, for both the new ufiddles and for teresa's ones. Then oi-peaks.sh makes graphs of these. The sky versions do look better than the original ones (although still a bit of variation).
Next thing to do is subtract the sky line from the spectra. This seems to work. I had to change from a gaussian to a hyper-gaussian to make it more square-topped. Otherwise, I was getting systematic residuals.
Used specplot2d to adjust the skyline height manually for each slit (edit jw4-sky.csv.sky_info) to give perfect subtraction. Sometimes a velocity shift was necessary. In such case, have to change ufiddle in jw4-sky.gnumeric (i.e., DON'T change peak velocity in jw4-sky.csv.sky_info since these are pegged to vhel). There are still some velocities that might be out a bit.
Command to produce image of each slit:
for i in $(seq -w 60) ; do printf '%s\n' SPMoi jw4-sky.csv $i y 3 \ '-80 80' '-200 200' 2 1 0 y n y '20 -0.2' n 2 n \ oi-sp2d-$i.png/png | src/wspecplot2d ; done