From d2bd345018e02c4ce5ac2b413c8efd7650e5d319 Mon Sep 17 00:00:00 2001 From: Robert Jedrzejewski Date: Mon, 8 May 2023 11:41:44 -0400 Subject: [PATCH 1/3] Add fix for HASP data sometimes creating a product with no good pixels, making the section that calculated the first and last good wavelengths throw an exception. Catch this case. Since this doesn't happen for Ullyses data, we don't need to catch this case in the Ullyses wrapper like we do in HASP --- ullyses/coadd.py | 51 ++++++++++++++++++++++++++---------------------- 1 file changed, 28 insertions(+), 23 deletions(-) diff --git a/ullyses/coadd.py b/ullyses/coadd.py index 3fc60c9..6279069 100644 --- a/ullyses/coadd.py +++ b/ullyses/coadd.py @@ -203,30 +203,35 @@ def coadd(self): if self.instrument == 'COS': self.output_varsum[indices[i]] = self.output_varsum[indices[i]] + variances[i] good_dq = np.where(self.output_exptime > 0.) - self.first_good_wavelength = self.output_wavelength[good_dq][0] - self.last_good_wavelength = self.output_wavelength[good_dq][-1] - nonzeros = np.where(self.output_sumweight != 0) - if self.instrument == 'COS': - # Using the variances (which only COS has) gives spikes in the error when the flux goes negative. - self.output_varsum[nonzeros] = np.where(self.output_varsum[nonzeros] < 0.5, 0.5, self.output_varsum[nonzeros]) - self.output_errors[nonzeros] = np.sqrt(self.output_varsum[nonzeros]) + if len(good_dq[0]) == 0: + self.first_good_wavelength = None + self.last_good_wavelength = None + print("Warning: No good pixels in product") else: - # For the moment calculate errors from the gross counts - self.output_errors[nonzeros] = np.sqrt(self.output_sumgcounts[nonzeros]) - self.output_flux[nonzeros] = self.output_sumflux[nonzeros] / self.output_sumweight[nonzeros] - # Calculate the conversion from counts to flux units for errors - conversion = self.output_flux[nonzeros] / self.sumnetcounts[nonzeros] - # Clean out NaNs from where flux and net are zero - good = np.where(~np.isnan(conversion)) - bad = np.where(np.isnan(conversion)) - lenwl = len(wavelength) - interpolated_values = np.interp(self.output_wavelength[bad], - self.output_wavelength[good], - conversion[good]) - conversion[bad] = interpolated_values - # Use conversion to calculate error in flux units (erg/cm^2/s/Ang) - self.output_errors[nonzeros] = self.output_errors[nonzeros] * conversion - self.signal_to_noise[nonzeros] = self.output_flux[nonzeros] / self.output_errors[nonzeros] + self.first_good_wavelength = self.output_wavelength[good_dq][0] + self.last_good_wavelength = self.output_wavelength[good_dq][-1] + nonzeros = np.where(self.output_sumweight != 0) + if self.instrument == 'COS': + # Using the variances (which only COS has) gives spikes in the error when the flux goes negative. + self.output_varsum[nonzeros] = np.where(self.output_varsum[nonzeros] < 0.5, 0.5, self.output_varsum[nonzeros]) + self.output_errors[nonzeros] = np.sqrt(self.output_varsum[nonzeros]) + else: + # For the moment calculate errors from the gross counts + self.output_errors[nonzeros] = np.sqrt(self.output_sumgcounts[nonzeros]) + self.output_flux[nonzeros] = self.output_sumflux[nonzeros] / self.output_sumweight[nonzeros] + # Calculate the conversion from counts to flux units for errors + conversion = self.output_flux[nonzeros] / self.sumnetcounts[nonzeros] + # Clean out NaNs from where flux and net are zero + good = np.where(~np.isnan(conversion)) + bad = np.where(np.isnan(conversion)) + lenwl = len(wavelength) + interpolated_values = np.interp(self.output_wavelength[bad], + self.output_wavelength[good], + conversion[good]) + conversion[bad] = interpolated_values + # Use conversion to calculate error in flux units (erg/cm^2/s/Ang) + self.output_errors[nonzeros] = self.output_errors[nonzeros] * conversion + self.signal_to_noise[nonzeros] = self.output_flux[nonzeros] / self.output_errors[nonzeros] return def get_targname(self): From 597d42e8e8a4772fae3e78d413f9ab7f34ec1e03 Mon Sep 17 00:00:00 2001 From: Robert Jedrzejewski Date: Mon, 5 Jun 2023 17:30:31 -0400 Subject: [PATCH 2/3] Added fix for STIS and CCD data too --- ullyses/coadd.py | 27 +++++++++++++++++++++------ 1 file changed, 21 insertions(+), 6 deletions(-) diff --git a/ullyses/coadd.py b/ullyses/coadd.py index 6279069..91f5642 100644 --- a/ullyses/coadd.py +++ b/ullyses/coadd.py @@ -344,8 +344,13 @@ def coadd(self): nonzeros = np.where(self.output_errors != 0.0) self.signal_to_noise[nonzeros] = self.output_flux[nonzeros] / self.output_errors[nonzeros] good_dq = np.where(self.output_exptime > 0.) - self.first_good_wavelength = self.output_wavelength[good_dq][0] - self.last_good_wavelength = self.output_wavelength[good_dq][-1] + if len(good_dq[0]) > 0: + self.first_good_wavelength = self.output_wavelength[good_dq][0] + self.last_good_wavelength = self.output_wavelength[good_dq][-1] + else: + self.first_good_wavelength = None + self.last_good_wavelength = None + print('No good data in product') return class CCDSegmentList(SegmentList): @@ -430,8 +435,13 @@ def coadd(self): nonzeros = np.where(self.output_errors != 0.0) self.signal_to_noise[nonzeros] = self.output_flux[nonzeros] / self.output_errors[nonzeros] good_dq = np.where(self.output_exptime > 0.) - self.first_good_wavelength = self.output_wavelength[good_dq][0] - self.last_good_wavelength = self.output_wavelength[good_dq][-1] + if len(good_dq[0]) > 0: + self.first_good_wavelength = self.output_wavelength[good_dq][0] + self.last_good_wavelength = self.output_wavelength[good_dq][-1] + else: + self.first_good_wavelength = None + self.last_good_wavelength = None + print('No good data in product') else: nelements = len(self.output_wavelength) self.output_sumsqerrors = np.zeros(nelements) @@ -449,8 +459,13 @@ def coadd(self): self.output_sumsqerrors[indices[i]] = self.output_sumsqerrors[indices[i]] + (err[i] * weight[i])**2 self.output_exptime[indices[i]] = self.output_exptime[indices[i]] + segment.exptime good_dq = np.where(self.output_exptime > 0.) - self.first_good_wavelength = self.output_wavelength[good_dq][0] - self.last_good_wavelength = self.output_wavelength[good_dq][-1] + if len(good_dq[0]) > 0: + self.first_good_wavelength = self.output_wavelength[good_dq][0] + self.last_good_wavelength = self.output_wavelength[good_dq][-1] + else: + self.first_good_wavelength = None + self.last_good_wavelength = None + print('No good data in product') nonzeros = np.where(self.output_sumweight != 0) self.output_flux[nonzeros] = self.output_sumflux[nonzeros] / self.output_sumweight[nonzeros] self.output_errors[nonzeros] = np.sqrt(self.output_sumsqerrors[nonzeros]) / self.output_sumweight[nonzeros] From 85f2ac43698a7d31f94ef86ced0d2f8657e18517 Mon Sep 17 00:00:00 2001 From: Robert Jedrzejewski Date: Fri, 28 Jul 2023 14:05:36 -0400 Subject: [PATCH 3/3] Fix a couple of issues that arose from HASP testing: Bug for STIS CCD data where input data maps to beyond the end of the product arrays Issue for some COS data that hasn't been processed with the latest calcos code so doesn't have variance columns in data table - use gcounts instead --- ullyses/coadd.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/ullyses/coadd.py b/ullyses/coadd.py index 91f5642..47bcffe 100644 --- a/ullyses/coadd.py +++ b/ullyses/coadd.py @@ -191,8 +191,11 @@ def coadd(self): weight = flux_weight[goodpixels] net_counts = segment.data['net'][goodpixels] * segment.exptime if self.instrument == 'COS': - variances = segment.data['variance_counts'][goodpixels] + segment.data['variance_bkg'][goodpixels] + segment.data['variance_flat'][goodpixels] - + try: + variances = segment.data['variance_counts'][goodpixels] + segment.data['variance_bkg'][goodpixels] + segment.data['variance_flat'][goodpixels] + except KeyError: + print('WARNING: Cant get variance keywords for COS data, using gross counts instead') + variances = gcounts flux = segment.data['flux'][goodpixels] for i in range(npts): self.sumnetcounts[indices[i]] = self.sumnetcounts[indices[i]] + net_counts[i] @@ -405,7 +408,7 @@ def create_output_wavelength_grid(self): self.delta_wavelength = max_delta_wavelength - wavegrid = np.arange(self.min_wavelength, self.max_wavelength, self.delta_wavelength) + wavegrid = np.arange(self.min_wavelength, self.max_wavelength+self.delta_wavelength, self.delta_wavelength) self.output_wavelength = wavegrid