Coverage for /wheeldirectory/casa-6.7.0-12-py3.10.el8/lib/py/lib/python3.10/site-packages/casatasks/private/task_sdfixscan.py: 6%

411 statements  

« prev     ^ index     » next       coverage.py v7.6.4, created at 2024-10-31 19:53 +0000

1# sd task for image processing (fft_mask or model) 

2import os 

3import time 

4 

5import numpy 

6import numpy.fft as npfft 

7 

8from casatasks import casalog 

9from casatools import ctsys 

10from casatools import image as iatool 

11from casatools import quanta 

12 

13from . import sdutil 

14 

15 

16def create_4d_image(infile, outfile): 

17 ia = iatool() 

18 ia.open(infile) 

19 image_shape = ia.shape() 

20 try: 

21 if len(image_shape) < 4: 

22 # add degenerate axes 

23 

24 cs = ia.coordsys() 

25 axistypes = cs.axiscoordinatetypes() 

26 no_stokes = 'Stokes' not in axistypes 

27 no_spectral = 'Spectral' not in axistypes 

28 stokes_axis = 'I' if no_stokes else '' 

29 outimage = ia.adddegaxes(outfile=outfile, spectral=no_spectral, 

30 stokes=stokes_axis) 

31 else: 

32 # generage complete copy of input image using subimage 

33 outimage = ia.subimage(outfile=outfile) 

34 finally: 

35 if len(image_shape) < 4: 

36 cs.done() 

37 ia.close() 

38 

39 return outimage 

40 

41 

42@sdutil.sdtask_decorator 

43def sdfixscan(infiles, mode, numpoly, beamsize, smoothsize, direction, maskwidth, 

44 tmax, tmin, outfile, overwrite): 

45 with sdfixscan_worker(**locals()) as worker: 

46 worker.initialize() 

47 worker.execute() 

48 worker.finalize() 

49 

50 

51class sdfixscan_worker(sdutil.sdtask_interface): 

52 def __init__(self, **kwargs): 

53 super(sdfixscan_worker, self).__init__(**kwargs) 

54 

55 def __del__(self, base=sdutil.sdtask_interface): 

56 # cleanup method must be called when the instance is 

57 # deleted 

58 self.cleanup() 

59 super(sdfixscan_worker, self).__del__() 

60 

61 def initialize(self): 

62 self.parameter_check() 

63 

64 # temporary filename 

65 tmpstr = time.ctime().replace(' ', '_').replace(':', '_') 

66 self.tmpmskname = 'masked.' + tmpstr + '.im' 

67 self.tmpconvname = 'convolve2d.' + tmpstr + '.im' 

68 self.tmppolyname = 'polyfit.' + tmpstr + '.im' 

69 # set tempolary filename 

70 self.tmprealname = [] 

71 self.tmpimagname = [] 

72 self.image = None 

73 self.convimage = None 

74 self.polyimage = None 

75 self.imageorg = None 

76 self.realimage = None 

77 self.imagimage = None 

78 if type(self.infiles) == str: 

79 self.tmprealname.append('fft.' + tmpstr + '.real..0.im') 

80 self.tmpimagname.append('fft.' + tmpstr + '.imag.0.im') 

81 else: 

82 for i in range(len(self.infiles)): 

83 self.tmprealname.append('fft.%s.%s.real.im' % (tmpstr, i)) 

84 self.tmpimagname.append('fft.%s.%s.imag.im' % (tmpstr, i)) 

85 

86 # default output filename 

87 if self.outfile == '': 

88 self.outfile = 'sdfixscan.out.im' 

89 casalog.post('outfile=%s' % self.outfile) 

90 

91 # threshold 

92 self.nolimit = 'nolimit' 

93 if self.tmin == 0.0 and self.tmax == 0.0: 

94 self.thresh = [] 

95 elif self.tmin > self.tmax: 

96 casalog.post('tmin > tmax. Swapped.') 

97 self.thresh = [self.tmax, self.tmin] 

98 elif self.tmin == self.tmax: 

99 if self.tmin > 0.0: 

100 casalog.post('tmin == tmax. Use tmin as minumum threshold.') 

