Coverage for /wheeldirectory/casa-6.7.0-12-py3.10.el8/lib/py/lib/python3.10/site-packages/casatasks/private/task

1################################################

2# Task to do wideband PB-corrections on images produced by MS-MFS.

4# v1.0: 2012.09.10, U.R.V.

6################################################

7import os

8import numpy as np

9import shutil

10from scipy import linalg

12from casatools import imager, image, quanta, measures

13from casatasks import casalog

15im = imager()

16ia = image()

17qa = quanta()

18me = measures()

21def widebandpbcor(

22 vis="",

23 imagename="mfim",

24 nterms=2,

25 threshold="1mJy",

26 action="pbcor",

27 reffreq="1.5GHz",

28 pbmin=0.001,

29 field="",

30 spwlist=[0],

31 chanlist=[0],

32 weightlist=[1],

33):

34 """

35 Wide-Band PB-correction. Specify a list of spwids and channel numbers at which

36 to compute primary beams. Supply weights to use for each beam. All three lists

37 spwlist, chanlist, weightlist must be of the same length. It is enough to

38 make one PB per spw (middle channel).

39 """

40 casalog.origin("widebandpbcor")

42 casalog.post(

43 "widebandpbcor is a temporary task, meant for use until a widebandpbcor option is enabled from within the tclean task.",

44 "WARN",

45 )

47 # Check nterms

48 if nterms < 1:

49 raise ValueError("Please specify a valid number of Taylor-terms (>0)")

51 # Check that all required input images exist on disk

52 taylorlist = []

53 residuallist = []

54 for i in range(0, nterms):

55 taylorlist.append(imagename + ".image.tt" + str(i))

56 residuallist.append(imagename + ".residual.tt" + str(i))

57 if not os.path.exists(taylorlist[i]):

58 raise RuntimeError(

59 "Taylor-coeff Restored Image : " + taylorlist[i] + " not found."

60 )

61 if not os.path.exists(residuallist[i]):

62 raise RuntimeError(

63 "Taylor-coeff Residual Image : " + residuallist[i] + " not found."

64 )

66 casalog.post(

67 "Using images " + str(taylorlist) + " and " + str(residuallist), "NORMAL"

68 )

70 ## Read imsize, cellsize, reffreq, phasecenter from the input images.

71 ia.open(taylorlist[0])

72 imsize = [ia.shape()[0], ia.shape()[1]]

73 iminfo = ia.summary()

74 ia.close()

76 ## Get cell size

77 cellx = (

78 str(

79 abs(

80 qa.convert(

81 qa.quantity(iminfo["incr"][0], iminfo["axisunits"][0]), "arcsec"

82 )["value"]

83 )

84 )

85 + " arcsec"

86 )

87 celly = (

88 str(

89 abs(

90 qa.convert(

91 qa.quantity(iminfo["incr"][1], iminfo["axisunits"][1]), "arcsec"

92 )["value"]

93 )

94 )

95 + " arcsec"

96 )

98 ## Get phasecenter

99 pcra = qa.quantity(iminfo["refval"][0], iminfo["axisunits"][0])

100 pcdec = qa.quantity(iminfo["refval"][1], iminfo["axisunits"][1])

101 pcdir = me.direction("J2000", pcra, pcdec)

102

103 phasecenter = (

104 "J2000 " + qa.formxxx(pcdir["m0"], "hms") + " " + qa.formxxx(pcdir["m1"], "dms")

105 )

106

107 ## Get reffreq, if not specified.

108 if len(reffreq) == 0:

109 rfreq = qa.convert(

110 qa.quantity(iminfo["refval"][3], iminfo["axisunits"][3]), "GHz"

111 )

112 reffreq = str(rfreq["value"]) + rfreq["unit"]

113

114 casalog.post(

115 "Using imsize : "

116 + str(imsize)

117 + " and cellsize : "

118 + str(cellx)

119 + ", "

120 + str(celly)

121 + " with phasecenter : "

122 + phasecenter

123 + " and reffreq : "

124 + reffreq,

125 "NORMAL",

126 )

127

128 if type(threshold) != str:

129 raise ValueError("Please specify imthreshold as a string with units")

130

131 imthreshold = qa.convert(qa.quantity(threshold), "Jy")["value"]

132

133 ret = False

134 if action == "pbcor":

135

136 # Extract thresholds

137 pbthreshold = pbmin

138

139 casalog.post("Using a pblimit of " + str(pbthreshold))

140

141 if len(chanlist) < nterms or len(spwlist) < nterms or len(weightlist) < nterms:

