New ELUT correction and propagation of live time uncertainty

[paolomassa]

We implemented a new ELUT correction which

• uses daily calibration files provided by @nicHoch
• applies spectral dependent correction based on coarse fit of the STIX count spectra

Correction has been implemented for imaging and for spectroscopy. However, in the case of spectroscopy it works only for CPD files and not for spectrograms (ELUT correction for spectrograms is meaningless).

We also added code for propagating live time uncertainty.

Finally, we use a calibrated version of the BKG detector transmission in stx_convert_pixel_data. This calibrated transmission has been provided by @MurielStiefel.

TODO list:

• Implement function for automatic download of the most appropriate daily calibration file from a time range provided as input
• Discuss what to do with spectrogram files: leave stx_convert_spectrogram as it is or modify it in a way similar to what done for stx_convert_pixel_data?
• Test that results are consistent with previous version based on daily calibration files provided by Olivier

[ennosigaeus]

A few observations on the new calibration file search logic and the ELUT correction flow — mostly clarifications and one question about intent. Nice work on the overall restructuring.

[ennosigaeus]

Two more observations — one on code duplication and one on an edge case in the spectral index calculation.

[ennosigaeus]

Observation (SRM/OSPEX duplication): The SRM construction and OSPEX integration (lines 552–700, ~150 lines) is now inlined here. The old stx_convert_science_data2ospex.pro still exists and was not modified in this PR.

This means there are now two copies of the SRM/OSPEX logic — the new inline version (with the updated ELUT correction, grid transmission, and livetime propagation) and the old standalone function. Fixes to one won't propagate to the other.

Not blocking, but worth discussing: should stx_convert_science_data2ospex be deprecated, or should the shared logic be extracted so both paths call the same code?

[paolomassa]

This is addressed in 1a24f45. The function stx_convert_science_data2ospex looks deprecated now, should we remove it?

[ennosigaeus]

Edge case (zero or negative spectrum values): If any spectrum bin is zero or negative after background subtraction, alog(spectrum(i)) produces -Inf or NaN. The cleanup at lines 107–108 (where(~finite(...))) catches bad values in the final index_final array, but the intermediate e_weighted values (computed from index_tmp in the first pass, lines 64–66) could already contain NaN. Since e_weighted feeds the second-pass recalculation (lines 87–101), NaN can propagate to adjacent good bins via alog(e_weighted(i-1)).

The flat-spectrum case (all equal values) is fine — produces all-zero indices as expected. The concern is specifically about bins with zero or negative counts.

[ennosigaeus]

Code review

Found 2 issues:

1. Time-rebin loop drops the trailing group from ELUT correction. In stx_convert_pixel_data.pro the outer while (i lt n_elements(duration)-1) exits before the final start index can be appended, and idx_time_min = iall[0:-2] / idx_time_max = iall[1:-1]-1 then treats the last iall entry purely as a boundary marker. Example: with N=10 and delta_time_min producing groups of 3, iall = [0,3,6] ⇒ groups (0,2) and (3,5) are processed but bins 6..9 are never visited by the loop at line 412. The corresponding count/error arrays for those bins are left as zeros from the fltarr(...) allocations at lines 407–408. Suggest appending the final i (or n_elements(duration)) to iall after the outer loop so the last group is included.

STIX-GSW/stix/idl/processing/spectrogram/stx_convert_pixel_data.pro

Lines 386 to 410 in edd1eb8

	i=0
	j=0
	total_time=0
	iall=[]
	while (i lt n_elements(duration)-1) do begin
	while (total_time lt delta_time_min) and (i+j le n_elements(duration)-1) do begin
	total_time = total(duration[i:i+j])
	j++
	endwhile
	iall = [iall,i]
	i = i+j
	j = 0
	total_time = 0
	endwhile
	idx_time_min = iall[0:-2]
	idx_time_max = iall[1:-1]-1
	count_rates_elut = fltarr(count_rates.dim)
	count_rates_error_elut = fltarr(count_rates_error.dim)
	n_times_rebinned = n_elements(idx_time_min)

2. NaN/Inf propagation in stx_estimate_spectral_index is still unaddressed (flagged in my 2026-03-30 review, no reply). When BKG subtraction yields zero or negative `spectrum` bins, `alog(spectrum(i))` at lines 61, 74, 89, 97 produces `NaN`/`-Inf`. The cleanup at lines 107–108 catches only `index_final`; the `e_weighted` values computed via `h0/h1` at lines 64–66 and 77–79 (where `h1` can also approach 0 when `index_tmp ≈ -1`) feed the second-pass loops at 87–101 and can contaminate adjacent good bins via `alog(e_weighted(i-1))`. Worth either guarding the inputs (`spectrum > 0`) before `alog`, or cleaning `e_weighted` between the two passes.

STIX-GSW/stix/idl/processing/spectrogram/stx_estimate_spectral_index.pro

Lines 57 to 102 in edd1eb8

	if idx_peak gt 0 then begin
	for i=0,idx_peak-1 do begin
	this_x=alog(e_mean(i+1))-alog(e_mean(i))
	this_y=alog(spectrum(i+1))-alog(spectrum(i))
	index_tmp(i)=this_y/this_x
	;use initial value to calculate weighted center
	h0=(e_high(i)^(index_tmp(i)+2)-e_low(i)^(index_tmp(i)+2))/(index_tmp(i)+2)
	h1=(e_high(i)^(index_tmp(i)+1)-e_low(i)^(index_tmp(i)+1))/(index_tmp(i)+1)
	e_weighted(i)=h0/h1
	endfor
	endif
	;above peak (negative index)
	for i=idx_peak+1,n_elements(spectrum)-1 do begin
	this_x=alog(e_mean(i))-alog(e_mean(i-1))
	this_y=alog(spectrum(i))-alog(spectrum(i-1))
	index_tmp(i)=this_y/this_x
	;use initial value to calculate weighted center
	h0=(e_high(i)^(index_tmp(i)+2)-e_low(i)^(index_tmp(i)+2))/(index_tmp(i)+2)
	h1=(e_high(i)^(index_tmp(i)+1)-e_low(i)^(index_tmp(i)+1))/(index_tmp(i)+1)
	e_weighted(i)=h0/h1
	endfor
	;for peak, no value is calculated
	;set e_weighted to center
	e_weighted(idx_peak)=e_mean(idx_peak)
	;recalculate slopes with weighted center
	for i=idx_peak+1,n_elements(spectrum)-1 do begin
	this_x=alog(e_weighted(i))-alog(e_weighted(i-1))
	this_y=alog(spectrum(i))-alog(spectrum(i-1))
	index_final(i)=this_y/this_x
	endfor
	if idx_peak gt 0 then begin
	for i=0,idx_peak-1 do begin
	this_x=alog(e_weighted(i+1))-alog(e_weighted(i))
	this_y=alog(spectrum(i+1))-alog(spectrum(i))
	index_final(i)=this_y/this_x
	endfor
	endif

Also: your 2026-05-08 question — "The function `stx_convert_science_data2ospex` looks deprecated now, should we remove it?" — is still open. The new `stx_convert_spectrogram2ospex` supersedes it; happy to drop the old one in this PR if you prefer.

(Reviewed at `edd1eb8`. The earlier suspicions about the `corr_low` formula in `stx_elut_correction.pro` and the BKG-subtraction order in `stx_convert_spectrogram.pro` turned out to be consistent with the documented intent on closer inspection.)

[Generated with Claude Code]

_{- If this code review was useful, please react with thumbs-up. Otherwise, react with thumbs-down.}