101 self.thresh = [self.tmin, self.nolimit] 

102 else: 

103 casalog.post('tmin == tmax. Use tmax as maximum threshold.') 

104 self.thresh = [self.nolimit, self.tmin] 

105 else: 

106 self.thresh = [self.tmin, self.tmax] 

107 

108 def parameter_check(self): 

109 if self.mode.lower() == 'model': 

110 # Pressed-out method 

111 # check input file 

112 if type(self.infiles) == list: 

113 if len(self.infiles) != 1: 

114 raise Exception("infiles allows only one input file for pressed-out method.") 

115 else: 

116 self.infiles = self.infiles[0] 

117 # check direction 

118 if type(self.direction) == list: 

119 if len(self.direction) != 1: 

120 raise Exception("direction allows only one direction for pressed-out method.") 

121 else: 

122 self.direction = self.direction[0] 

123 elif self.mode.lower() == 'fft_mask': 

124 # FFT-based basket-weaving method 

125 # check input file 

126 if type(self.infiles) == str or \ 

127 (type(self.infiles) == list and len(self.infiles) < 2): 

128 raise Exception("infiles should be a list of input images for Basket-Weaving.") 

129 

130 # check direction 

131 if type(self.direction) == float: 

132 raise Exception('direction must have at least two different direction.') 

133 else: 

134 if len(self.direction) < 2: 

135 raise Exception('direction must have at least two different direction.') 

136 else: 

137 raise Exception('Unsupported processing mode: %s' % (self.mode)) 

138 

139 def execute(self): 

140 if self.mode.lower() == 'model': 

141 self.__execute_press() 

142 elif self.mode.lower() == 'fft_mask': 

143 self.__execute_basket_weaving() 

144 

145 def __execute_press(self): 

146 ### 

147 # Pressed-out method (Sofue & Reich 1979) 

148 ### 

149 casalog.post('Apply Pressed-out method') 

150 

151 # CAS-5410 Use private tools inside task scripts 

152 ia = iatool() 

153 

154 # mask 

155 self.image = ia.newimagefromimage(infile=self.infiles, outfile=self.tmpmskname) 

156 # back-up original mask name 

157 is_initial_mask = (self.image.maskhandler('default')[0] != '') 

158 temp_maskname = "temporal" 

159 # imshape = self.image.shape() 

160 # ndim = len(imshape) 

161 # nx = imshape[0] 

162 # ny = imshape[1] 

163 if len(self.thresh) == 0: 

164 casalog.post('Use whole region') 

165 else: 

166 # mask pixels beyond thresholds 

167 maskstr = ("mask('%s')" % self.tmpmskname) 

168 if self.thresh[0] != self.nolimit: 

169 maskstr += (" && '%s'>=%f" % (self.tmpmskname, self.thresh[0])) 

170 if self.thresh[1] != self.nolimit: 

171 maskstr += (" && '%s'<=%f" % (self.tmpmskname, self.thresh[1])) 

172 # Need to flush to image once to calcmask ... sigh 

173 self.image.done() 

174 self.image = ia.newimage(self.tmpmskname) 

175 self.image.calcmask(mask=maskstr, name=temp_maskname, asdefault=True) 

176 

177 # smoothing 

178 # bmajor = 0.0 

179 # bminor = 0.0 

180 # CAS-5410 Use private tools inside task scripts 

181 qa = quanta() 

182 if type(self.beamsize) == str: 

183 qbeamsize = qa.quantity(self.beamsize) 

184 else: 

185 qbeamsize = qa.quantity(self.beamsize, 'arcsec') 

186 if type(self.smoothsize) == str: 

187 # bmajor = smoothsize 

188 # bminor = smoothsize 

189 qsmoothsize = qa.quantity(self.smoothsize) 

190 else: 

191 # bmajor = '%sarcsec' % (beamsize*smoothsize) 

192 # bminor = '%sarcsec' % (beamsize*smoothsize) 

193 qsmoothsize = qa.mul(qbeamsize, self.smoothsize) 

194 bmajor = qsmoothsize 

195 bminor = qsmoothsize 

196 pa = qa.quantity(0.0, 'deg') 