142 raise ValueError(

143 "Please specify channel/spw/weight lists of lengths >= nterms"

144 )

145

146 # Make a list of PBs from all the specifield spw/chan pairs.

147 if len(chanlist) != len(spwlist):

148 raise ValueError("Spwlist and Chanlist must be the same length")

149

150 if len(spwlist) != len(weightlist):

151 raise ValueError("Spwlist and Weightlist must be the same length")

152

153 casalog.post(

154 "Calculating PBs for spws : "

155 + str(spwlist)

156 + " weightlist : "

157 + str(weightlist)

158 + " chanlist : "

159 + str(chanlist),

160 "NORMAL",

161 )

162

163 pbdirname = imagename + ".pbcor.workdirectory"

164 if not os.path.exists(pbdirname):

165 os.system("mkdir " + pbdirname)

166

167 pblist = _makePBList(

168 msname=vis,

169 pbprefix=pbdirname + "/" + imagename + ".pb.",

170 field=field,

171 spwlist=spwlist,

172 chanlist=chanlist,

173 imsize=imsize,

174 cellx=cellx,

175 celly=celly,

176 phasecenter=phasecenter,

177 )

178

179 pbcubename = pbdirname + "/" + imagename + ".pb.cube"

180

181 casalog.post(

182 "Concatenating PBs at all specified frequencies into a cube", "NORMAL"

183 )

184 newia = ia.imageconcat(

185 outfile=pbcubename, infiles=pblist, relax=True, overwrite=True

186 )

187 newia.close()

188

189 # Delete individual pb images

190 for pbim in pblist:

191 shutil.rmtree(pbim)

192

193 # Make lists of image names

194 pblist = []

195 imlist = []

196 imlistpbcor = []

197 for i in range(0, nterms):

198 imlist.append(imagename + ".image.tt" + str(i))

199 pblist.append(pbdirname + "/" + imagename + ".pb.tt" + str(i))

200 imlistpbcor.append(imagename + ".pbcor.image.tt" + str(i))

201

202 # Calculate Taylor Coeffs from this cube

203 ret = _calcTaylorFromCube(

204 imtemplate=imagename + ".image.tt0",

205 reffreq=reffreq,

206 cubename=pbcubename,

207 newtay=pblist,

208 pbthreshold=pbthreshold,

209 weightlist=weightlist,

210 )

211

212 # Calculate PB alpha and beta ( just for information )

213 pbalphaname = pbdirname + "/" + imagename + ".pb.alpha"

214 ret = _calcPBAlpha(

215 pbtay=pblist, pbthreshold=pbthreshold, pbalphaname=pbalphaname

216 )

217

218 # Divide out the PB polynomial

219 ret = _dividePBTaylor(imlist, pblist, imlistpbcor, pbthreshold)

220

221 # Recalculate Alpha.

222 if ret == True or action == "calcalpha":

223 if nterms == 1:

224 casalog.post("Cannot compute spectral index with only one image", "WARN")

225 return

226

227 if ret == True:

228 imname = imagename + ".pbcor"

229 else:

230 imname = imagename

231 residuallist = []

232 imagelist = []

233 for ii in range(0, nterms):

234 residuallist.append(imagename + ".residual.tt" + str(ii))

235 imagelist.append(imname + ".image.tt" + str(ii))

236 _compute_alpha_beta(

237 imname, nterms, imagelist, residuallist, imthreshold, [], True

238 )

239

240

241###############################################

242

243

244def _calcPBAlpha(pbtay=[], pbthreshold=0.1, pbalphaname="pbalpha.im"):

245 nterms = len(pbtay)

246 if nterms < 2:

247 return False

248

249 if os.path.exists(pbalphaname):

250 rmcmd = "rm -rf " + pbalphaname

251 os.system(rmcmd)

252 if not os.path.exists(pbalphaname):

253 cpcmd = "cp -r " + pbtay[0] + " " + pbalphaname

254 os.system(cpcmd)

255

256 ia.open(pbalphaname)

257 alpha = ia.getchunk()

258 ia.close()

259

260 ptay = []

261 ia.open(pbtay[0])

262 ptay.append(ia.getchunk())

263 ia.close()

264 ia.open(pbtay[1])

265 ptay.append(ia.getchunk())

266 ia.close()

267

268 ptay[0][ptay[0] < pbthreshold] = 1.0

269 ptay[1][ptay[0] < pbthreshold] = 0.0

270

271 alpha = ptay[1] / ptay[0]

