File: pywsclean.py

package info (click to toggle)
wsclean 3.6-3
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid
  • size: 16,296 kB
  • sloc: cpp: 129,246; python: 22,066; sh: 360; ansic: 230; makefile: 185
file content (276 lines) | stat: -rw-r--r-- 9,024 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
"""
A Python interface to WSClean
This wrapper can be used to call the (C++) WSClean imager.
"""

import _wsclean
import numpy


class ImagingParameters(object):
    """Parameters for imaging"""

    msPath = ""
    imageWidth = 512
    imageHeight = 512
    pixelScaleX = "1asec"
    pixelScaleY = "1asec"
    doNormalize = 1
    extraParameters = ""


class ImagingData(object):
    """Information about the imaging operation"""

    dataSize = 0


class Operator(object):
    """Class that wraps WSClean as an operator, so that it is easy
    to get an image from data 'in memory' (and the inverse). Currently, the
    operator will write that data to the MODEL_DATA Measurement Set before imaging.
    The read/write/backward/forward methods can only be used within a "with" context.
    """

    _userdata = None
    _parameters = None
    _imagingdata = None

    def __init__(self, parameters):
        """Constructor: only sets parameters"""
        self._parameters = parameters

    def __enter__(self):
        """Context manager entrance: initialize WSClean"""
        self._userdata, self._imagingdata = _wsclean.initialize(
            self._parameters
        )
        return self

    def __exit__(self, type, value, traceback):
        """Destructor: release WSClean resources"""
        if self._userdata != None:
            print("Releasing resources for WSClean...")
            _wsclean.deinitialize(self._userdata)
            self._userdata = None
            self._imagingdata = None

    def data_size(self):
        """Get the number of visibilities"""
        if self._userdata == None:
            raise RuntimeError(
                'Operator.data_size() was called outside "with" block'
            )
        return self._imagingdata.dataSize

    def image_size(self):
        """Get the number of pixels in the image"""
        if self._userdata == None:
            raise RuntimeError(
                'Operator.image_size() was called outside "with" block'
            )
        return self._parameters.imageWidth * self._parameters.imageHeight

    def read(self):
        """Read the visibilities and return as a (data,weight) tuple."""
        if self._userdata == None:
            raise RuntimeError(
                'Operator.read() was called outside "with" block'
            )
        data = numpy.ascontiguousarray(
            numpy.zeros(self._imagingdata.dataSize, dtype=numpy.complex128)
        )
        weights = numpy.ascontiguousarray(
            numpy.zeros(self._imagingdata.dataSize, dtype=numpy.float64)
        )
        _wsclean.read(self._userdata, data, weights)
        return data, weights

    def write(self, filename, data):
        """Write a FITS image with the correct keywords etc."""
        if self._userdata == None:
            raise RuntimeError(
                'Operator.write() was called outside "with" block'
            )
        dataCont = numpy.ascontiguousarray(data)
        _wsclean.write(self._userdata, filename, data)

    def forward(self, dataOut, dataIn):
        """Perform the forward operation. This is 'prediction': convert
        an image into visibilities. dataOut should be an complex double array
        that will be filled with visibilities, dataIn should be an array
        of doubles, representing the image for the operator input."""

        if self._userdata == None:
            raise RuntimeError(
                'Operator.forward() was called outside "with" block'
            )

        if numpy.shape(dataOut)[0] != self.data_size():
            raise RuntimeError(
                "Size of output argument ("
                + str(numpy.shape(dataOut)[0])
                + ") does not match the image size ("
                + str(self.data_size())
                + ")"
            )

        if numpy.shape(dataIn)[0] != self.image_size():
            raise RuntimeError(
                "Shape of input argument ("
                + str(numpy.shape(dataIn)[0])
                + ") does not match the number of visibilities ("
                + str(self.image_size())
                + ")"
            )

        if dataOut.dtype.name != "complex128":
            raise RuntimeError(
                "The dataOut parameter of forward() should be of type complex128, but was "
                + dataOut.dtype.name
            )

        if dataIn.dtype.name != "float64":
            raise RuntimeError(
                "The dataIn parameter of forward() should be of type float64, but was "
                + dataIn.dtype.name
            )

        dataOutCont = numpy.ascontiguousarray(dataOut)
        dataInCont = numpy.ascontiguousarray(dataIn)
        _wsclean.operator_A(self._userdata, dataOut, dataIn)

    def backward(self, dataOut, dataIn):
        """Perform the backward operation. This is the 'imaging' step:
        convert visibilities into an image. dataOut should be an array
        of doubles, which will be filled with the image, dataOut should be an array
        of complex doubles, representing the visibilities for the operator input.
        """

        if self._userdata == None:
            raise RuntimeError(
                'Operator.backward() was called outside "with" block'
            )

        if numpy.shape(dataOut)[0] != self.image_size():
            raise RuntimeError(
                "Size of output argument ("
                + str(numpy.shape(dataOut)[0])
                + ") does not match the image size ("
                + str(self.image_size())
                + ")"
            )

        if numpy.shape(dataIn)[0] != self.data_size():
            raise RuntimeError(
                "Shape of input argument ("
                + str(numpy.shape(dataIn)[0])
                + ") does not match the number of visibilities ("
                + str(self.data_size())
                + ")"
            )

        if dataOut.dtype.name != "float64":
            raise RuntimeError(
                "The dataOut parameter of forward() should be of type float64, but was "
                + dataOut.dtype.name
            )

        if dataIn.dtype.name != "complex128":
            raise RuntimeError(
                "The dataIn parameter of forward() should be of type complex128, but was "
                + dataIn.dtype.name
            )

        dataOutCont = numpy.ascontiguousarray(dataOut)
        dataInCont = numpy.ascontiguousarray(dataIn)
        _wsclean.operator_At(self._userdata, dataOut, dataIn)


class WSClean(object):
    """The Python interface to WSClean"""

    datacolumn = ""
    """The column used for imaging; empty means CORRECTED_DATA if it exists, otherwise
	use DATA."""

    width = 1024
    """Image width"""

    height = 1024
    """Image height"""

    scale = "1asec"
    """Pixel scale of image. Units can e.g. be deg, amin, asec, masec. There should not
	be a space between the number and its unit."""

    niter = 0
    """Number of clean or moresane iterations"""

    gain = -1
    """Gain per minor iteration. -1 means use WSClean's default."""

    mgain = -1
    """Gain per major iteration. -1 means use WSClean's default."""

    __weightpar = ""

    def __init__(self):
        return

    def image(self, msnames, nameprefix):
        """Run WSClean to make an image on the specified list of measurement sets"""
        plist = self.__get_parameterlist(nameprefix)
        import os

        msnamelist = " ".join(msnames)
        cmd = "wsclean " + str(plist) + " " + msnamelist
        print(cmd)
        os.system(cmd)
        return

    def predict(self, msnames, nameprefix):
        """Run WSClean to predict"""
        plist = self.__get_parameterlist(nameprefix)
        import os

        msnamelist = " ".join(msnames)
        cmd = "wsclean -predict " + str(plist) + " " + msnamelist
        print(cmd)
        os.system(cmd)
        return

    def set_uniform_weighting(self):
        """Enable uniform weighting"""
        self.__weightpar = "-weight uniform"

    def set_natural_weighting(self):
        """Enable natural weighting"""
        self.__weightpar = "-weight natural"

    def set_briggs_weighting(self, robustness):
        """Enable Briggs' weighting with a given robustness"""
        self.__weightpar = "-weight briggs " + str(robustness)

    def __get_parameterlist(self, prefixname):
        plist = (
            "-size "
            + str(self.width)
            + " "
            + str(self.height)
            + " -scale "
            + str(self.scale)
        )
        if self.datacolumn != "":
            plist += " -data-column " + self.datacolumn
        if self.__weightpar != "":
            plist += " " + self.__weightpar
        if self.niter != 0:
            plist += " -niter " + str(self.niter)
        if self.gain != -1:
            plist += " -gain " + str(self.gain)
        if self.mgain != -1:
            plist += " -mgain " + str(self.mgain)
        if prefixname != "":
            plist += " -name " + prefixname
        return plist