Add airmass dependent Bandpasses #444

jmeyers314 · 2024-01-17T07:18:09Z

Work in progress, looking for feedback.

Adds

an airmass kwarg to RubinBandpass. Performs linear interpolation of log(throughput) vs airmass X, which seems to be pretty good from investigation in a separate notebook.
EnvelopeRatio photon operator. This takes an "envelope" (X=1.0) bandpass object, and reads the target (X > 1.0) bandpass function from the base config. It then sets photon fluxes to zero with probability 1 - bandpass(wavelength)/envelope(wavelength), effectively changing the filter from X=1 to X>1 and still getting the Poisson statistics right.

I'm still not sure how to:

Ensure that any precomputed skyCatalog fluxes are using the X=1.0 throughputs.
Coerce galsim.drawImage() to use the envelope bandpass when initially sampling photons behind-the-scenes with WavelengthSampler.

I'll dig in tomorrow, but any advice appreciated.

jmeyers314 · 2024-01-18T20:58:44Z

Thinking more about how to coerce .drawImage to use the X=1.0 "envelope" bandpass instead of the actual @image.bandpass when using the EnvelopeRatio photon operator: we could place that logic directly into LSST_SiliconBuilder. That'd make this work for our most important case I think, but errors might creep in if using a different stamp type.

Alternatively, we could refactor to require that image.bandpass be the X=1 envelope, and only build the X>1 actual bandpass when configuring LSST_SiliconBuilder. That'd have the downside though that image.bandpass is slightly a lie.

@rmjarvis , do you have an opinion here?

cwwalter · 2024-01-24T18:47:22Z

Sorry have been busy with first trip to summit after moving. @rmjarvis did you see Josh's question above?

rmjarvis · 2024-01-24T18:55:50Z

I thought you would just specify the config bandpass as being the one without the airmass factor and let the effective bandpass with the atmosphere be emergent from the combination of that with the EnvelopeRatio. But I admit I haven't thought about this enough to understand why that doesn't work for wavelength sampling as Josh implied above.

cwwalter · 2024-01-29T17:51:28Z

Hi Josh,

What's your current status on this and what do you still need help/input on?

I think we wanted to start the rehearsal production really soon so we should figure out how to do it:

You asked above for advice on how to:

Ensure that any precomputed skyCatalog fluxes are using the X=1.0 throughputs.
Coerce galsim.drawImage() to use the envelope bandpass when initially sampling photons behind-the-scenes with WavelengthSampler.

Did you get any out of this issue?

jmeyers314 · 2024-01-29T20:15:58Z

Jim, Eli and I chatted a bit about this last Friday. I think we decided that instead of initially generating photons from the outer envelope bandpass and then probabilistically culling them to match the realized bandpass, we'd instead generate from the good 'ole X=1.2 throughputs and just modify photon fluxes (weights) to get to the appropriate bandpass. That's much easier than trying to coordinate envelope+realized bandpass. The downside is the Poisson statistics don't end up quite right, but since the airmass and QE variations are pretty small, I bet we can live with that. (And it's still better than what we had before).

I can probably deliver the above by ~tomorrow.

cwwalter · 2024-01-29T20:23:00Z

Intersting. So, the interpolation would return a weights table?

I thought I remembered the Galsim photon list can take a weight for entry, but can't find that now. Will GalSim itself apply the weight?

jmeyers314 · 2024-02-01T22:35:40Z

Okay, I think this one is ready to look at now. Tests pass for me locally (with branch tickets/DM-42513 of obs_lsst_data).

jmeyers314 · 2024-02-01T23:22:31Z

Okay, looks like the one CI failure is where it can't find the file that lives on an unmerged branch of obs_lsst_data. Good.

jchiang87

Thanks, Josh. I have some minor comments/requests having to do with handling the per-CCD QE curve side of things, otherwise lgtm.

imsim/bandpass.py

config/imsim-config.yaml

tests/test_stamp.py

jmeyers314 · 2024-02-06T23:13:49Z

@jchiang87 , I think I've addressed all your comments now. Let me know if there's anything else I should look at.

@rmjarvis , could you take a peek at this, especially the tests to make sure they address some of the concerns you raised at the last imsim technical meeting?

jchiang87

Thanks, Josh. LGTM.

rmjarvis

I have a few questions, but this looks fine to me too.

rmjarvis · 2024-02-08T17:46:26Z

imsim/bandpass.py

-        logger.warning("Warning: Using the old bandpass files from GalSim, not lsst.sims")
-        file_name = f"LSST_{band}.dat"
-    bp = galsim.Bandpass(file_name, wave_type='nm')
+    if sum([camera is not None, det_name is not None]) == 1:


Interested implementation of this. I generally prefer

if (camera is not None) != (det_name is not None):

Not sure if that's more or less obvious as to intent, but I figured I put it out there.

Yeah, I like that. Shortened it to if (camera is None) != (det_name is None):.

rmjarvis · 2024-02-08T17:51:33Z

imsim/bandpass.py

+        )
+
+        if camera is not None:
+            old_path = Path(os.getenv("OBS_LSST_DATA_DIR"))


Is this necessarily in the environment? If not, should maybe do something graceful if it isn't. Otherwise I think getenv returns None, and then Path probably takes this to be the current directory, which might lead to something odd happening downstream.

Turns out Path(None) raises a TypeError. I think it's okay to panic here; user asked for a specific detector bandpass but we're unable to comply. Reraised a ValueError with a more helpful message.

rmjarvis · 2024-02-08T17:55:51Z

imsim/stamp.py

@@ -720,16 +743,32 @@ def draw(self, prof, image, method, offset, config, base, logger):
            # In case we had to make a bigger image, just copy the part we need.
            image += fft_image[image.bounds]
            base['realized_flux'] = fft_image.added_flux