197 # masked channels are replaced by zero and convolved here. 

198 self.convimage = self.image.convolve2d(outfile=self.tmppolyname, 

199 major=bmajor, minor=bminor, pa=pa, 

200 overwrite=True) 

201 self.convimage.done() 

202 

203 # get dTij (original - smoothed) 

204 self.convimage = ia.imagecalc(outfile=self.tmpconvname, 

205 pixels='"{org}" - "{conv}"'.format(org=self.tmpmskname, 

206 conv=self.tmppolyname), 

207 overwrite=True) 

208 

209 # polynomial fit 

210 fitaxis = 0 

211 if self.direction == 0.0: 

212 fitaxis = 0 

213 elif self.direction == 90.0: 

214 fitaxis = 1 

215 else: 

216 raise Exception( 

217 "Sorry, the task don't support inclined scan " 

218 "with respect to horizontal or vertical axis, right now.") 

219 if os.path.exists(self.tmppolyname): 

220 # CAS-5410 Use private tools inside task scripts 

221 ctsys.removetable([self.tmppolyname]) 

222 self.convimage.setbrightnessunit('K') 

223 # Unfortunately, ia.fitprofile is very fragile. 

224 # Using numpy instead for fitting with masked pixels (KS, 2014/07/02) 

225 # resultdic = self.convimage.fitprofile(model=self.tmppolyname, axis=fitaxis, 

226 # poly=self.numpoly, ngauss=0, multifit=True, gmncomps=0 ) 

227 self.__polynomial_fit_model(image=self.tmpmskname, 

228 model=self.tmppolyname, axis=fitaxis, order=self.numpoly) 

229 polyimage = ia.newimage(self.tmppolyname) 

230 # set back defalut mask (need to get from self.image) 

231 avail_mask = polyimage.maskhandler('get') 

232 if is_initial_mask: 

233 casalog.post("copying mask from %s" % (self.infiles)) 

234 polyimage.calcmask("mask('%s')" % self.infiles, asdefault=True) 

235 else: # no mask in the original image 

236 polyimage.calcmask('T', asdefault=True) 

237 if temp_maskname in avail_mask: 

238 polyimage.maskhandler('delete', name=temp_maskname) 

239 

240 # subtract fitted image from original map 

241 subtracted = ia.imagecalc(outfile=self.outfile, 

242 pixels='"{org}" - "{fit}"'.format(org=self.infiles, 

243 fit=self.tmppolyname), 

244 overwrite=self.overwrite) 

245 subtracted.done() 

246 

247 # finalization 

248 polyimage.done(remove=True) 

249 self.convimage.done(remove=True) 

250 self.image.done() 

251 

252 def __polynomial_fit_model(self, image=None, model=None, axis=0, order=2): 

253 if not image or not os.path.exists(image): 

254 raise RuntimeError("No image found to fit.") 

255 if os.path.exists(model): 

256 # CAS-5410 Use private tools inside task scripts 

257 ctsys.removetable([model]) 

258 tmpia = iatool() 

259 modelimg = tmpia.newimagefromimage(infile=image, outfile=model) 

260 try: 

261 if tmpia.isopen(): 

262 tmpia.close() 

263 imshape = modelimg.shape() 

264 # the axis order of [ra, dec, chan(, pol)] is assumed throughout the task. 

265 ndim = len(imshape) 

266 nx = imshape[0] 

267 ny = imshape[1] 

268 # an xy-plane can be fit simultaneously (doing per plane to save memory) 

269 if ndim == 3: 

270 def get_blc(i, j): 

271 return [0, 0, i] 

272 

273 def get_trc(i, j): 

274 return [nx - 1, ny - 1, i] 

275 

276 imshape2 = imshape[2] 

277 imshape3 = 1 

278 elif ndim == 4: 

279 def get_blc(i, j): 

280 return [0, 0, i, j] 

281 

282 def get_trc(i, j): 

283 return [nx - 1, ny - 1, i, j] 

284 

285 imshape2 = imshape[2] 

286 imshape3 = imshape[3] 

287 else: # ndim == 2 

288 def get_blc(i, j): 

289 return [0, 0] 

