Feature #1530
Python Function to Convert GObservations to GCTABackground3D
Status: | Feedback | Start date: | |||
---|---|---|---|---|---|
Priority: | Normal | Due date: | |||
Assigned To: | - | % Done: | 90% | ||
Category: | - | ||||
Target version: | - | ||||
Duration: |
Description
I’ve been working on veritas backgrounds for a bit, and ended up developing this function (attached python script, function write_background_to_fits()
) that creates a background fits file for GCTABackground3D from a GObservations object. The function can optionally apply variable-width gaussian smoothing to each bin.
While trying to construct backgrounds with a limited number of events, I found a decent denoising algorithm that acts like a low-band-pass filter, smearing out large single-bin variations, but leaving multi-bin variations intact.
The denoising is beneficial, since applying a gaussian smoothing tends to alter the shape of the background, and when you make a background from a limited number of events, statistical fluctuations tend to make the background worse (i.e. less representative of the actual background shape). So, instead of gaussian smoothing, the attached script can instead do 'totalvariation’ denoising, through the denoise_tv_chambolle_2d()
function I “borrowed” from the skimage python module .
You can see a examples of the algorithm here and here .
The python script is attached, along with two examples where large relative fluctuations in the raw events histogram (the left plots) are smoothed out with the denoise_tv_chambolle_2d(...,weight=10) function (the right plots), while normalising for the total number of events.
Note that the python script does require numpy, scipy, and astropy, but these are easy to get through pip. Running the script should explain how to use it.
For this feature: I’d like the attached python script to be turned a cscript, but I’m not sure I could easily figure out how to write a fits file using GFits structures, or how to deal with the numpy-dependent parts, since I think we’re trying to avoid relying on external libraries? It might also be nice to add more advanced image filtering, like full fourier cleaning to get rid of high-frequency noise. So, I’m just attaching it here for now.
Recurrence
No recurrence.
Related issues
History
#1 Updated by Kelley-Hoskins Nathan about 9 years ago
- File bck_0.46TeV.png added
- Description updated (diff)
attached a second example denoising plot.
#2 Updated by Kelley-Hoskins Nathan about 9 years ago
- File obs2background.py added
and of course I forget the reason for the post.
#3 Updated by Deil Christoph about 9 years ago
Jürgen, I’d be very interested to hear what the plan is for such fancy algorithms (image processing in this case) that are non-trivial to re-implement.
Will Gammalib / ctools stick with the no external dependency approach or will dependencies like Scipy or scikit-image or Astropy be allowed for certain functionality like background model production here?
#4 Updated by Mayer Michael almost 9 years ago
- Start date deleted (
10/02/2015)
I currently see two options in implementing such a function.
1. Write a cscript
that depends on an external package (no unit testing possible then, I guess). I basically would have such a script ready for this. If we want to go this way I can commit it.
2. Implement some smoothing algorithms in gammalib
for this and write a cscript
or ctool
that uses them. Judging from the python code Nathan provided it doesn’t seem too hard to translate it into C++. In particular, we have GMatrix
classes, that mimc somewhat the numpy functionality.
For the latter one, we could then simply add nice functions to GSkyMap
, e.g. GSkyMap::smooth_gauss()
, GSkyMap::tv_denoise()
or similar. However, I have no idea about the license, are we allowed to simply copy and translate such things to C++?
#5 Updated by Mayer Michael almost 9 years ago
- Status changed from New to Feedback
- % Done changed from 0 to 90
I guess the need for such a script is quite high. So my proposal for the moment is to make it an example script.
I have committed a cscript
-like script csbkgmodel.py
to the examples folder on branch 1530-csbkgmodel. Jürgen, do you think we could let this life in the examples folder until we have come up with a better solution?
#6 Updated by Knödlseder Jürgen almost 9 years ago
- Related to Feature #1729: Add support to smooth sky maps added