ratt-ru / pfb-imaging

Preconditioned forward/backward clean algorithm
MIT License
7 stars 5 forks source link

Trouble using ms with multiple spws #82

Closed y-mhiri closed 1 year ago

y-mhiri commented 1 year ago

The grid worker results with the following error when trying to image the multi spws measurement set from the VLA Casa Imaging guide (http://casa.nrao.edu/Data/EVLA/SNRG55/SNR_G55_10s.calib.tar.gz) :

Traceback (most recent call last): File "/home/vagrant/.local/bin/pfb", line 33, in sys.exit(load_entry_point('pfb-clean', 'console_scripts', 'pfb')()) File "/home/vagrant/.local/lib/python3.8/site-packages/click/core.py", line 1134, in call return self.main(args, kwargs) File "/home/vagrant/.local/lib/python3.8/site-packages/click/core.py", line 1059, in main rv = self.invoke(ctx) File "/home/vagrant/.local/lib/python3.8/site-packages/click/core.py", line 1665, in invoke return _process_result(sub_ctx.command.invoke(sub_ctx)) File "/home/vagrant/.local/lib/python3.8/site-packages/click/core.py", line 1401, in invoke return ctx.invoke(self.callback, ctx.params) File "/home/vagrant/.local/lib/python3.8/site-packages/click/core.py", line 767, in invoke return __callback(args, kwargs) File "/vagrant/pfb-clean/pfb/workers/grid.py", line 50, in grid return _grid(opts) File "/vagrant/pfb-clean/pfb/workers/grid.py", line 135, in _grid xdso = xr.concat(xdsb, dim='row', File "/home/vagrant/.local/lib/python3.8/site-packages/xarray/core/concat.py", line 243, in concat return _dataset_concat( File "/home/vagrant/.local/lib/python3.8/site-packages/xarray/core/concat.py", line 466, in _dataset_concat align(*datasets, join=join, copy=False, exclude=[dim], fill_value=fill_value) File "/home/vagrant/.local/lib/python3.8/site-packages/xarray/core/alignment.py", line 772, in align aligner.align() File "/home/vagrant/.local/lib/python3.8/site-packages/xarray/core/alignment.py", line 565, in align self.reindex_all() File "/home/vagrant/.local/lib/python3.8/site-packages/xarray/core/alignment.py", line 542, in reindex_all self.results = tuple( File "/home/vagrant/.local/lib/python3.8/site-packages/xarray/core/alignment.py", line 543, in self._reindex_one(obj, matching_indexes) File "/home/vagrant/.local/lib/python3.8/site-packages/xarray/core/alignment.py", line 527, in _reindex_one dim_pos_indexers = self._get_dim_pos_indexers(matching_indexes) File "/home/vagrant/.local/lib/python3.8/site-packages/xarray/core/alignment.py", line 493, in _get_dim_pos_indexers indexers = obj_idx.reindex_like(aligned_idx, **self.reindex_kwargs) File "/home/vagrant/.local/lib/python3.8/site-packages/xarray/core/indexes.py", line 514, in reindex_like raise ValueError( ValueError: cannot reindex or align along dimension 'chan' because the (pandas) index has duplicate values

The error is raised when the datasets of the xds are concatenated to produce nband datasets. The error comes from the fact that datasets from multiple scans with the same frequency channels are concatenated.

landmanbester commented 1 year ago

Thanks for reporting @y-mhiri . Indeed, overlapping SPW's will currently cause some pain. I am looking into this, I think there are two avenues that can solve the problem

i) Don't --concat by default. This has the drawback that we have to grid data mapping to the same imaging band independently which means duplicating FFT's but it should work as it stands. You can combine the dirty images after the fact.

ii) Take the weighted sum of the data where they overlap on the time frequency grid. This happens on what is effectively "corrected data" with the correct weights so it should work. But I would need to write a custom concatenation routine for that.

Can you try the grid worker with --no-concat and let me know if you run into trouble please? You would have to combine the dirty images and PSF's but I think this is what happens by default if you use the clean or spotless deconvolution routines. Let's keep this open as it really needs a more elegant solution...

y-mhiri commented 1 year ago

I tried the grid worker with --no-concat and got the following traceback : Traceback (most recent call last): File "/home/vagrant/.local/bin/pfb", line 33, in sys.exit(load_entry_point('pfb-clean', 'console_scripts', 'pfb')()) File "/home/vagrant/.local/lib/python3.8/site-packages/click/core.py", line 1134, in call return self.main(args, kwargs) File "/home/vagrant/.local/lib/python3.8/site-packages/click/core.py", line 1059, in main rv = self.invoke(ctx) File "/home/vagrant/.local/lib/python3.8/site-packages/click/core.py", line 1665, in invoke return _process_result(sub_ctx.command.invoke(sub_ctx)) File "/home/vagrant/.local/lib/python3.8/site-packages/click/core.py", line 1401, in invoke return ctx.invoke(self.callback, ctx.params) File "/home/vagrant/.local/lib/python3.8/site-packages/click/core.py", line 767, in invoke return __callback(args, kwargs) File "/vagrant/pfb-clean/pfb/workers/grid.py", line 50, in grid return _grid(opts) File "/vagrant/pfb-clean/pfb/workers/grid.py", line 446, in _grid fitsout.append(dds2fits_mfs(dds, 'DIRTY', basename, norm_wsum=True)) File "/vagrant/pfb-clean/pfb/utils/fits.py", line 191, in dds2fits_mfs datas[t] += ds.get(column).data IndexError: list index out of range

I am looking into it in the next few days. Also when using --no-concat the dask.compute() function is quite slow to run. Do you have an idea on why ?

landmanbester commented 1 year ago

Ah that error has been fixed in the solarkat branch. I have merged spwindexing into solarkat so you may have some luck with the solarkat branch. I can't currently test it though as I have very limited access to internet. Give me until Monday then I can run my tests.

I suspect the reason it is slower is because it will now do the gridding on a per scan and spectral window basis, increasing the number of required FFT's. That is really why we need a better solution for this problem

landmanbester commented 1 year ago

There was indeed a bug there that I had to hunt down. I have now changed how the mapping logic works. The init worker now takes a channels_per_image argument which controls how many datasets a spectral window is split into (note there is no longer an nband argument, this is set during gridding). By default this is set to -1 which means make one dataset per spectral window. It should now be possible to coarsen the spectral resolution to an arbitrary number of imaging bands during gridding. The grid worker will try to do a weighted sum of the "corrected data" where there is an overlap. It should now be possible to combine arbitrary overlapping spectral windows as long as all the spectral windows have the same shape. I am working on making this more general but @y-mhiri please have a go with the solarkat branch and let me know if you run into issues

y-mhiri commented 1 year ago

Hi Landman, Many thanks for your work ! I'm trying that this afternoon and let you know what i get

landmanbester commented 1 year ago

I should add that you can now do --concat-row (the default) and it should work. I will eventually add similar functionality for arbitrarily coarsening the time axis during gridding but that is a WIP. For now you either stick to the time resolution specified at the init stage or you combine all of them into a single output image

y-mhiri commented 1 year ago

So I went with the solarkat branch and got the right dirty image with the grid worker using the default parameters, so it seems to work perfectly !

landmanbester commented 1 year ago

Good to know. Note there was a small typo which may have resulted in inconsistent flags. Fixed in the solarkat branch now. Please reopen this one if you run into further issues