272

273 ia.open(pbalphaname)

274 ia.putchunk(alpha)

275 ia.calcmask(mask='"' + pbtay[0] + '"' + ">" + str(pbthreshold))

276 ia.close()

277

278

279#################################################

280def _makePBList(

281 msname="",

282 pbprefix="",

283 field="",

284 spwlist=[],

285 chanlist=[],

286 imsize=[],

287 cellx="10.0arcsec",

288 celly="10.0arcsec",

289 phasecenter="",

290):

291 # casalog.post('Making PB List from the following spw,chan pairs')

292 pblist = []

293 try:

294 for aspw in range(0, len(spwlist)):

295 # casalog.post('spw=', spwlist[aspw], ' chan=',chanlist[aspw])

296 im.open(msname)

297 sspw = str(spwlist[aspw]) + ":" + str(chanlist[aspw])

298 ret = im.selectvis(field=field, spw=sspw)

299 if ret == False:

300 msg = "Error in constructing PBs at the specified spw:chan of " + sspw

301 raise RuntimeError(msg)

302 im.defineimage(

303 nx=imsize[0],

304 ny=imsize[1],

305 cellx=cellx,

306 celly=celly,

307 nchan=1,

308 start=chanlist[aspw],

309 stokes="I",

310 mode="channel",

311 spw=[spwlist[aspw]],

312 phasecenter=phasecenter,

313 )

314 im.setvp(dovp=True)

315 pbname = pbprefix + str(spwlist[aspw]) + "." + str(chanlist[aspw])

316 im.makeimage(type="pb", image=pbname, compleximage="", verbose=False)

317 pblist.append(pbname)

318 im.close()

319 # shutil.rmtree('evlavp.tab')

320 except Exception as exc:

321 msg = (

322 "Error in constructing PBs at the specified frequencies. Please check spwlist"

323 "/chanlist. Exception: {}".format(exc)

324 )

325 raise RuntimeError(msg)

326

327 return pblist

328

329

330##############################################################################

331def _calcTaylorFromCube(

332 imtemplate="",

333 reffreq="1.42GHz",

334 cubename="sim.pb",

335 newtay=[],

336 pbthreshold=0.0001,

337 weightlist=[],

338):

339 for tay in range(0, len(newtay)):

340 if os.path.exists(newtay[tay]):

341 rmcmd = "rm -rf " + newtay[tay]

342 os.system(rmcmd)

343 if not os.path.exists(newtay[tay]):

344 cpcmd = "cp -r " + imtemplate + " " + newtay[tay]

345 os.system(cpcmd)

346

347 ia.open(cubename)

348 pbcube = ia.getchunk()

349 csys = ia.coordsys()

350 shp = ia.shape()

351 ia.close()

352

353 if csys.axiscoordinatetypes()[3] == "Spectral":

354 # restfreq = csys.restfrequency()['value'][0]/1e+09; # convert more generally..

355 restfreq = csys.referencevalue()["numeric"][3] / 1e09

356 # convert more generally..

357 freqincrement = csys.increment()["numeric"][3] / 1e09

358 freqlist = []

359 for chan in range(0, shp[3]):

360 freqlist.append(restfreq + chan * freqincrement)

361 elif csys.axiscoordinatetypes()[3] == "Tabular":

362 freqlist = (csys.torecord()["tabular2"]["worldvalues"]) / 1e09

363 else:

364 raise RuntimeError("Unknown frequency axis. Exiting.")

365

366 reffreqGHz = qa.convert(qa.quantity(reffreq), "GHz")["value"]

367 freqs = (np.array(freqlist, "f") - reffreqGHz) / reffreqGHz

368

369 casalog.post(

370 "Using PBs at "

371 + str(freqlist)

372 + " GHz, to compute Taylor coefficients about "

373 + str(reffreqGHz)

374 + " GHz",

375 "NORMAL",

376 )

377 # casalog.post(freqs)

378

379 nfreqs = len(freqlist)

380 if len(weightlist) == 0:

381 weightarr = np.ones(freqs.shape)

382 else:

383 if len(weightlist) != nfreqs:

384 raise RuntimeError(

385 "Weight list must be the same length as nFreq : " + nfreqs

386 )

387 else:

388 weightarr = np.array(weightlist)

389

390 nterms = len(newtay)

391 if nterms > 5:

392 raise RuntimeError("Cannot handle more than 5 terms for PB computation")

393

394 ptays = []

395 if nterms > 0:

396 ia.open(newtay[0])