+            # print(f"{fft_image.added_flux =                          }")
+            # print(f"{fft_image.array.sum() =                         }")


Probably could remove most of these debugging print statements now.

rmjarvis · 2024-02-08T18:04:56Z

tests/test_stamp.py


    # There is also a code path where the xsize,ysize is dictated by the calling routine.
    del config['stamp']['size']
    image, _ = galsim.config.BuildStamp(config, obj_num, logger=logger, do_noise=False,
                                        xsize=128, ysize=128)
    print(obj_num, image.center, image.array.shape, image.array.sum())
    assert image.array.shape == (128,128)
-    assert 290_000 < image.array.sum() < 320_000  # 306675
+    assert 320_000 < image.array.sum() < 350_000  # 338_265


These all changed because of something in the sims code I suppose? New version?

Not the sims code, this is b/c I added kwargs to imsim-config.yaml to use the measured QE curves in the obs_lsst_data repo instead of the fiducial QE curve from the throughputs repo. Here's the difference (for R22):

rmjarvis · 2024-02-08T18:05:58Z

tests/test_stamp.py

@@ -446,7 +448,127 @@ def new_toplevel_only(self, object_types):
    image, _ = galsim.config.BuildStamp(config, obj_num, logger=logger, do_noise=False)
    print(obj_num, image.center, image.array.shape, image.array.sum())
    assert image.array.shape == (32, 32)
-    assert 5 < image.array.sum() < 10  # 8
+    assert 7 < image.array.sum() < 13  # 10


How about this? Do you understand why this changed?

I think same as above: real QE measurements are better than the baseline.

jmeyers314 added 2 commits January 16, 2024 23:02

Add airmass throughput interpolator

0ec8d15

Add EnvelopeRatio

6b4faaa

jmeyers314 added 13 commits January 30, 2024 15:03

Switch to reweighting photons

176a261

Do bandpass reweighting in LSST_SiliconBuilder

2af952b

Keep bandpass finite

f8be010

Start with linear Bandpass interpolation

09ea0ad

bandpass_fiducial -> fiducial_bandpass

17d547e

whitespace

1f3962c

checkpoint

deaa5fc

get fft-mode working with airmass-dependent bandpasses

81f5d2f

Add QE-dependent bandpasses

587eebf

Cleanup

e116d03

Cleanup

5c98a6a

Make compatible with skyCatalogs v1.6.0rc2

c3cbe8d

Use camera instead of instrument_name

4273c84

jmeyers314 marked this pull request as ready for review February 1, 2024 22:35

jmeyers314 requested review from rmjarvis and jchiang87 February 1, 2024 22:35

jchiang87 requested changes Feb 2, 2024

View reviewed changes

imsim/bandpass.py Show resolved Hide resolved

imsim/bandpass.py Outdated Show resolved Hide resolved

imsim/bandpass.py Outdated Show resolved Hide resolved

config/imsim-config.yaml Show resolved Hide resolved

tests/test_stamp.py Outdated Show resolved Hide resolved

jmeyers314 added 3 commits February 2, 2024 12:24

Document AtmInterpolator

ce585b1

Warning points to rubin_sim not lsst.sims

ecc5fc4

Always compute fiducial_bandpass

b80347c

jmeyers314 added 3 commits February 2, 2024 12:35

test against obs_lsst_data LsstCam

779ca69

use initial_flux_bandpass

94faa0c

Use CCD-dependent QE by default

7cf23a1

jchiang87 approved these changes Feb 8, 2024

View reviewed changes

rmjarvis approved these changes Feb 8, 2024

View reviewed changes

jmeyers314 added 3 commits February 8, 2024 17:11

Message for missing OBS_LSST_DATA_DIR

286bf8c

Simplify double None check

603febc

remove debugging prints

351177c

jmeyers314 merged commit 1a05c40 into main Feb 9, 2024
3 checks passed

jmeyers314 deleted the airmass branch February 9, 2024 01:48

jchiang87 mentioned this pull request Mar 7, 2024

Allow the atmospheric transmission curves to vary by airmass #179

Closed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Add airmass dependent Bandpasses #444

Add airmass dependent Bandpasses #444

jmeyers314 commented Jan 17, 2024

jmeyers314 commented Jan 18, 2024

cwwalter commented Jan 24, 2024

rmjarvis commented Jan 24, 2024

cwwalter commented Jan 29, 2024

jmeyers314 commented Jan 29, 2024

cwwalter commented Jan 29, 2024

jmeyers314 commented Feb 1, 2024

jmeyers314 commented Feb 1, 2024

jchiang87 left a comment

jmeyers314 commented Feb 6, 2024

jchiang87 left a comment

rmjarvis left a comment

rmjarvis Feb 8, 2024

jmeyers314 Feb 9, 2024

rmjarvis Feb 8, 2024

jmeyers314 Feb 9, 2024

rmjarvis Feb 8, 2024

rmjarvis Feb 8, 2024

jmeyers314 Feb 9, 2024

rmjarvis Feb 8, 2024

jmeyers314 Feb 9, 2024

Add airmass dependent Bandpasses #444

Add airmass dependent Bandpasses #444

Conversation

jmeyers314 commented Jan 17, 2024

jmeyers314 commented Jan 18, 2024

cwwalter commented Jan 24, 2024

rmjarvis commented Jan 24, 2024

cwwalter commented Jan 29, 2024

jmeyers314 commented Jan 29, 2024

cwwalter commented Jan 29, 2024

jmeyers314 commented Feb 1, 2024

jmeyers314 commented Feb 1, 2024

jchiang87 left a comment

Choose a reason for hiding this comment

jmeyers314 commented Feb 6, 2024

jchiang87 left a comment

Choose a reason for hiding this comment

rmjarvis left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment