No Photons Produced if sp.min() or sp.max() fails

[FredJenningsAstro]

Hi,

I've come across something weird that's probably a very niche case, but it took a while for me to narrow it down and was something behaving not the way I expected, so I though I'd bring it up.

I have a sphere object created from a dataset which I filter with a particle cut using a function very similar to the one given as an example here; https://hea-www.cfa.harvard.edu/~jzuhone/pyxsim/source_models/thermal_sources.html?highlight=filtering

I then do the following because I want to record the filtered field extremas;

self._logger.info(f"ds filtered field min {('filtered_gas',cut_field[1])} = {self.sp.min([('filtered_gas',cut_field[1] )])}")

where _cutfield is looped through ('gas','density'), ('PartType0', 'StarFormationRate') etc...

This is in a try/except because some of the fields (associated with PartType0) do not belong to the derived_field since the filtered_type='gas' in the particle cut. This is a lazy way of doing it, but I thought no biggie if it fails for a given field, since this line is only logging and is not actually doing anything.

However, it turns out that if this line fails, then if pyxsim.make_photons() is called with the sphere as the data source after this, 0 photons will be generated. The only way I've found to fix this apart from being more careful in the logging line is to do something like sp2 = ds.sphere(sp.center, sp.radius) and then calling make_photons on sp2. I assume this is messing something up in the yT cache which then causes the photon generation to fail in pyxsim?

[jzuhone]

What this probably means is that the failure of the logging line is signaling that something has gone wrong in the setup of the particle filter. Could you paste a full code example?

[FredJenningsAstro]

Hi John, thanks very much as always for the help.

I'm sorry, now looking back at my particle filter I have found where it is going wrong. Embarrassingly a wrong index into (field_type, field) meant that the requires argument in yt.add_particle_filter was being passed an empty list. Fixing this to pass the gas field names properly seems to fix the whole issue.

I am a little confused around the particle filter though; firstly, passing an empty list to requires didn't stop the filtered field seemingly being set-up correctly, and when I checked, say, the density and temperature extrema of the filtered field they were as expected. This is why I didn't spot this indexing typo first time around, and didn't suspect the filter.

Also, like in the cookbook example, I filter both gas and PartType0 in the filter, with the filtered_type = "gas". How does this work, because PartType0 contains some fields, like StarFormationRate, not present in gas, and which are not present in the derived field list for the resultant filtered field, yet they do make an obvious difference to the emission and evidently do need to be in the filter?

I attach a MWE below (with the typo corrected), although it is very similar to the cookbook;

import yt
import pyxsim

cuts_dict_arr = [{'field':("gas", "density"), 'less_than':simba_rho_th},
                {'field':("gas", "temperature"), 'gtr_than': 2e5*unyt.K},
                {'field':("PartType0", "density"), 'less_than':simba_rho_th},
                {'field':("PartType0", "temperature"), 'gtr_than': 2e5*unyt.K},
                {'field':("PartType0", "DelayTime"), 'equals':0},
                {'field':("PartType0", "StarFormationRate"), 'equals':0}]

ds = yt.load( "path/to/box", default_species_fields="ionized")

def _filtered_gas(pfilter, data):
    pfilter = True
    for i, cut in enumerate(cuts_dict_arr):
        if cut["equals"] != None:
            pfilter &= data[cut["field"]] == cut["equals"]
        else:
            if cut["gtr_than"] != None:
                pfilter &= data[cut["field"]] > cut["gtr_than"]
            if cut["less_than"] != None:
                pfilter &= data[cut["field"]] < cut["less_than"]        
    return pfilter

yt.add_particle_filter("filtered_gas", function=_filtered_gas, filtered_type="gas", requires=[x["field"][1] for x in cuts_dict_arr if  x["field"][0] == "gas" ])
ds.add_particle_filter("filtered_gas")  

[jzuhone]

Hi @FredJenningsAstro,

I am a little confused around the particle filter though; firstly, passing an empty list to requires didn't stop the filtered field seemingly being set-up correctly, and when I checked, say, the density and temperature extrema of the filtered field they were as expected. This is why I didn't spot this indexing typo first time around, and didn't suspect the filter.

I'm not sure about this one--I'll have to investigate. It is rather odd that this occurred.

Also, like in the cookbook example, I filter both gas and PartType0 in the filter, with the filtered_type = "gas". How does this work, because PartType0 contains some fields, like StarFormationRate, not present in gas, and which are not present in the derived field list for the resultant filtered field, yet they do make an obvious difference to the emission and evidently do need to be in the filter?

This one is easier to explain. "gas" is simply an alias for "PartType0". This does not mean that every "PartType0" field gets mapped to "gas", but it works because they are referring to the same set of particles. We should probably create a mapping from "PartType0", "StarFormationRate" to "gas","star_formation_rate" for consistency, however.

[FredJenningsAstro]

Okay, thanks so much for clearing that up!