397 ptays.append(ia.getchunk())

398 ptays[0].fill(0.0)

399 ia.close()

400 pstart = [0.0]

401 if nterms > 1:

402 ia.open(newtay[1])

403 ptays.append(ia.getchunk())

404 ptays[1].fill(0.0)

405 ia.close()

406 pstart = [0.0, 0.0]

407 if nterms > 2:

408 ia.open(newtay[2])

409 ptays.append(ia.getchunk())

410 ptays[2].fill(0.0)

411 ia.close()

412 pstart = [0.0, 0.0, 0.0]

413 if nterms > 3:

414 ia.open(newtay[3])

415 ptays.append(ia.getchunk())

416 ptays[3].fill(0.0)

417 ia.close()

418 pstart = [0.0, 0.0, 0.0, 0.0]

419 if nterms > 4:

420 ia.open(newtay[4])

421 ptays.append(ia.getchunk())

422 ptays[4].fill(0.0)

423 ia.close()

424 pstart = [0.0, 0.0, 0.0, 0.0, 0.0]

425

426 ##### Fit a nterms-term polynomial to each point !!!

427

428 # Linear fit

429 ptays = _linfit(ptays, freqs, pbcube[:, :, 0, :], weightarr, pbthreshold)

430

431 # Set all values below the pbthreshold, to zero

432 for tt in range(0, nterms):

433 ptays[tt][ptays[0] < pbthreshold] = 0.0

434

435 # Write to disk.

436 if nterms > 0:

437 ia.open(newtay[0])

438 ia.putchunk(ptays[0])

439 ia.close()

440 if nterms > 1:

441 ia.open(newtay[1])

442 ia.putchunk(ptays[1])

443 ia.close()

444 if nterms > 2:

445 ia.open(newtay[2])

446 ia.putchunk(ptays[2])

447 ia.close()

448 if nterms > 3:

449 ia.open(newtay[3])

450 ia.putchunk(ptays[3])

451 ia.close()

452 if nterms > 4:

453 ia.open(newtay[4])

454 ia.putchunk(ptays[4])

455 ia.close()

456

457

458####################################################

459def _linfit(ptays, freqs, pcube, wts, pbthresh):

460 # casalog.post 'Calculating PB Taylor Coefficients by applying Inv Hessian to Taylor-weighted sums')

461 nterms = len(ptays)

462 hess = np.zeros((nterms, nterms))

463 rhs = np.zeros((nterms, 1))

464 soln = np.zeros((nterms, 1))

465 shp = ptays[0].shape

466 if len(freqs) != pcube.shape[2]:

467 casalog.post(

468 "Mismatch in frequency axes : " + len(freqs) + " and " + pcube.shape[2],

469 "SEVERE",

470 )

471 return ptays

472 if len(freqs) != len(wts):

473 casalog.post(

474 "Mismatch in lengths of freqs : " + len(freqs) + " and wts : " + len(wts),

475 "SEVERE",

476 )

477 return ptays

478

479 for ii in range(0, nterms):

480 for jj in range(0, nterms):

481 hess[ii, jj] = np.mean(freqs ** (ii + jj) * wts)

482

483 normval = hess[0, 0]

484 hess = hess / normval

485

486 invhess = linalg.inv(hess)

487 casalog.post("Hessian : " + str(hess), "NORMAL")

488 casalog.post("Inv Hess : " + str(invhess), "NORMAL")

489

490 casalog.post("Calculating Taylor-coefficients for the PB spectrum", "NORMAL")

491

492 for x in range(0, shp[0]):

493 for y in range(0, shp[1]):

494 if (

495 pcube[x, y, 0] > pbthresh

496 ): # Calculate coeffs only where the largest beam is above thresh

497 for ii in range(0, nterms):

498 rhs[ii, 0] = (

499 np.mean((freqs ** (ii)) * pcube[x, y, :] * wts) / normval

500 )

501 soln = np.dot(invhess, rhs)

502 for ii in range(0, nterms):

503 ptays[ii][x, y] = soln[ii, 0]

504 if x % 100 == 0:

505 casalog.post("--- finished rows " + str(x) + " to " + str(x + 99), "NORMAL")

507 return ptays

510#############

513#############

514def _dividePB(nterms, pbcoeffs, targetpbs):

515 if len(pbcoeffs) != nterms or len(targetpbs) != nterms:

516 casalog.post(

517 "To divide out the PB spectrum, PB coeffs and target images must have same nterms",

518 "SEVERE",

519 )

520 return []

521 correctedpbs = []

