Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Possible error in weighting scheme for calwf3 #569

Open
mgennaro opened this issue Jun 21, 2022 · 3 comments
Open

Possible error in weighting scheme for calwf3 #569

mgennaro opened this issue Jun 21, 2022 · 3 comments

Comments

@mgennaro
Copy link

mgennaro commented Jun 21, 2022

I believe the weighting scheme in cridcalc (for slope fitting calculations) might be erroneous. The code seems to imply that there is an implicit assumption of EQUAL sampling intervals. This is not the case in general for WFC3 up the ramp samples. I believe the error comes from using Mike Regan optimal weighting scheme which was developed for JWST, where the equal intervals assumption is correct, but was ported to WFC3 without paying attention to this detail. But maybe I am wrong.

Tagging a few WFC3 folks. I brought this up a while back, but I see it has not been fixed. I was reminded about it by reading a Tech memo for Roman!
@Vb2341 @protonchain

Line 1182 under https://github.com/spacetelescope/hstcal/blob/master/pkg/wfc3/calwf3/wf3ir/cridcalc.c

i.e.:

       wt=fabs(pow(fabs(k-((ndata-1)/2.))/((ndata-1)/2.),power));

The changes are the following:

  1. instead of k (the read sample index), we need x[k], i.e. the TIME of the k-th read

  2. instead of (ndata-1)/2 (the index of the “central” or middle read, or, for even number of total reads, the mid-point between the two central reads), we want ( x[ndata-1] + x[0] )/2, i.e. the mid-point of the interval IN TIME (not in read-index), i.e. the time of the last read minus the time of the first read in that interval divided by two

So I think the new line should be (you can and should make it in one line, because those variables I am wriritng are not defined, but for clarity in the email I will spell it out):

       xmid = (x[ndata-1] + x[0]) / 2;
       half_interval =  (x[ndata-1] - x[0]) / 2;
       deltas = (x - xmid)/half_interval;

       wt=fabs(pow(fabs(deltas),power));

And deltas can be written in one line as

       deltas = (2*x-x[ndata-1]-x[0]) / (x[ndata-1]-x[0] );

I think the external fabs could be removed because a positive number to any power stays a positive number, but there might be a numeric reason for that fabs to be there.

Just note that these “intervals” are not necessarily 1:1 w.r.t. the seqeunce of reads, as, at this stage, the pipeline may have already split the full interval in multiple sub-intervals due to cosmic ray identification.

@mdlpstsci
Copy link
Contributor

@mgennaro Does this ticket need further comment from other WFC3 team members?

Aside: I do not see any other ticket discussing this issue in Git. I will note the code associated with this bug was last touched approximately ten years ago.

@mgennaro
Copy link
Author

@mdlpstsci I confirm there is no other git issue taht I know of on this topic. This was part of an old email thread between you, me, Varun and Sylvia. The email I sent is from April 15 2019, and you replied on April 18 2019. You mentioned fixing this into this branch:

https://github.com/mdlpstsci/hstcal/blob/wfc3_experiment/pkg/wfc3/calwf3/wf3ir/cridcalc.c

That branch may have never been merged, I guess.
I just thought that opening an issue would be beneficial for tracking.

@mdlpstsci
Copy link
Contributor

Form Sylvia via email on 24 June 2022:
Hi Mario – thank you submitting that ticket. It’s on our wish list of things to do though unfortunately we haven’t had the resources to take it on yet. Thanks for thinking of us! -Sylvia

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

2 participants