ratt-ru / owlcat

miscellaneous utility scripts for manipulating radio interferometry data
4 stars 5 forks source link

add support for zeroing specified region of FITS image #64

Closed o-smirnov closed 3 years ago

o-smirnov commented 4 years ago

Propose --null-region NXxNY[@X0,Y0].

IanHeywood commented 4 years ago

If this is for what I think it's for then reading in positions from a LSM (and possibly selecting by tags) and zeroing the image according to those might be a friendly addition.

o-smirnov commented 4 years ago

I can see why you'd want that. :) Hmmm, I wonder if a separate script in the tigger-lsm package would not be more appropriate (owlcat has no Tigger dependency at present, so I hesitate to add it for just one feature...)

o-smirnov commented 4 years ago

My motivation in this particular case is to subtract peripheral sources as best as possible, so that people with weird imaging algorithms can apply them to the complex source in the middle of the field, while oversampling the f*ck out of the PSF.

If I ever get time, this might also allow me to resuscitate the crazy wide-band image project. Remember the showstopper there? Can't image an L-band FoV at X-band resolution...

bennahugo commented 4 years ago

@IanHeywood I've been down this particular road before and lsms do not model off axis sources correctly for wide bandwidths (unless we can plug in a polynomial fit similar to wsclean's catalog format). In fact I already wrote this functionality into catdagger to blank pixels of tagged lsms

--remove-tagged-dE-components-from-model-images

The only way I ever got de calibration working is through degridding. In this case all the functionality to subtract tagged sources is already part of cubical. You just need to specify MODEL_COLUMN+-my.DicoModel@tags.reg and save it into another model column? There is an option to not write corrected data / flags back, so you can use this to simulate

IanHeywood commented 4 years ago

Here's MeerKAT-16 data (full bandwidth) with a partial LSM model:

dE_LADUMA16

Here's VLA data (1-2 GHz) with a partial LSM model for the shown sources (L-R: 1GC, 2GC, 3GC):

Screen Shot 2020-01-30 at 09 47 52

It might not work for MeerKAT-64 data (although I see no reason why it shouldn't and would love to know why it didn't), but I don't agree that the technique does not work for wide bandwidths.

bennahugo commented 4 years ago

The last time I tried this with MK 16 data (1300 and TriA) off axis sources it actually worsened the image for all choices of intervals. Each time bdsf fitted completely unphysical spectra off intrinsic cubes and the calibration was bound by the modelling errors. If you have a better source finder I would love to know what it is?

When I in turn calibrate with the ddfacet model it actually works, so the problem points at the bdsf source finding and spectrum fitting code).

IanHeywood commented 4 years ago

Ah there's the difference. I make per-sub-band LSMs (4-8 of them) and calibrate the sub-bands (G+dE) independently, and all in the "apparent" domain, assuming nothing about the primary beam. The LSM models are spectrally flat. I've never used the PyBDSF spectral measurement functionality. Fitted in-band spectral indices (for VLA L-band and MeerKAT) tend to be a bit ropey unless the source is very bright. Fitted curvature measurements I would never trust for anything other than things that were calibrator-bright. I think even your 1934 model suffers from this? The two close-in FR-IIs have wild spectral indices for what is surely just synchrotron emission.

bennahugo commented 4 years ago

But this would surely leave jump discontinuities in the data? I prefer to use clean component models and set the number of degridding bands much higher to mitigate this problem.

bennahugo commented 4 years ago

Yes the models are really meant for phaseup @IanHeywood. I'm busy working on deep WSClean polynomial models for both fields in L and UHF bands. I will send these when I'm happy with the dynamic range on the maps

IanHeywood commented 4 years ago

But this would surely leave jump discontinuities in the data?

Undoubtedly, but if it works and doesn't bias unmodelled sources then is that an issue?

If you've got an interpolative de-gridder then most of the sky model will be smooth (although I thought the reason crystalball was written was because people weren't happy with that aspect of wsclean). And I guess you'll end up with discontinuities at the time/frequency boundaries of the solution intervals no matter how the model looks, unless you enforce smoothness on them somehow.

bennahugo commented 4 years ago

No the degridder is not quite interpolative, but because I'm predicting small patches at a time I can set the degridding bands to however many I need without blowing the memory budget. It is true that the solution intervals will do this as well, but typcally one would not have only 4 or so bands to cover the entire bandwidth?

bennahugo commented 4 years ago

What I do think would be great to have is online prediction of wsclean models within cubical?

IanHeywood commented 4 years ago

For killMS I'm solving in 4 bands, but that was to dodge the crash that I was getting when using 8 and it tried to solve for a completely flagged chunk. This is having set up model with DDFacet with --Freq-NDegridBand=8.

For DI-phasecal with CubiCal I was using 8 solution intervals. The model at this stage is typically from a wsclean run with -channels-out 8, but since BDA is enabled I have to follow up with a predictions, which then uses -predict-channels 64. But since CubiCal stopped working on IDIA I reverted to CASA gaincal which of course uses 1 interval for the entire band. It gives very similar results to CubiCal at this stage though.