522

523 if nterms == 1:

524 det = pbcoeffs[0]

525 det[abs(det) == 0.0] = 1.0

526 correctedpbs.append(targetpbs[0] / det)

527

528 if nterms == 2:

529 det = pbcoeffs[0] ** 2

530 det[abs(det) == 0.0] = 1.0

531 correctedpbs.append(pbcoeffs[0] * targetpbs[0] / det)

532 correctedpbs.append(

533 (-1 * pbcoeffs[1] * targetpbs[0] + pbcoeffs[0] * targetpbs[1]) / det

534 )

535

536 if nterms == 3:

537 det = pbcoeffs[0] ** 3

538 det[abs(det) == 0.0] = 1.0

539 correctedpbs.append((pbcoeffs[0] ** 2) * targetpbs[0] / det)

540 correctedpbs.append(

541 (

542 -1 * pbcoeffs[0] * pbcoeffs[1] * targetpbs[0]

543 + (pbcoeffs[0] ** 2) * targetpbs[1]

544 )

545 / det

546 )

547 correctedpbs.append(

548 (

549 (pbcoeffs[1] ** 2 - pbcoeffs[0] * pbcoeffs[2]) * targetpbs[0]

550 + (-1 * pbcoeffs[0] * pbcoeffs[1]) * targetpbs[1]

551 + (pbcoeffs[0] ** 2) * targetpbs[2]

552 )

553 / det

554 )

555

556 return correctedpbs

557

558

559#####

560##############################################################################

561def _dividePBTaylor(imlist=[], pblist=[], imlistpbcor=[], pbthreshold=0.1):

562 casalog.post("Dividing the Image polynomial by the PB polynomial", "NORMAL")

563

564 if len(imlist) != len(pblist):

565 raise RuntimeError("Image list and PB list must be of the same length")

566

567 if len(imlist) != len(imlistpbcor):

568 raise RuntimeError("Image list and imlistpbcor must be of the same length")

569

570 nterms = len(imlist)

571

572 # Read PB coefficient images

573 pbcoeffs = []

574 for tay in range(0, nterms):

575 if not os.path.exists(pblist[tay]):

576 raise RuntimeError("PB Coeff " + pblist[tay] + " does not exist ")

577

578 ia.open(pblist[tay])

579 pbcoeffs.append(ia.getchunk())

580 ia.close()

581

582 # Read Images to normalize

583 inpimages = []

584 for tay in range(0, nterms):

585 if not os.path.exists(imlist[tay]):

586 raise RuntimeError(

587 "Image Coeff " + imlist[tay] + " does not exist ", "SEVERE"

588 )

589

590 ia.open(imlist[tay])

591 inpimages.append(ia.getchunk())

592 ia.close()

593

594 # Divide the two polynomials.

595 normedims = _dividePB(nterms, pbcoeffs, inpimages)

596

597 if len(normedims) == 0:

598 raise RuntimeError("Could not divide the beam")

599

600 for tay in range(0, nterms):

601 imtemp = imlist[0]

602 normname = imlistpbcor[tay]

603 casalog.post("Writing PB-corrected images " + normname)

604 if os.path.exists(normname):

605 rmcmd = "rm -rf " + normname

606 os.system(rmcmd)

607 cpcmd = "cp -r " + imtemp + " " + normname

608 os.system(cpcmd)

609 ia.open(normname)

610 ia.putchunk(normedims[tay])

611 ia.calcmask(mask='"' + pblist[0] + '"' + ">" + str(pbthreshold))

612 ia.close()

613

614 return True

615

616

617###################################################

618

619

620####################################################

621def _compute_alpha_beta(

622 imagename, nterms, taylorlist, residuallist, threshold, beamshape, calcerror

623):

624 imtemplate = imagename + ".image.tt0"

625 nameintensity = imagename + ".image.tt0"

626 namealpha = imagename + ".image.alpha"

627 nameerror = namealpha + ".error"

628 namebeta = imagename + ".image.beta"

629

630 # if(os.path.exists(nameintensity)):

631 # rmcmd = 'rm -rf ' + nameintensity;

632 # os.system(rmcmd);

633 # casalog.post( 'Creating new image : ', nameintensity)

634 # cpcmd = 'cp -r ' + imtemplate + ' ' + nameintensity;

635 # os.system(cpcmd);

636

637 if os.path.exists(namealpha):

638 rmcmd = "rm -rf " + namealpha

639 os.system(rmcmd)

640 casalog.post("Creating new image : " + namealpha, "NORMAL")

