Analytical Gaussian fit for 5 and 7 points profiles

[detoma]

The new wavelength positions for the 5 and 7 points are not equally spaced.

The code needs to handle it when it computes the d_lambda for the analytical gaussian fit and we need to implement a least-square gaussian fit to take advantage of the extra points.

The 3 wavelengths to use for the analytical gaussian fits needs to be specified. For the 5-points, we will use the 3 center wavelengths. For 1074 and 1079 they will be:

• 1074.59, 1074.70, 1074.81
• 1079.69, 1079.80, 1079.91

For the 7-points the we will use the wavelength position 2, 4 and 6, i.e., we will skip the two farther and the two closer points to the center wavelength. For 1074 and 1079 they will be:

• 1074.59, 1074.70, 1074.81
• 1079.69, 1079.80, 1079.91

that are the same wavelengths as in the 3-point wave program.

Tasks

• make wave region option specifying these 3 wavelengths
• if above option is present, use those 3 wavelengths; if not, use center 3 wavelengths
• split wing_offset option into {red,blue}_reference_wavelength options
• record wavelengths used for "fit": REF_RED and REF_BLUE
• record method used for "fit" in FITMETHD: "analytic", "3-term", "4-term"
• reprocess 20240409 to compare to 3-term and 4-term results
• Mike to debug code since something isnt working.

[jburkepile]

Since observing time is limited by the helicopter fly times, we have agreed that one high quality 5-point observation for the other 3 coronal lines served to the community (637 nm, 706 nm and 789 nm) are much better than two 3-point observations that do not capture the info needed to fit the line profile except under optimal conditions.

The observing program will scan the lines at these locations:

1. Fe X 637 nm: 637.33, 637.365, 637.40, 637.435, 637.47 nm (wavelength stepsize = 0.035 nm)
   Ben has captured this in recipe: 637_05wave_0.035step_2beam_16sums_4reps_BOTH.rcp

2. Fe XV 706 nm: 706.13, 706.165, 706.20, 706.235, 706.27 (wavelength stepsize = 0.035 nm )
   Ben has captured this in recipe: 706_05wave_0.035step_2beam_16sums_4reps_BOTH.rcp

3. Fe XI 789 nm: 789.315, 789.35, 789.40, 789.45, 789.485 (wavelength stepsize = 0.035, 0.05, 0.05, 0.035 nm)
   Ben has captured this in recipe: 789_05wave_0.05step_2beam_16sums_4reps_BOTH.rcp

[detoma]

We need to perform an analytical gaussian fit for all UCoMP data.

If only 3 wavelength position are present, this will be the end of it. If more than 3 wavelength positions are present, this analytical fit will be used as first guess for the least square fit.

The three wavelengths to be used for the analytical gaussian fit are not necessarily the three center wavelengths. We will use the closest wavelengths to the "nominal" reference wavelengths used for the 3-point data:

ref_blue = 1074.59
ref_red = 1074.81

ref_blue = 1079.69
ref_red = 1079.91

ref_blue = 789.33
ref_red = 789.47

ref_blue = 637.35
ref_red = 637.45*

ref_blue = 706.14
red_red = 706.26

This is what I use in my code to find the closest wavelength position:

wave_blue = nearest_element(ref_blue, wave, p_blue) 
wave_red = nearest_element(ref_red, wave, p_red)

where wave is the array with the wavelength positions and nearest_element.pro is the routine written by Mankoff

*reference lines updated per discussion below

[mgalloy]

Making a new option, wing_offset, for each wave region, which will be optional and undefined for wave regions not specified below. In those cases, the center three wavelengths will be used.

Wave region	center_wavelength	wing_offset
637	637.4	0.05
706	706.2	0.06
789	789.4	0.07
1074	1074.7	0.11
1079	1079.8	0.11

[mgalloy]

OK, here is potential problem for the current method. Center wavelength is 637.4 and offset is +/- 0.06, so 637.34 and 637.46. For the below file, the method uses extensions 4, 5, and 7:

level1$ ucomp cat 20240330.211133.ucomp.637.l1.p9.fts
  1  Corona Stokes IQUV [637.24 nm]  science  -----   637.24  -----
  2  Corona Stokes IQUV [637.28 nm]  science  -----   637.28  -----
  3  Corona Stokes IQUV [637.32 nm]  science  -----   637.32  -----
  4  Corona Stokes IQUV [637.36 nm]  science  -----   637.36  -----
  5  Corona Stokes IQUV [637.40 nm]  science  -----   637.40  -----
  6  Corona Stokes IQUV [637.44 nm]  science  -----   637.44  -----
  7  Corona Stokes IQUV [637.48 nm]  science  -----   637.48  -----
  8  Corona Stokes IQUV [637.52 nm]  science  -----   637.52  -----
  9  Corona Stokes IQUV [637.56 nm]  science  -----   637.56  -----

637.34 and 637.46 are both directly in the middle of two wavelengths and round-off is giving me the higher one in both cases.

[detoma]

This was a special scan. It can be fixed by selecting 637.40 +/- 0.05.

However, I see an issue on 2024 04 09, when the red and blue wavelengths were not symmetric
20240409.182654.38.ucomp.637.l0.fts 40 exts sci 637.33, 637.36, 637.40, 637.43, 637.47
20240409.213526.90.ucomp.637.l0.fts 40 exts sci 637.33, 637.36, 637.40, 637.43, 637.47

