ctools

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":https://en.wikipedia.org/wiki/Total_variation_denoising 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":https://github.com/scikit-image/scikit-image/blob/master/skimage/restoration/_denoise.py .

You can see a examples of the algorithm "here":http://www.tp.umu.se/~nylen/fnm/pylect/advanced/image_processing/auto_examples/plot_lena_tv_denoise.html#example-plot-lena-tv-denoise-py and "here":http://scikit-image.org/docs/dev/auto_examples/plot_denoise.html .

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.