641 cpcmd = "cp -r " + imtemplate + " " + namealpha

642 os.system(cpcmd)

643

644 if calcerror == True:

645 if os.path.exists(nameerror):

646 rmcmd = "rm -rf " + nameerror

647 os.system(rmcmd)

648 casalog.post("Creating new image : " + nameerror, "NORMAL")

649 cpcmd = "cp -r " + imtemplate + " " + nameerror

650 os.system(cpcmd)

651

652 if nterms > 2:

653 if os.path.exists(namebeta):

654 rmcmd = "rm -rf " + namebeta

655 os.system(rmcmd)

656 casalog.post("Creating new image : " + namebeta, "NORMAL")

657 cpcmd = "cp -r " + imtemplate + " " + namebeta

658 os.system(cpcmd)

659

660 ## Open and read the images to compute alpha/beta with

661 ptay = []

662 for i in range(0, nterms):

663 ia.open(taylorlist[i])

664 ptay.append(ia.getchunk())

665 ia.close()

666

667 ## If calc error, open residual images too

668 if calcerror == True:

669 pres = []

670 for i in range(0, nterms):

671 ia.open(residuallist[i])

672 pres.append(ia.getchunk())

673 ia.close()

674

675 ## Initialize arrays of empty images

676 # ia.open(nameintensity);

677 # intensity = ia.getchunk();

678 # intensity.fill(0.0);

679 # ia.close();

680 ia.open(namealpha)

681 alpha = ia.getchunk()

682 alpha.fill(0.0)

683 ia.close()

684 if calcerror:

685 ia.open(nameerror)

686 aerror = ia.getchunk()

687 aerror.fill(0.0)

688 ia.close()

689 if nterms > 2:

690 ia.open(namebeta)

691 beta = ia.getchunk()

692 beta.fill(0.0)

693 ia.close()

694

695 ## Calc alpha,beta from ptay0,ptay1,ptay2

696 # intensity = ptay[0];

697 # ia.open(nameintensity);

698 # ia.putchunk(intensity);

699 # ia.close();

700

701 ptay[0][ptay[0] < 1e-06] = 1.0

702 ptay[0][ptay[0] < threshold] = 1.0

703 ptay[1][ptay[0] < threshold] = 0.0

704 if nterms > 2:

705 ptay[2][ptay[0] < threshold] = 0.0

706

707 alpha = ptay[1] / ptay[0]

708

709 if nterms > 2:

710 beta = (ptay[2] / ptay[0]) - 0.5 * alpha * (alpha - 1)

711

712 ia.open(namealpha)

713 ia.putchunk(alpha)

714 ia.calcmask(mask='"' + nameintensity + '"' + ">" + str(threshold))

715 ia.setbrightnessunit("")

716 ia.close()

717 if nterms > 2:

718 ia.open(namebeta)

719 ia.putchunk(beta)

720 ia.calcmask(mask='"' + nameintensity + '"' + ">" + str(threshold))

721 ia.setbrightnessunit("")

722 ia.close()

723

724 # calc error

725 if calcerror:

726

727 pres[1][ptay[1] == 0.0] = 0.0

728 ptay[1][pres[1] == 0.0] = 1.0

729

730 aerror = np.abs(alpha) * np.sqrt(

731 (pres[0] / ptay[0]) ** 2 + (pres[1] / ptay[1]) ** 2

732 )

733 ia.open(nameerror)

734 ia.putchunk(aerror)

735 ia.calcmask(mask='"' + nameintensity + '"' + ">" + str(threshold))

736 ia.setbrightnessunit("")

737 ia.close()

738

739 # Set the new restoring beam, if beamshape was used

740 if beamshape != []:

741 if _set_clean_beam(nameintensity, beamshape) == False:

742 return False

743 if _set_clean_beam(namealpha, beamshape) == False:

744 return False

745 if nterms > 2:

746 if _set_clean_beam(namebeta, beamshape) == False:

747 return False

748 if calcerror == True:

749 if _set_clean_beam(nameerror, beamshape) == False:

750 return False

751

752

753####################################################

754# Set the restoring beam to the new one.

755def _set_clean_beam(imname, beamshape):

756 ia.open(imname)

757 try:

758 ia.setrestoringbeam(major=beamshape[0], minor=beamshape[1], pa=beamshape[2])

759 except Exception as e:

760 casalog.post("Error setting restoring beam : " + e, "WARN")

761 ia.close()

762 return False

763 ia.close()

764 return True

765

766

767####################################################

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

413 statements