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.

Sign up for GitHub

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

HLA-1362: Modified the minimum RMS computation used as a discriminant for choos… #1908

Merged
merged 4 commits into from
Nov 21, 2024
Merged

HLA-1362: Modified the minimum RMS computation used as a discriminant for choos… #1908

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
9e809c4
Modified the minimum RMS computation used as a discriminant for choos…
mdlpstsci Nov 4, 2024
8f33148
Include a minimum RMS value for the SBC detector as there is often
mdlpstsci Nov 21, 2024
eebac3c
Update the CHANGELOG
mdlpstsci Nov 21, 2024
2bba5f6
Update latest change with GitHub PR numbe
mdlpstsci Nov 21, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 4 additions & 0 deletions CHANGELOG.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,10 @@ number of the code change for that issue. These PRs can be viewed at:

3.7.2 (unreleased)
==================
- Include a minimum RMS value for the SBC detector, as is done for the other
detectors, as there seems to be a lot of noise in the source catalogs due to
a low detection threshold. [#1908]

- Force an exit with a return code, KEYWORD_UPDATE_PROBLEM, in try/exception block
when invoking refine_product_headers in hapsequencer.py and hapmultisequencer.py.
If the FITS header keywords are not properly updated, this can cause errors during
Expand Down
33 changes: 22 additions & 11 deletions drizzlepac/haputils/catalog_utils.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -289,6 +289,15 @@
is_zero_background_defined = True
log.info(f"Input image contains excessive zero values in the background. Median: {self.bkg_median:.6f} RMS: {self.bkg_rms_median:.6f}")

# Compute a minimum rms value based upon information directly from the data
minimum_rms = 1.0 / self.keyword_dict['texpo_time']
if np.nan_to_num(self.bkg_rms_median) < minimum_rms:
self.bkg_rms_median = minimum_rms
log.info("")
log.info(f"Minimum RMS of input based upon the total exposure time: {minimum_rms:.6f}")
log.info(f"Median RMS has been updated - Median: {self.bkg_median:.6f} RMS: {self.bkg_rms_median:.6f}")
log.info("")

Check warning on line 299 in drizzlepac/haputils/catalog_utils.py

View check run for this annotation

Codecov / codecov/patch

drizzlepac/haputils/catalog_utils.py#L295-L299 

Added lines #L295 - L299 were not covered by tests

# BACKGROUND COMPUTATION 2 (sigma_clipped_stats)
# If the input data is not the unusual case of SBC "excessive zero background", compute
# a sigma-clipped background which returns only single values for mean,
Expand Down Expand Up @@ -343,17 +352,18 @@
log.info("")

# Compute a minimum rms value based upon information directly from the data
if self.keyword_dict["detector"].upper() != 'SBC':
minimum_rms = self.keyword_dict['atodgn'] * self.keyword_dict['readnse'] \
* self.keyword_dict['ndrizim'] / self.keyword_dict['texpo_time']

# Compare a minimum rms based upon input characteristics versus the one computed and use
# the larger of the two values.
if (bkg_rms < minimum_rms):
bkg_rms = minimum_rms
log.info(f"Mimimum RMS of input based upon the readnoise, gain, number of exposures, and total exposure time: {minimum_rms:.6f}")
log.info(f"Sigma-clipped RMS has been updated - Background mean: {bkg_mean:.6f} median: {bkg_median:.6f} rms: {bkg_rms:.6f}")
log.info("")
if self.keyword_dict['detector'].upper() != 'SBC':
minimum_rms = np.sqrt(self.keyword_dict['numexp']) * self.keyword_dict['readnse'] / self.keyword_dict['texpo_time']
else:
minimum_rms = np.sqrt(np.clip(bkg_median * self.keyword_dict['texpo_time'], a_min=1.0, a_max=None)) / self.keyword_dict['texpo_time']
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
minimum_rms = np.sqrt(np.clip(bkg_median * self.keyword_dict['texpo_time'], a_min=1.0, a_max=None)) / self.keyword_dict['texpo_time']
minimum_rms = np.sqrt(
(bkg_median * self.keyword_dict['texpo_time']).clip(min=1.0)
) / self.keyword_dict['texpo_time']

I like using the ndarray.clip() method rather than the np.clip function. But if you have some reason to prefer the np.clip function version, your code will also work. In either case, you could omit the max value (or a_max, which seems the same as max), since it defaults to None.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I have no preference for either function/method, but since the function is already coded I will continue to use it :). I will note that one does have to specify both a_min and a_max (TypeError: clip() missing 1 required positional argument: 'a_max') as I discovered when I initially implemented the code!

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OK, sounds OK to me. (You can omit the max value with ndarray.clip.)

I have one other suggestion. I think the version of the rms that is used in the "excessive zero background" case should also have the minimum applied. The code is a little simpler because the median value is known to be zero. Here's a patch:

index e9ceb8fe..dbadf50f 100755
--- a/drizzlepac/haputils/catalog_utils.py
+++ b/drizzlepac/haputils/catalog_utils.py
@@ -289,6 +289,12 @@ class CatalogImage:
 
                 is_zero_background_defined = True
                 log.info(f"Input image contains excessive zero values in the background. Median: {self.bkg_median:.6f} RMS: {self.bkg_rms_median:.6f}")
+                # Compute a minimum rms value based upon information directly from the data
+                minimum_rms = 1. / self.keyword_dict['texpo_time']
+                if np.nan_to_num(self.bkg_rms_median) < minimum_rms:
+                    self.bkg_rms_median = minimum_rms
+                    log.info(f"Minimum RMS of input based upon the total exposure time: {minimum_rms:.6f}")
+                    log.info(f"Median RMS has been updated - Median: {self.bkg_median:.6f} RMS: {self.bkg_rms_median:.6f}")
 
         # BACKGROUND COMPUTATION 2 (sigma_clipped_stats)
         # If the input data is not the unusual case of SBC "excessive zero background", compute

I included code that handles the issue you ran across where the computed rms is NaN.


# Compare a minimum rms based upon input characteristics versus the rms computed and use
# the larger of the two values.
if (bkg_rms < minimum_rms):
bkg_rms = minimum_rms
log.info(f"Minimum RMS of input based upon the readnoise, number of exposures, and total exposure time: {minimum_rms:.6f}")
log.info(f"Sigma-clipped RMS has been updated - Background mean: {bkg_mean:.6f} median: {bkg_median:.6f} rms: {bkg_rms:.6f}")
log.info("")

# Generate two-dimensional background and rms images with the attributes of
# the input data, but the content based on the sigma-clipped statistics.
Expand Down Expand Up @@ -490,6 +500,7 @@
keyword_dict["texpo_time"] = self.imghdu[0].header["TEXPTIME"]
keyword_dict["exptime"] = self.imghdu[0].header["EXPTIME"]
keyword_dict["ndrizim"] = self.imghdu[0].header["NDRIZIM"]
keyword_dict["numexp"] = self.imghdu[0].header["NUMEXP"]
if keyword_dict["detector"].upper() != "SBC":
if keyword_dict["instrument"].upper() == 'WFPC2':
atodgn = self._get_max_key_value(self.imghdu[0].header, 'ATODGAIN')
Expand Down
Loading