290 

291 def get_trc(i, j): 

292 return [nx - 1, ny - 1] 

293 

294 imshape2 = 1 

295 imshape3 = 1 

296 

297 for i3 in range(imshape3): 

298 for i2 in range(imshape2): 

299 blc = get_blc(i2, i3) 

300 trc = get_trc(i2, i3) 

301 dslice = modelimg.getchunk(blc, trc) 

302 mslice = modelimg.getchunk(blc, trc, getmask=True) 

303 model = self._get_polyfit_model_array(dslice.reshape(nx, ny), 

304 mslice.reshape(nx, ny), 

305 axis, order) 

306 modelimg.putchunk(model, blc) 

307 

308 # the fit model image itself is free from invalid pixels 

309 modelimg.calcmask('T', asdefault=True) 

310 except Exception: 

311 raise 

312 finally: 

313 modelimg.close() 

314 

315 def _get_polyfit_model_array(self, data, mask, axis, order): 

316 if axis == 1: 

317 tmp = data.transpose() 

318 data = tmp 

319 tmp = mask.transpose() 

320 mask = tmp 

321 del tmp 

322 nx = data.shape[0] 

323 ny = data.shape[1] 

324 x = range(nx) 

325 flag = mask ^ True # invert mask for masked array 

326 mdata = numpy.ma.masked_array(data, flag) 

327 retc = numpy.ma.polyfit(x, mdata, order) 

328 del flag 

329 coeffs = retc.transpose() 

330 tmpmodel = numpy.zeros([nx, ny]) 

331 for iy in range(ny): 

332 tmpmodel[:, iy] = numpy.poly1d(coeffs[iy])(x) 

333 if axis == 1: 

334 return tmpmodel.transpose() 

335 return tmpmodel 

336 

337 def __execute_basket_weaving(self): 

338 ### 

339 # Basket-Weaving (Emerson & Grave 1988) 

340 ### 

341 casalog.post('Apply Basket-Weaving') 

342 

343 # CAS-5410 Use private tools inside task scripts 

344 ia = iatool() 

345 

346 # initial setup 

347 outimage = ia.newimagefromimage( 

348 infile=self.infiles[0], outfile=self.outfile, overwrite=self.overwrite) 

349 imshape_out = outimage.shape() 

350 # ndim_out = len(imshape_out) 

351 coordsys = outimage.coordsys() 

352 axis_types = coordsys.axiscoordinatetypes() 

353 # direction axis should always exist 

354 try: 

355 direction_axis0 = axis_types.index('Direction') 

356 direction_axis1 = axis_types[direction_axis0 + 1:].index('Direction') \ 

357 + direction_axis0 + 1 

358 except IndexError: 

359 raise RuntimeError('Direction axes don\'t exist.') 

360 finally: 

361 coordsys.done() 

362 nx = imshape_out[direction_axis0] 

363 ny = imshape_out[direction_axis1] 

364 tmp = [] 

365 nfile = len(self.infiles) 

366 for i in range(nfile): 

367 tmp.append(numpy.zeros(imshape_out, dtype=float)) 

368 maskedpixel = numpy.array(tmp) 

369 del tmp 

370 

371 # direction 

372 dirs = [] 

373 if len(self.direction) == nfile: 

374 dirs = self.direction 

375 else: 

376 casalog.post('direction information is extrapolated.') 

377 for i in range(nfile): 

378 dirs.append(self.direction[i % len(self.direction)]) 

379 

380 # maskwidth 

381 masks = [] 

382 if isinstance(self.maskwidth, int) or isinstance(self.maskwidth, float): 

383 for i in range(nfile): 

384 masks.append(self.maskwidth) 

385 elif isinstance(self.maskwidth, list): # and nfile != len(self.maskwidth): 

386 for i in range(nfile): 

387 masks.append(self.maskwidth[i % len(self.maskwidth)]) 

388 for i in range(len(masks)): 

389 masks[i] = 0.01 * masks[i] 

390 

391 # mask 

392 for i in range(nfile): 

393 self.realimage = create_4d_image(self.infiles[i], self.tmprealname[i]) 

