Coverage for /wheeldirectory/casa-6.7.0-12-py3.10.el8/lib/py/lib/python3.10/site-packages/casatasks/private/imagerhelpers/imager_base.py: 74%
405 statements
« prev ^ index » next coverage.py v7.6.4, created at 2024-11-01 07:19 +0000
« prev ^ index » next coverage.py v7.6.4, created at 2024-11-01 07:19 +0000
1import os
2import math
3import shutil
4import string
5import time
6import re
7import copy
8from typing import TYPE_CHECKING
10from casatools import (
11 synthesisimager,
12 synthesisdeconvolver,
13 synthesisnormalizer,
14 iterbotsink,
15 ctsys,
16 table,
17 image,
18)
19from casatasks import casalog
20from casatasks.private.imagerhelpers.summary_minor import SummaryMinor
21if TYPE_CHECKING:
22 from casatasks.private.imagerhelpers.input_parameters import ImagerParameters
24ctsys_hostinfo = ctsys.hostinfo
25_tb = table()
26_ia = image()
28"""
29A set of helper functions for tclean.
31Summary...
33"""
35#############################################
36class PySynthesisImager:
38 def __init__(self,params: 'ImagerParameters'):
39 ################ Tools
40 self.initDefaults()
42 # Check all input parameters, after partitioning setup.
44 # Selection Parameters. Dictionary of dictionaries, indexed by 'ms0','ms1',...
45 self.allselpars = params.getSelPars()
46 # Imaging/Deconvolution parameters. Same for serial and parallel runs
47 self.alldecpars = params.getDecPars()
48 self.allimpars = params.getImagePars()
49 self.allgridpars = params.getGridPars()
50 self.allnormpars = params.getNormPars()
51 self.weightpars = params.getWeightPars()
52 # Iteration parameters
53 self.iterpars = params.getIterPars() ## Or just params.iterpars
55 # CFCache params
56 self.cfcachepars = params.getCFCachePars()
57 ## Number of fields ( main + outliers )
58 self.NF = len(self.allimpars.keys())
59 self.stopMinor = {} ##[0]*self.NF
60 for immod in range(0, self.NF):
61 self.stopMinor[str(immod)] = 1.0
62 ## Number of nodes. This gets set for parallel runs
63 ## It can also be used serially to process the major cycle in pieces.
64 self.NN = 1
65 ## for debug mode automask incrementation only
66 self.ncycle = 0
67 # For the nmajor parameter
68 self.majorCnt = 0
70 # isvalid = self.checkParameters()
71 # if isvalid==False:
72 # casalog.post('Invalid parameters')
74 #############################################
75 # def checkParameters(self):
76 # # Copy the imagename from impars to decpars, for each field.
77 # for immod in range(0,self.NF):
78 # self.alldecpars[str(immod)]['imagename'] = self.allimpars[str(immod)]['imagename']
79 # return True
81 #############################################
82 def makeCFCache(self, exists):
83 # Make the CFCache and re-load it. The following calls become
84 # NoOps (in SynthesisImager.cc) if the gridder is not one
85 # which uses CFCache.
86 if exists:
87 casalog.post("CFCache already exists")
88 else:
89 self.dryGridding();
90 self.fillCFCache();
91 self.reloadCFCache();
94 def initializeImagers(self):
96 ## Initialize the tool for the current node
97 self.SItool = synthesisimager()
99 ## casalog.post('impars ', self.allimpars['0']['specmode'], 'frame', self.allimpars['0']['outframe'])
100 ## Send in selection parameters for all MSs in the list.
101 for mss in sorted((self.allselpars).keys()):
102 # if(self.allimpars['0']['specmode']=='cubedata'):
103 # self.allselpars[mss]['outframe']='Undefined'
104 self.SItool.selectdata(self.allselpars[mss])
105 # self.SItool.selectdata( **(self.allselpars[mss]) )
107 ## For each image-field, define imaging parameters
108 # nimpars = copy.deepcopy(self.allimpars)
109 # for fld in range(0,self.NF):
110 # self.SItool.defineimage( **( nimpars[str(fld)] ) )
112 # If cfcache directory already exists, assume that it is
113 # usable and is correct. makeCFCache call then becomes a
114 # NoOp.
116 cfCacheName=''
117 exists=False
118 if(self.allgridpars['0']['gridder'].startswith('awpr') or self.allgridpars['0']['gridder'].startswith('awph') ):
119 cfCacheName=self.allgridpars['0']['cfcache'];
120 if (cfCacheName == ''):
121 cfCacheName = self.allimpars['0']['imagename'] + '.cf'
122 self.allgridpars['0']['cfcache']= cfCacheName
123 exists = (os.path.exists(cfCacheName) and os.path.isdir(cfCacheName));
124 else:
125 cfCacheName=''
126 exists=True
128 for fld in range(0,self.NF):
129 # casalog.post("self.allimpars=",self.allimpars,"\n")
131 # print(f'####allimpars={self.allimpars[str(fld)]} \n allgridpars={self.allgridpars[str(fld)]}')
132 self.SItool.defineimage(
133 self.allimpars[str(fld)], self.allgridpars[str(fld)]
134 )
136 ###for cases when synthesisnormalizer is setup in c++ send the normalizer info
137 ###all images have the same normtype etc..so first one is good enough
138 self.SItool.normalizerinfo(self.allnormpars["0"])
139 ###commenting this out so that tuneSelect is done after weighting
140 ###CAS-11687
141 # For cube imaging: align the data selections and image setup
143 #if self.allimpars['0']['specmode'] != 'mfs' and self.allimpars['0']['specmode'] != 'cubedata':
144 # self.SItool.tuneselectdata()
145 ###For cubes create cfcache ahead of each partition trying
146 ### to create it as it is not multiprocess safe
147 if("cube" in self.allimpars['0']['specmode']):
148 self.makeCFCache(exists);
150 #############################################
152 def initializeDeconvolvers(self):
153 for immod in range(0, self.NF):
154 self.SDtools.append(synthesisdeconvolver())
155 self.SDtools[immod].setupdeconvolution(decpars=self.alldecpars[str(immod)])
157 #############################################
158 ## Overloaded by ParallelCont
159 def initializeNormalizers(self):
160 for immod in range(0, self.NF):
161 self.PStools.append(synthesisnormalizer())
162 normpars = self.allnormpars[str(immod)]
163 self.PStools[immod].setupnormalizer(normpars=normpars)
165 #############################################
167 def initializeIterationControl(self):
168 # note that in CASA5 this is casac.synthesisiterbot
169 self.IBtool = iterbotsink()
170 itbot = self.IBtool.setupiteration(iterpars=self.iterpars)
172 #############################################
173 def estimatememory(self):
174 # casalog.post("MEMORY usage ", self.SItool.estimatememory(), type(self.SItool.estimatememory()))
175 # griddermem=0
176 if self.SItool != None:
177 griddermem = self.SItool.estimatememory()
178 deconmem = 0
179 for immod in range(0, self.NF):
180 ims = self.allimpars[str(immod)]["imsize"]
181 if type(ims) == int:
182 ims = [ims, ims]
183 if len(ims) == 1:
184 ims.append(ims[0])
185 # casalog.post('shape', self.allimpars[str(immod)]['imsize'], len(ims) )
186 # casalog.post("DECON mem usage ", self.SDtools[immod].estimatememory(ims))
187 if len(self.SDtools) > immod:
188 if self.SDtools != None:
189 deconmem += self.SDtools[immod].estimatememory(ims)
190 availmem = ctsys_hostinfo()["memory"]["available"]
191 if (deconmem + griddermem) > 0.8 * availmem:
192 casalog.post(
193 "Memory available "
194 + str(availmem)
195 + " kB is very close to amount of required memory "
196 + str(deconmem + griddermem)
197 + " kB",
198 "WARN",
199 )
200 else:
201 casalog.post(
202 "Memory available "
203 + str(availmem)
204 + " kB and required memory "
205 + str(deconmem + griddermem)
206 + " kB",
207 "INFO2",
208 )
210 ############################################
211 def restoreImages(self):
212 for immod in range(0, self.NF):
213 self.SDtools[immod].restore()
215 #############################################
216 def pbcorImages(self):
217 for immod in range(0, self.NF):
218 self.SDtools[immod].pbcor()
220 def getSummary(self,fullsummary,fignum=1):
221 summ = self.IBtool.getiterationsummary()
222 casalog.post('getSummary call: fullsummary='+str(fullsummary))
223 if ('stopcode' in summ):
224 summ['stopDescription'] = self.getStopDescription(summ['stopcode'])
225 if ('summaryminor' in summ):
226 summ['summaryminor'] = SummaryMinor.convertMatrix(summ['summaryminor'],fullsummary)
227 #self.plotReport( summ, fignum )
229 return summ
231 #############################################
232 def deleteImagers(self):
233 if self.SItool != None:
234 self.SItool.done()
235 self.SItool = None
237 def deleteDeconvolvers(self):
238 for immod in range(0, len(self.SDtools)):
239 self.SDtools[immod].done()
241 def deleteNormalizers(self):
242 for immod in range(0, len(self.PStools)):
243 self.PStools[immod].done()
245 def deleteIterBot(self):
246 if self.IBtool != None:
247 self.IBtool.done()
249 def deleteCluster(self):
250 # casalog.post('no cluster to delete')
251 return
253 def deleteWorkDir(self):
254 # No .workdirectory to delete
255 return
257 def initDefaults(self):
258 # Reset globals/members
259 self.NF = 1
260 self.stopMinor = {"0": 1.0} # Flag to call minor cycle for this field or not.
261 self.NN = 1
262 self.SItool = None
263 self.SDtools = []
264 self.PStools = []
265 self.IBtool = None
267 #############################################
269 def deleteTools(self):
270 self.deleteImagers()
271 self.deleteDeconvolvers()
272 self.deleteNormalizers()
273 self.deleteIterBot()
274 self.deleteWorkDir()
275 self.initDefaults()
276 self.deleteCluster()
278 #############################################
280 def getStopDescription(self, stopflag):
281 stopreasons = [
282 "iteration limit", # 1
283 "threshold", # 2
284 "force stop", # 3
285 "no change in peak residual across two major cycles", # 4
286 "peak residual increased by more than 3 times from the previous major cycle", # 5
287 "peak residual increased by more than 3 times from the minimum reached", # 6
288 "zero mask", # 7
289 "any combination of n-sigma and other valid exit criterion", # 8
290 "reached nmajor", # 9
291 ]
292 if stopflag > 0:
293 return stopreasons[stopflag - 1]
294 return None
296 def hasConverged(self):
297 # Merge peak-res info from all fields to decide iteration parameters
298 time0 = time.time()
299 self.IBtool.resetminorcycleinfo()
300 for immod in range(0, self.NF):
301 initrec = self.SDtools[immod].initminorcycle()
302 # print('INIT Minor cycle dict {}'.format(initrec))
303 self.IBtool.mergeinitrecord(initrec)
305 # # Run interactive masking (and threshold/niter editors)
306 # self.runInteractiveGUI2()
308 # Check with the iteration controller about convergence.
309 reachedNmajor = (
310 self.iterpars["nmajor"] >= 0 and self.majorCnt >= self.iterpars["nmajor"]
311 )
312 stopflag = self.IBtool.cleanComplete(reachedMajorLimit=reachedNmajor)
313 if stopflag > 0:
314 casalog.post(
315 "Reached global stopping criterion : "
316 + self.getStopDescription(stopflag),
317 "INFO",
318 )
320 # revert the current automask to the previous one
321 # if self.iterpars['interactive']:
322 for immod in range(0, self.NF):
323 if self.alldecpars[str(immod)]["usemask"].count("auto") > 0:
324 prevmask = self.allimpars[str(immod)]["imagename"] + ".prev.mask"
325 if os.path.isdir(prevmask):
326 # Try to force rmtree even with an error as an nfs mounted disk gives an error
327 # shutil.rmtree(self.allimpars[str(immod)]['imagename']+'.mask')
328 shutil.rmtree(
329 self.allimpars[str(immod)]["imagename"] + ".mask",
330 ignore_errors=True,
331 )
332 # For NFS mounted disk it still leave .nfs* file(s)
333 if os.path.isdir(
334 self.allimpars[str(immod)]["imagename"] + ".mask"
335 ):
336 import glob
338 if glob.glob(
339 self.allimpars[str(immod)]["imagename"] + ".mask/.nfs*"
340 ):
341 for item in os.listdir(prevmask):
342 src = os.path.join(prevmask, item)
343 dst = os.path.join(
344 self.allimpars[str(immod)]["imagename"]
345 + ".mask",
346 item,
347 )
348 if os.path.isdir(src):
349 shutil.move(src, dst)
350 else:
351 shutil.copy2(src, dst)
352 shutil.rmtree(prevmask)
353 else:
354 shutil.move(
355 prevmask,
356 self.allimpars[str(immod)]["imagename"] + ".mask",
357 )
358 casalog.post(
359 "["
360 + str(self.allimpars[str(immod)]["imagename"])
361 + "] : Reverting output mask to one that was last used ",
362 "INFO",
363 )
364 casalog.post(
365 "***Time taken in checking hasConverged " + str(time.time() - time0),
366 "INFO3",
367 )
368 return stopflag > 0
370 #############################################
371 def updateMask(self):
372 # Setup mask for each field ( input mask, and automask )
373 maskchanged = False
374 time0 = time.time()
375 for immod in range(0, self.NF):
376 maskchanged = maskchanged | self.SDtools[immod].setupmask()
378 # Run interactive masking (and threshold/niter editors), if interactive=True
379 maskchanged = maskchanged | self.runInteractiveGUI2()
381 time1 = time.time()
382 casalog.post("Time to update mask " + str(time1 - time0) + "s", "INFO3")
383 ## Return a flag to say that the mask has changed or not.
384 return maskchanged
386 #############################################
387 def runInteractiveGUI2(self):
388 maskchanged = False
389 forcestop = True
390 if self.iterpars["interactive"] == True:
391 self.stopMinor = self.IBtool.pauseforinteraction()
392 # casalog.post("Actioncodes in python : " , self.stopMinor)
394 for akey in self.stopMinor:
395 if self.stopMinor[akey] < 0:
396 maskchanged = True
397 self.stopMinor[akey] = abs(self.stopMinor[akey])
399 # Check if force-stop has happened while savemodel != "none".
400 # If so, warn the user that unless the Last major cycle has happened,
401 # the model won't have been written into the MS, and to do a 'predict' run.
402 forcestop = True
403 for akey in self.stopMinor:
404 forcestop = forcestop and self.stopMinor[akey] == 3
406 if self.iterpars["savemodel"] != "none":
407 if forcestop == True:
408 self.predictModel()
409 # if self.iterpars['savemodel'] == "modelcolumn":
410 # wstr = "Saving model column"
411 # else:
412 # wstr = "Saving virtual model"
413 # casalog.post("Model visibilities may not have been saved in the MS even though you have asked for it. Please check the logger for the phrases 'Run (Last) Major Cycle' and '" + wstr +"'. If these do not appear, then please save the model via a separate tclean run with niter=0,calcres=F,calcpsf=F. It will pick up the existing model from disk and save/predict it. Reason for this : For performance reasons model visibilities are saved only in the last major cycle. If the X button on the interactive GUI is used to terminate a run before this automatically detected 'last' major cycle, the model isn't written. However, a subsequent tclean run as described above will predict and save the model. ","WARN")
415 # casalog.post('Mask changed during interaction : ', maskchanged)
416 return maskchanged or forcestop
418 #############################################
419 def makePSF(self):
420 self.makePSFCore()
421 divideInPython = (
422 self.allimpars["0"]["specmode"] == "mfs"
423 or self.allimpars["0"]["deconvolver"] == "mtmfs"
424 )
426 ### Gather PSFs (if needed) and normalize by weight
427 for immod in range(0, self.NF):
428 # for cube normalization is done in C++
429 if divideInPython:
430 self.PStools[immod].gatherpsfweight()
431 self.PStools[immod].dividepsfbyweight()
432 self.check_psf(immod)
434 def check_psf(self, immod):
435 if self.SDtools != []:
436 if immod <= len(self.SDtools) - 1:
437 self.SDtools[immod].checkrestoringbeam()
439 #############################################
440 def calcVisAppSens(self):
442 return self.SItool.apparentsens()
444 #############################################
446 def runMajorCycle(self, isCleanCycle=True):
447 # @param isCleanCycle: true if this is part of the major/minor cleaning loop, false if this is being used for some other purpose (such as generating the residual image)
448 if self.IBtool != None:
449 lastcycle = self.IBtool.cleanComplete(lastcyclecheck=True) > 0
450 else:
451 lastcycle = True
453 divideInPython = (
454 self.allimpars["0"]["specmode"] == "mfs"
455 or self.allimpars["0"]["deconvolver"] == "mtmfs"
456 )
457 ##norm is done in C++ for cubes
458 if not divideInPython:
459 self.runMajorCycleCore(lastcycle)
460 if isCleanCycle:
461 self.majorCnt += 1
462 if self.IBtool != None:
463 self.IBtool.endmajorcycle()
464 return
466 for immod in range(0, self.NF):
467 self.PStools[immod].dividemodelbyweight()
468 self.PStools[immod].scattermodel()
469 self.runMajorCycleCore(lastcycle)
470 if isCleanCycle:
471 self.majorCnt += 1
473 if self.IBtool != None:
474 self.IBtool.endmajorcycle()
475 ### Gather residuals (if needed) and normalize by weight
476 for immod in range(0, self.NF):
477 self.PStools[immod].gatherresidual()
478 self.PStools[immod].divideresidualbyweight()
479 self.PStools[immod].multiplymodelbyweight()
481 #############################################
482 def predictModel(self):
483 for immod in range(0, self.NF):
484 self.PStools[immod].dividemodelbyweight()
485 self.PStools[immod].scattermodel()
487 self.predictModelCore()
488 ###return the model images back to whatever state they were
489 for immod in range(0, self.NF):
490 self.PStools[immod].multiplymodelbyweight()
492 ## self.SItool.predictmodel();
494 #############################################
495 def dryGridding(self):
496 self.SItool.drygridding(**(self.cfcachepars))
498 #############################################
499 ## Overloaded for parallel runs
500 def fillCFCache(self):
501 cfcName = self.allgridpars["0"]["cfcache"]
502 cflist = []
503 if not (cfcName == ""):
504 cflist = [f for f in os.listdir(cfcName) if re.match(r"CFS*", f)]
505 # cflist = ["CFS_0_0_CF_0_0_0.im"];
506 self.cfcachepars["cflist"] = cflist
508 # self.SItool.fillcfcache(**(self.cfcachepars), self.allgridpars['0']['gridder'],cfcName);
510 self.SItool.fillcfcache(
511 cflist,
512 self.allgridpars["0"]["gridder"],
513 cfcName,
514 self.allgridpars["0"]["psterm"],
515 self.allgridpars["0"]["aterm"],
516 self.allgridpars["0"]["conjbeams"],
517 )
519 #############################################
520 def reloadCFCache(self):
521 self.SItool.reloadcfcache()
523 #############################################
524 def makePB(self):
525 ###for cube standard gridder pb is made in c++ with psf
526 if not (
527 "stand" in self.allgridpars["0"]["gridder"]
528 and "cube" in self.allimpars["0"]["specmode"]
529 ):
530 self.makePBCore()
531 for immod in range(0, self.NF):
532 self.PStools[immod].normalizeprimarybeam()
534 #############################################
535 def makePBCore(self):
536 self.SItool.makepb()
538 #############################################
539 def checkPB(self):
540 """Checks for common problem cases in the .pb image"""
541 if self.SItool is None:
542 # Seems to be None for specmode='mfs', parallel=True
543 return
545 import numpy as np
547 facetIdx = 0 # TODO iterate over facets
548 imagename = self.SItool.getImageName(facetIdx, "PB")
549 _ia.open(imagename)
550 # Case 1: non-zeroes on edge of .pb
551 pixelVals = _ia.getregion().copy()
552 pixelVals[1:-2][
553 1:-2
554 ] = 0 # zero out everything that isn't at the edge of 'right ascension' and 'declination' indexes
555 if pixelVals.max() > 0:
556 idx = np.unravel_index([pixelVals.argmax()], pixelVals.shape)
557 idx = [
558 x[0] for x in idx
559 ] # (array([296]), array([147]), array([0]), array([0])) --> [296, 147, 0, 0]
560 casalog.post(
561 f"Warning! Non-zero values at the edge of the .pb image can cause unexpected aliasing effects! (found value {pixelVals.max()} at index {idx})",
562 "WARN",
563 )
564 # release the image
565 _ia.close()
566 _ia.done()
568 #############################################
569 def makeSdImage(self):
570 self.makeSdImageCore()
571 for immod in range(0, self.NF):
572 self.PStools[immod].gatherresidual()
573 self.PStools[immod].divideresidualbyweight(singledish=True)
575 #############################################
576 def makeSdPSF(self):
577 self.makeSdPSFCore()
578 for immod in range(0, self.NF):
579 self.PStools[immod].gatherresidual()
580 self.PStools[immod].dividepsfbyweight()
582 #############################################
583 def makeSdImageCore(self):
584 self.SItool.makesdimage()
586 #############################################
587 def makeSdPSFCore(self):
588 self.SItool.makesdpsf()
590 #############################################
591 def makeImage(
592 self, imagetype="observed", image="", compleximage="", imagefieldid=0
593 ):
594 """
595 This should replace makeSDImage, makeSDPSF and makePSF
596 etc in the long run
597 But for now you can do the following images i.e string recognized by type
598 "observed", "model", "corrected", "psf", "residual", "singledish-observed",
599 "singledish", "coverage", "holography", "holography-observed"
600 For holography the FTmachine should be SDGrid and the baselines
601 selected should be those that are pointed up with the antenna which is rastering.
602 """
603 self.SItool.makeimage(imagetype, image, compleximage, imagefieldid)
605 #############################################
607 ## Overloaded for parallel runs
608 def setWeighting(self):
609 ## Set weighting parameters, and all pars common to all fields.
610 self.SItool.setweighting(**(self.weightpars))
611 ##moved the tuneselect after weighting so that
612 ##the weight densities use all the data selected CAS-11687
613 ###For cube imaging: align the data selections and image setup
614 ### the tuneSelect is now done in C++ CubeMajorCycleAlgorith.cc
615 # if self.allimpars['0']['specmode'] != 'mfs' and self.allimpars['0']['specmode'] != 'cubedata':
616 # self.SItool.tuneselectdata()
618 # casalog.post("get set density from python")
619 # self.SItool.getweightdensity()
620 # self.SItool.setweightdensity()
622 #############################################
623 ## Overloaded for parallel runs
624 def makePSFCore(self):
625 self.SItool.makepsf()
627 #############################################
628 ## Overloaded for parallel runs
629 def runMajorCycleCore(self, lastcycle):
630 controldict = {"lastcycle": lastcycle}
631 if ("0" in self.alldecpars) and ("usemask" in self.alldecpars["0"]):
632 controldict["usemask"] = self.alldecpars["0"]["usemask"]
633 self.SItool.executemajorcycle(controls=controldict)
635 #############################################
636 ## Overloaded for parallel runs
637 def predictModelCore(self):
638 self.SItool.predictmodel()
640 #############################################
642 def runMinorCycle(self):
643 return self.runMinorCycleCore()
645 #############################################
647 def runMinorCycleCore(self):
649 # Set False for release packages.
650 # Only set this to True for testing and debugging automask in parallel mode
651 # since in parallel mode, runtime setting of the enviroment variable
652 # currently does not work.
653 # False = disable always save intermediate images mode
654 alwaysSaveIntermediateImages = False
656 # Get iteration control parameters
657 iterbotrec = self.IBtool.getminorcyclecontrols()
659 # TT debug - comment out after debugging....
660 casalog.post("Minor Cycle controls : " + str(iterbotrec))
662 self.IBtool.resetminorcycleinfo()
664 #
665 # Run minor cycle
666 self.ncycle += 1
667 retval = False
668 for immod in range(0, self.NF):
669 if self.stopMinor[str(immod)] < 3:
671 # temporarily disable the check (=> always save the intermediate images
672 if alwaysSaveIntermediateImages or (
673 "SAVE_ALL_RESIMS" in os.environ
674 and os.environ["SAVE_ALL_RESIMS"] == "true"
675 ):
676 resname = self.allimpars[str(immod)]["imagename"] + ".residual"
677 tempresname = (
678 self.allimpars[str(immod)]["imagename"]
679 + ".inputres"
680 + str(self.ncycle)
681 )
682 if os.path.isdir(resname):
683 shutil.copytree(resname, tempresname)
684 modname = self.allimpars[str(immod)]["imagename"] + ".model"
685 tempmodname = (
686 self.allimpars[str(immod)]["imagename"]
687 + ".inputmod"
688 + str(self.ncycle)
689 )
690 if os.path.isdir(modname):
691 shutil.copytree(modname, tempmodname)
693 exrec = self.SDtools[immod].executeminorcycle(iterbotrecord=iterbotrec)
695 # casalog.post('.... iterdone for ', immod, ' : ' , exrec['iterdone'])
696 retval = retval or exrec["iterdone"] > 0
697 self.IBtool.mergeexecrecord(exrec, immod)
698 if alwaysSaveIntermediateImages or (
699 "SAVE_ALL_AUTOMASKS" in os.environ
700 and os.environ["SAVE_ALL_AUTOMASKS"] == "true"
701 ):
702 maskname = self.allimpars[str(immod)]["imagename"] + ".mask"
703 tempmaskname = (
704 self.allimpars[str(immod)]["imagename"]
705 + ".autothresh"
706 + str(self.ncycle)
707 )
708 if os.path.isdir(maskname):
709 shutil.copytree(maskname, tempmaskname)
711 # Some what duplicated as above but keep a copy of the previous mask
712 # for interactive automask to revert to it if the current mask
713 # is not used (i.e. reached deconvolution stopping condition).
714 ## no longer needed as of CAS-9386 for cubes.
715 # if self.iterpars['interactive'] and self.alldecpars[str(immod)]['usemask']=='auto-thresh':
716 if self.alldecpars[str(immod)]["usemask"].count("auto") > 0:
717 maskname = self.allimpars[str(immod)]["imagename"] + ".mask"
718 prevmaskname = (
719 self.allimpars[str(immod)]["imagename"] + ".prev.mask"
720 )
721 if os.path.isdir(maskname):
722 if os.path.isdir(prevmaskname):
723 shutil.rmtree(prevmaskname)
724 shutil.copytree(maskname, prevmaskname)
725 return retval
727 #############################################
728 def runMajorMinorLoops(self):
729 self.runMajorCycle(isCleanCycle=False)
730 while not self.hasConverged():
731 self.runMinorCycle()
732 self.runMajorCycle()
734 #############################################
736 def plotReport(self, summ={}, fignum=1):
738 if not (
739 "summaryminor" in summ
740 and "summarymajor" in summ
741 and "threshold" in summ
742 and summ["summaryminor"].shape[0] == 6
743 ):
744 casalog.post(
745 "Cannot make summary plot. Please check contents of the output dictionary from tclean."
746 )
747 return summ
749 import matplotlib.pyplot as pl
750 from numpy import max as amax
752 # 0 : iteration number (within deconvolver, per cycle)
753 # 1 : peak residual
754 # 2 : model flux
755 # 3 : cyclethreshold
756 # 4 : deconvolver id
757 # 5 : subimage id (channel, stokes..)
759 pl.ioff()
761 fig, ax = pl.subplots(nrows=1, ncols=1, num=fignum)
762 pl.clf()
763 minarr = summ["summaryminor"]
764 if minarr.size == 0:
765 casalog.post("Zero iteration: no summary plot is generated.", "WARN")
766 else:
767 iterlist = minarr[0, :]
768 eps = 0.0
769 peakresstart = []
770 peakresend = []
771 modfluxstart = []
772 modfluxend = []
773 itercountstart = []
774 itercountend = []
775 peakresstart.append(minarr[1, :][0])
776 modfluxstart.append(minarr[2, :][0])
777 itercountstart.append(minarr[0, :][0] + eps)
778 peakresend.append(minarr[1, :][0])
779 modfluxend.append(minarr[2, :][0])
780 itercountend.append(minarr[0, :][0] + eps)
781 for ii in range(0, len(iterlist) - 1):
782 if iterlist[ii] == iterlist[ii + 1]:
783 peakresend.append(minarr[1, :][ii])
784 peakresstart.append(minarr[1, :][ii + 1])
785 modfluxend.append(minarr[2, :][ii])
786 modfluxstart.append(minarr[2, :][ii + 1])
787 itercountend.append(iterlist[ii] - eps)
788 itercountstart.append(iterlist[ii] + eps)
790 peakresend.append(minarr[1, :][len(iterlist) - 1])
791 modfluxend.append(minarr[2, :][len(iterlist) - 1])
792 itercountend.append(minarr[0, :][len(iterlist) - 1] + eps)
794 # pl.plot( iterlist , minarr[1,:] , 'r.-' , label='peak residual' , linewidth=1.5, markersize=8.0)
795 # pl.plot( iterlist , minarr[2,:] , 'b.-' , label='model flux' )
796 # pl.plot( iterlist , minarr[3,:] , 'g--' , label='cycle threshold' )
798 pl.plot(itercountstart, peakresstart, "r.--", label="peak residual (start)")
799 pl.plot(
800 itercountend,
801 peakresend,
802 "r.-",
803 label="peak residual (end)",
804 linewidth=2.5,
805 )
806 pl.plot(itercountstart, modfluxstart, "b.--", label="model flux (start)")
807 pl.plot(
808 itercountend, modfluxend, "b.-", label="model flux (end)", linewidth=2.5
809 )
810 pl.plot(
811 iterlist, minarr[3, :], "g--", label="cycle threshold", linewidth=2.5
812 )
813 maxval = amax(minarr[1, :])
814 maxval = max(maxval, amax(minarr[2, :]))
816 bcols = ["b", "g", "r", "y", "c"]
817 minv = 1
818 niterdone = len(minarr[4, :])
820 if (
821 len(summ["summarymajor"].shape) == 1
822 and summ["summarymajor"].shape[0] > 0
823 ):
824 pl.vlines(
825 summ["summarymajor"], 0, maxval, label="major cycles", linewidth=2.0
826 )
828 pl.hlines(
829 summ["threshold"],
830 0,
831 summ["iterdone"],
832 linestyle="dashed",
833 label="threshold",
834 )
836 pl.xlabel("Iteration Count")
837 pl.ylabel("Peak Residual (red), Model Flux (blue)")
839 box = ax.get_position()
840 ax.set_position([box.x0, box.y0, box.width, box.height * 0.8])
842 pl.legend(
843 loc="lower center",
844 bbox_to_anchor=(0.5, 1.05),
845 ncol=3,
846 fancybox=True,
847 shadow=True,
848 )
850 pl.savefig("summaryplot_" + str(fignum) + ".png")
851 pl.ion()
853 return summ
855 #############################################
857 def unlockimages(self, imageid=0):
858 """
859 Will try to unlock images attached for the image or outlier field id
860 in this instance
861 """
862 retval = False
863 if len(self.PStools) > imageid:
864 retval = self.SItool.unlockimages(imageid)
865 retval = retval and self.PStools[imageid].unlockimages()
866 return retval
869#######################################################
870#######################################################