RFE: fluximage create psfmaps

[kglotfelty]

We've talked about this before but don't think it ever made it to "the" list. I'm doing celldetect which now needs the psfmap ... so having them created automatically in fluximage would be helpful.

To run mkpsfmap we need the energy -- same mono energy used to make expmaps, and the ECF -- so that would probably need to be a new hidden fluximage parameter. Default probably the 1sigma, 0.393 value -- though I now see that celldetect has been using 0.8 with the old psftable cal files.

[DougBurke]

Yes, it's been on my back burner so long I've lost site of it...

If we use a spectrum weights file for the bands parameter then I think it would make sense just to skip the psf map creation. I don't think we want to try and "guess" the energy from this (eg we could take the position of the peak value in the spectrum file but I think at this point we have gone way beyond the use case we want for this script).

I assume the thought is to always create the psfmap (modulo the weight issue above), rather than to being optional?

Doug

[DougBurke]

Also, would it make sense to combine psfmaps in merge_obs/flux_obs - presumably providing a choice of methods based on your memo/thread/whatever?

[kglotfelty]

I would think optional -- since it can be slow for large images (eg HRC)

mkpsfmap does allow for a spectrum -- not sure it'll be the same format as mkinstmap -- that'd be too easy. [I just don't remember where we landed with that.]

The merge is tricky . Yes, maybe options, but probably the "ugly" one is the best (exposure weighted min being careful to only consider pixels in the per-obi fov).

[DougBurke]

Should the PSF map be filtered by the SKY FOV file (that is, to make sure the outside is 0/not in the filtered subspace)? I think so, since we want to be able - in the merging case - to only "combine" pixels that are actually on chip. This filtering could be done at the merging stage, but I think I'd like to do it earlier.

EDITED TO ADD an alternative/corollary - in fluximage we allow the data to be thresholded; should this be applied to the PSF map too (I would think so for consistency)

[kglotfelty]

That works -- that is I can't think of a reason not to.

I tested celldetect w/ both 0's and NaN's outside the FOV. I know it's OK w/ wavdetect since that's how L3 does it.

The only thing is just make sure to use [opt full] to keep the same dimensions.

[DougBurke]

See #169 for what I've got so far.

[kglotfelty]

Quick test looks mostly good. I tried with an obs with 4 chips and got

% fluximage 1522/primary/\*evt2.fits.gz"[ccd_id=0:3]" \
    out=trymultichippsfmap xygrid=3072.5:5120.5:1,3072.5:5120.5:1 \
    mode=h clob+ band=broad psfecf=0.8 
...
Creating 4 PSF maps for obsid 1522

I'm only seeing 1 psfmap -- just not sure if it tried to make the same thing 4 times? The final file looks good.

I also tried with a single chip and in the case the psfmap did already include the FOV boundary

% fluximage \
  "578/primary/acisf00578N003_evt2.fits[ccd_id=7,sky=region(578/primary/acisf00578_000N003_fov1.fits[ccd_id=7])]"\
    binsize=1 bands=broad outroot=trypsfmap psfecf=0.8 clob+

% dmstat trypsfmap_broad.psfmap cen- 
PSFMAP[arcsec]
    min:    0.62097513676         @:    ( 4109 4124 )
    max:    7.0572948456          @:    ( 3920 5065 )
   mean:    2.4446281347 
  sigma:    1.4849433555 
    sum:    2830952.7188 
   good:    1158030 
   null:    772680 

Ahh ... it probably has the boundary because the event file I supplied had the fov boundary in it.

[DougBurke]

The screen message is definitely wrong. I'd originally thought that I had to do this per-chip then realised that was wrong and it just needed to be done per band, but obviously didn't get the counting right in the label (it should have only been run once on this case).

[DougBurke]

I've fixed the screen message. I have also implemented an initial attempt at filtering the output. I have not tried it in combination with a pre-filtered event file - e.g. the second example in https://github.com/cxcsds/ciao-contrib/issues/164#issuecomment-440273817 - which could be interesting (since the masked PSF map has 0 for "outside" whereas the mkpsfmap output for "outside" the subspace of the input image is NaN).

[kglotfelty]

I'm not saying to do this ... just a couple of notes ... [I haven't looked at latest changes so feel free to skip if irrelevant.]

Instead of filtering the output, you could filter the input to mkpsfmap and it'll set the pixels outside of the subspace=NaN.

As for say applying threshold via dmimgthresh , you can set val=INDEF to set output pixels to NaN instead of 0 (I think you may also be able to just use nan too). And|or if you want to make the NaN's be 0, then use dmimgthresh with cut=INDEF val=0

[DougBurke]

When I get home I shall have to play with

% mkspfmap "broad.img[sky=mask(broad.expmap)]" broad.psf ...

as that would save an extra step / inventing a temporary file name.

[kglotfelty]

[DougBurke]

I've just pushed a commit which adds basic support to flux_obs - the psfmerge parameter.

The implementation is a bit up in the air (eg expmap weighting is not supproted although advertised), but I'd be interested in whether it works. Also, the exposure-time weighting is a bit different from the thread (as it now accounts for the total observation time per pixel, or at least it should).

[DougBurke]

The fluximage command has the psfecf parameter, and the flux_obs and merge_obs have this and the psfmerge parameter to control how the PSF maps can be combined:

% cat `paccess flux_obs` | grep psf
psfecf,r,h,INDEF,0,1,"If set, create PSF map with this ECF"
psfmerge,s,h,"min",exptime|expmap|min|max|mean|median|mid,,"How are the PSF maps combined?"

There's some "fun" behind the scenes for the exptime and expmap options, in particular in emulating the DataModel's merging rules.

Also, I do use a mask to filter the per-obsid PSF maps, but then immediately remove this MASK from the subspace (since mkpsfmap sets invalid pixels to NaN we do not lose this filter, and having a "clean" SKY subspace makes merging a bit easier).

[DougBurke]

One thing I'm unhappy with the current code - e.g. as of https://github.com/cxcsds/ciao-contrib/pull/169/commits/7562b1b6dcd5694dbb450ca994d8ab77e6858cda - is that the file naming doesn't follow the other "thresholded" files. Since we use the thresholded exposure map to mask the PSF map, it probably should.

This is not a large or hard change to make, just a bit messy.

[DougBurke]

Fixed in #169