394 self.imagimage = self.realimage.subimage(outfile=self.tmpimagname[i]) 

395 

396 # replace masked pixels with 0.0 

397 if not self.realimage.getchunk(getmask=True).all(): 

398 casalog.post("Replacing masked pixels with 0.0 in %d-th image" % (i)) 

399 self.realimage.replacemaskedpixels(0.0) 

400 self.realimage.close() 

401 self.imagimage.close() 

402 

403 # Below working images are all 4D regardless of dimension of input images 

404 # image shape for temporary images (always 4D) 

405 ia.open(self.tmprealname[0]) 

406 imshape = ia.shape() 

407 # ndim = len(imshape) 

408 ia.close() 

409 

410 if len(self.thresh) == 0: 

411 casalog.post('Use whole region') 

412 else: 

413 for i in range(nfile): 

414 self.realimage = ia.newimage(self.tmprealname[i]) 

415 for iaxis2 in range(imshape[2]): 

416 for iaxis3 in range(imshape[3]): 

417 pixmsk = self.realimage.getchunk([0, 0, iaxis2, iaxis3], 

418 [nx - 1, ny - 1, iaxis2, iaxis3]) 

419 for ix in range(pixmsk.shape[0]): 

420 for iy in range(pixmsk.shape[1]): 

421 if self.thresh[0] == self.nolimit: 

422 if pixmsk[ix][iy] > self.thresh[1]: 

423 maskedpixel[i][ix][iy][iaxis2][iaxis3] = pixmsk[ix][iy] 

424 pixmsk[ix][iy] = 0.0 

425 elif self.thresh[1] == self.nolimit: 

426 if pixmsk[ix][iy] < self.thresh[0]: 

427 maskedpixel[i][ix][iy][iaxis2][iaxis3] = pixmsk[ix][iy] 

428 pixmsk[ix][iy] = 0.0 

429 else: 

430 if pixmsk[ix][iy] < self.thresh[0] \ 

431 or pixmsk[ix][iy] > self.thresh[1]: 

432 maskedpixel[i][ix][iy][iaxis2][iaxis3] = pixmsk[ix][iy] 

433 pixmsk[ix][iy] = 0.0 

434 self.realimage.putchunk(pixmsk, [0, 0, iaxis2, iaxis3]) 

435 self.realimage.close() 

436 maskedvalue = None 

437 if any(maskedpixel.flatten() != 0.0): 

438 maskedvalue = maskedpixel.mean(axis=0) 

439 del maskedpixel 

440 

441 # set weight factor 

442 weights = numpy.ones(shape=(nfile, nx, ny), dtype=float) 

443 eps = 1.0e-5 

444 dtor = numpy.pi / 180.0 

445 for i in range(nfile): 

446 scan_direction = '' 

447 if abs(numpy.sin(dirs[i] * dtor)) < eps: # direction is around 0 deg 

448 maskw = 0.5 * nx * masks[i] 

449 scan_direction = 'horizontal' 

450 elif abs(numpy.cos(dirs[i] * dtor)) < eps: # direction is around 90 deg 

451 maskw = 0.5 * ny * masks[i] 

452 scan_direction = 'vertical' 

453 else: 

454 maskw = 0.5 * numpy.sqrt(nx * ny) * masks[i] 

455 for ix in range(nx): 

456 halfwx = (nx - 1) // 2 

457 for iy in range(ny): 

458 halfwy = (ny - 1) // 2 

459 if scan_direction == 'horizontal': 

460 # dd = abs(float(ix) - 0.5*(nx-1)) 

461 dd = abs(float(ix) - halfwx) # for CAS-9434 

462 elif scan_direction == 'vertical': 

463 # dd = abs(float(iy) - 0.5*(ny-1)) 

464 dd = abs(float(iy) - halfwy) # for CAS-9434 

465 else: 

466 tand = numpy.tan((dirs[i] - 90.0) * dtor) 

467 # dd = abs((float(ix) - 0.5*(nx-1)) * tand - (float(iy) - 0.5*(ny-1))) 

468 dd = abs((float(ix) - halfwx) * tand - (float(iy) - halfwy)) # for CAS-9434 