There is a not a number that works well in this case. We need to do the following:
We find the nearest wavelength using +/- 0.05 and if there are two equidistant values, we select the larger value for the blue wing and the smaller value for the red wing.

This should select:
637.36, 637.40, 637.44 for 2024 03 30
and
637.36, 637.40, 637.43, for 2024 04 09

Unfortunately, uneven spacing like on 2024 04 09 will happen again by mistake or on purpose to avoid certain wavelength positions.

[bberkeyU]

@detoma What are the "correct" human-selected wavelengths from the 20240330.211133 file?

Based on the recommendation that +/-.05 would work better, it sounds like the correct answer is to round toward the central wavelength. Is that the correct interpretation?

[mgalloy]

Can I say, round towards the central wavelength when within 0.001? Floating point numbers can't be reliably compared for equality/non-equality in these edge cases:

IDL> print, (0.55 - 0.45) gt (0.45 - 0.35)
   1

[detoma]

The scan for 2024 03 30 was correct. We scanned a wide wavelength range to find the line with equally spaced steps.

There was a human error in selecting wavelengths for 2024 04 09 which made the steps NOT equally spaced. This is likely to happen again.

We cannot round towards the central wavelength. In order for the analytic gaussian fit to work, we want the red and blue wings to be ~0.05-0.06 nm apart from the center wavelength for 637.40nm.

Try to implement what I suggested above. It should work.

We need to be careful going forward about adopting new wavelength scans.

[mgalloy]

It can be fixed by selecting 637.40 +/- 0.05.

So, 637 should always use +/- 0.05? or just for this day? Doesn't this just move the problem to different wavelengths?

I don't know how to reconcile the following two statements:

We find the nearest wavelength using +/- 0.05 and if there are two equidistant values, we select the larger value for the blue wing and the smaller value for the red wing.

and

We cannot round towards the central wavelength.

[detoma]

It is not always the three center wavelengths that are the best. See my updated comments to below. I think the best approach is to define independently the blue and red reference wavelengths.

The spacing adopted on 2024 04 09 for 637.40nm was a mistake. Likewise on 2024 04 09, we took an extra point in the 1079nm scan by mistake. We need to account for human errors and there may be cases when we want them unequally spaced. I am trying to find a strategy that works for what was done in 2021-2022 and 2024.

[detoma]

We need to define the blue and red reference wavelength independently for flexibility. This will allow us to handle wavelength scans that are intentionally or unintentionally non evenly spaced. With a proper epoch, this will fix the issue for 2024 04 09.

I will test again which are the best three wavelength for 637.40 the day after the eclipse and update the ticket later.

[detoma]

It seems the wavelengths are listed incorrectly in the summary file but are correct in the 637 file.

Mike: can you please look into this? Is that a rounding error of some sort?

The wavelength in the 637 file are the following:
637.330 637.365 637.400 637.435 637.470

For the day after the eclipse, we want to use these wavelengths for the analytical fit: 637.365 637.400 637.435

Our analytical gaussian fit expects the 3 wavelengths to be at the same distance, so I am glad the data is correct, otherwise the fit would have failed.

I updated the reference red and blue lines for 637.

[mgalloy]

@detoma Are you saying the wavelengths are not correct in the 20240409 637 summary file, i.e., 20240409.ucomp.637.files.txt? I have just update the pipeline to always give 3 decimal places, so now that file reads:

20240409.182654.38.ucomp.637.l0.fts        40 exts    sci  637.330, 637.365, 637.400, 637.435, 637.470
20240409.213526.90.ucomp.637.l0.fts        40 exts    sci  637.330, 637.365, 637.400, 637.435, 637.470

[detoma]

Thank you Mike. This looks correct.

[mgalloy]

The velocity was coming out with opposite sign compared to what expected. This was because of a bug, now fixed, that flipped the red and blue reference wings. Running a test in process.analytic-gaussian to have results to compare to process.

[mgalloy]

The test run is complete for 20240409 in process.analytic-gaussian.

[detoma]

This seems to work fine based on data in 2022 but cannot compare velocities in 2024 because the files in process/ and the ones in process.analytic-gaussian/ are processed with two different distortions and are very different.

I need two files in 2024 where the only difference is the analytical algorithm.

[mgalloy]

OK, process.old-dist has a version with the old distortion and the new analytic Gaussian calculation with reference wavelengths.

[detoma]

I cannot find two images with the same distortion: process and process.old-dist use a different distortion. The velocities are very different, but I cannot tell how much is from the distortion and how much from the wavelength selection.

I am going to pass this issue to Mike. He is the best one to test it and make sure it works since he wrote the code.

[mgalloy]

Those runs should both use the old distortion. From the process.old-dist run:

level1$ ucomp cat -r 20240409.214658.ucomp.1074.l1.p5.fts | grep DIST
DISTORTC=                    T / whether distortion corrected
DISTORTF= 'dist_coeff_20210708_1074.sav' / created Fri Jul 30 16:10:59 2021

The distortion file was not recorded previously, i.e., in process, because there was only one distortion file, dist_coeff_20210708_1074.sav.

[detoma]

Mike is better if you do this.