469 dd = dd / numpy.sqrt(1.0 + tand * tand) 

470 if dd < maskw: 

471 cosd = numpy.cos(0.5 * numpy.pi * dd / maskw) 

472 weights[i][ix][iy] = 1.0 - cosd * cosd 

473 if weights[i][ix][iy] == 0.0: 

474 weights[i][ix][iy] += eps * 0.01 

475 """ 

476 if abs(numpy.sin(dirs[i]*dtor)) < eps: 

477 # direction is around 0 deg 

478 maskw = 0.5 * nx * masks[i] 

479 for ix in range(nx): 

480 for iy in range(ny): 

481 dd = abs( float(ix) - 0.5 * (nx-1) ) 

482 if dd < maskw: 

483 cosd = numpy.cos(0.5*numpy.pi*dd/maskw) 

484 weights[i][ix][iy] = 1.0 - cosd * cosd 

485 if weights[i][ix][iy] == 0.0: 

486 weights[i][ix][iy] += eps*0.01 

487 elif abs(numpy.cos(dirs[i]*dtor)) < eps: 

488 # direction is around 90 deg 

489 maskw = 0.5 * ny * masks[i] 

490 for ix in range(nx): 

491 for iy in range(ny): 

492 dd = abs( float(iy) - 0.5 * (ny-1) ) 

493 if dd < maskw: 

494 cosd = numpy.cos(0.5*numpy.pi*dd/maskw) 

495 weights[i][ix][iy] = 1.0 - cosd * cosd 

496 if weights[i][ix][iy] == 0.0: 

497 weights[i][ix][iy] += eps*0.01 

498 else: 

499 maskw = 0.5 * numpy.sqrt( nx * ny ) * masks[i] 

500 for ix in range(nx): 

501 for iy in range(ny): 

502 tand = numpy.tan((dirs[i]-90.0)*dtor) 

503 dd = abs( ix * tand - iy - 0.5 * (nx-1) * tand + 0.5 * (ny-1) ) 

504 dd = dd / numpy.sqrt( 1.0 + tand * tand ) 

505 if dd < maskw: 

506 cosd = numpy.cos(0.5*numpy.pi*dd/maskw) 

507 weights[i][ix][iy] = 1.0 - cosd * cosd 

508 if weights[i][ix][iy] == 0.0: 

509 weights[i][ix][iy] += eps*0.01 

510 """ 

511 # shift 

512 xshift = -((ny - 1) // 2) 

513 yshift = -((nx - 1) // 2) 

514 for ix in range(int(xshift), 0, 1): 

515 tmp = weights[i, :, 0].copy() 

516 weights[i, :, 0:ny - 1] = weights[i, :, 1:ny].copy() 

517 weights[i, :, ny - 1] = tmp 

518 for iy in range(int(yshift), 0, 1): 

519 tmp = weights[i, 0:1].copy() 

520 weights[i, 0:nx - 1] = weights[i, 1:nx].copy() 

521 weights[i, nx - 1:nx] = tmp 

522 

523 # FFT 

524 for i in range(nfile): 

525 self.realimage = ia.newimage(self.tmprealname[i]) 

526 self.imagimage = ia.newimage(self.tmpimagname[i]) 

527 for iaxis2 in range(imshape[2]): 

528 for iaxis3 in range(imshape[3]): 

529 pixval = self.realimage.getchunk([0, 0, iaxis2, iaxis3], 

530 [nx - 1, ny - 1, iaxis2, iaxis3]) 

531 pixval = pixval.reshape((nx, ny)) 

532 pixfft = npfft.fft2(pixval) 

533 pixfft = pixfft.reshape((nx, ny, 1, 1)) 

534 self.realimage.putchunk(pixfft.real, [0, 0, iaxis2, iaxis3]) 

535 self.imagimage.putchunk(pixfft.imag, [0, 0, iaxis2, iaxis3]) 

536 del pixval, pixfft 

537 self.realimage.close() 

538 self.imagimage.close() 

539 

540 # weighted mean 

541 for ichan in range(imshape[2]): 

542 for iaxis3 in range(imshape[3]): 

543 pixout = numpy.zeros(shape=(nx, ny), dtype=complex) 

544 denom = numpy.zeros(shape=(nx, ny), dtype=float) 

545 for i in range(nfile): 

546 self.realimage = ia.newimage(self.tmprealname[i]) 

547 self.imagimage = ia.newimage(self.tmpimagname[i]) 

548 pixval = self.realimage.getchunk( 

549 [0, 0, ichan, iaxis3], [nx - 1, ny - 1, ichan, iaxis3]) \ 

550 + self.imagimage.getchunk([0, 0, ichan, iaxis3], 

551 [nx - 1, ny - 1, ichan, iaxis3]) * 1.0j 

552 pixval = pixval.reshape((nx, ny)) 

553 pixout = pixout + pixval * weights[i] 

554 denom = denom + weights[i] 

555 self.realimage.close() 

556 self.imagimage.close() 

557 pixout = pixout / denom 

558 pixout = pixout.reshape((nx, ny, 1, 1)) 

559 self.realimage = ia.newimage(self.tmprealname[0]) 

560 self.imagimage = ia.newimage(self.tmpimagname[0]) 

561 self.realimage.putchunk(pixout.real, [0, 0, ichan, iaxis3]) 

562 self.imagimage.putchunk(pixout.imag, [0, 0, ichan, iaxis3]) 

563 self.realimage.close() 

564 self.imagimage.close() 

565 

566 # inverse FFT 

567 self.realimage = ia.newimage(self.tmprealname[0]) 

568 self.imagimage = ia.newimage(self.tmpimagname[0]) 

569 for ichan in range(imshape[2]): 

570 for iaxis3 in range(imshape[3]): 

571 pixval = self.realimage.getchunk([0, 0, ichan, iaxis3], 

572 [nx - 1, ny - 1, ichan, iaxis3]) \ 

573 + self.imagimage.getchunk([0, 0, ichan, iaxis3], 

574 [nx - 1, ny - 1, ichan, iaxis3]) * 1.0j 

575 pixval = pixval.reshape((nx, ny)) 

576 pixifft = npfft.ifft2(pixval) 

577 pixifft = pixifft.reshape((nx, ny, 1, 1)) 

578 self.realimage.putchunk(pixifft.real, blc=[0, 0, ichan, iaxis3]) 

579 del pixval, pixifft 

580 if maskedvalue is not None: 

581 self.realimage.putchunk(self.realimage.getchunk() + maskedvalue) 

582 

583 # put result into outimage 

584 chunk = self.realimage.getchunk() 

585 outimage.putchunk(chunk.reshape(imshape_out)) 

586 # handling of output image mask 

587 maskstr = "" 

588 for name in self.infiles: 

589 if len(maskstr) > 0: 

590 maskstr += " || " 

591 maskstr += ("mask('%s')" % (name)) 

592 outimage.calcmask(maskstr, name="basketweaving", asdefault=True) 

593 

594 self.realimage.close() 

595 self.imagimage.close() 

596 outimage.close() 

597 

598 def finalize(self): 

599 pass 

600 

601 def cleanup(self): 

602 # finalize image analysis tool 

603 if hasattr(self, 'image') and self.image is not None: 

604 if self.image.isopen(): 

605 self.image.done() 

606 tools = ['convimage', 'imageorg', 'realimage', 'imagimage'] 

607 for t in tools: 

608 if hasattr(self, t): 

609 v = getattr(self, t) 

610 if v and v.isopen(): 

611 v.done(remove=True) 

612 

613 # remove tempolary files 

614 filelist = ['tmpmskname', 'tmpconvname', 'tmppolyname', 

615 'tmprealname', 'tmpimagname'] 

616 existing_files = [] 

617 for s in filelist: 

618 if hasattr(self, s): 

619 f = getattr(self, s) 

620 if isinstance(f, list): 

621 for g in f: 

622 if os.path.exists(g): 

623 existing_files.append(g) 

624 else: 

625 if os.path.exists(f): 

626 existing_files.append(f) 

627 # CAS-5410 Use private tools inside task scripts 

628 if len(existing_files) > 0: 

629 ctsys.removetable(existing_files)