def predict_vis (msname, image, column="MODEL_DATA", 
        chanchunk=None, chanstart=0, chanstep=1):
  """Converts image into predicted visibilities"""

  # CASA to convert them
  casaimage = '{0}/{1}.img'.format(OUTPUT, os.path.basename(image))
  
  # convert to CASA image
  img = pyrap.images.image(image)
  img.saveas(casaimage)

  imgshp = img.shape()
  ftab = table(msname+'/SPECTRAL_WINDOW')
  numchans = ftab.getcol('NUM_CHAN')[0]
  # default chunk list is entire chanel range. Update this if needed
  chunklist = [ (0, numchans, None, None) ]
  if len(imgshp) == 4 and imgshp[0] > 1:
      nimgchan = imgshp[0]

      print("image cube has {0} channels, MS has {1} channels".format(nimgchan, numchans))

      imgchansize = imgshp[1]*imgshp[2]*imgshp[3]*4  # size of an image channel in bytes
      if chanchunk is None:
          mem_bytes = os.sysconf('SC_PAGE_SIZE') * os.sysconf('SC_PHYS_PAGES')  # e.g. 4015976448
          chanchunk = max((mem_bytes/20)/imgchansize, 1)
          print("based on available memory ({0}), max image chunk is {1} channels".format(mem_bytes, chanchunk))
      if chanchunk < nimgchan:
          mschanstep = numchans*chanstep/nimgchan
          if numchans%nimgchan:
              warn("MS channels not evenly divisible into $nimgchan image channels, chunking may be incorrect")
          chunklist = []
          for chan0 in range(0, nimgchan, chanchunk):
              imch0, imch1 = chan0, (min(chan0+chanchunk, nimgchan)-1)
              msch0 = chanstart + imch0*mschanstep
              msnch = (imch1-imch0+1)*mschanstep/chanstep
              # overlap each chunk from 1 onwards by a half-chunk back to take care of extrapolated visibilties
              # from previous channel
              if imch0:
                  imch0 -= 1
                  msch0 -= mschanstep/2
                  msnch += mschanstep/2
              print("image chunk {0}~{1} corresponds to MS chunk {2}~{3}".format(imch0, imch1, msch0, msch0+msnch-1))
              chunklist.append((msch0, msnch, imch0, imch1))
              
  # even in fill-model mode where it claims to ignore image parameters, the image channelization
  # arguments need to be "just so" as per below, otherwise it gives a GridFT: weights all zero message
  kw0 = {}
  kw0.update(ms=msname, model=casaimage, 
      niter=0, fixed=1, mode="channel", operation="csclean",
      img_nchan=1, img_chanstart=chanstart, img_chanstep=numchans*chanstep)
  kw0['fillmodel'] = 1
                   
  blc = [0]*len(imgshp)
  trc = [ x-1 for x in imgshp ]

  # now loop over image frequency chunks
  for ichunk, (mschanstart, msnumchans, imgch0, imgch1) in enumerate(chunklist):
      if len(chunklist) > 1:
          blc[0], trc[0] = imgch0, imgch1
          print("writing CASA image for slice {0} {1}".format(blc, trc))
          casaimage1 = II("{0}.{1}.img".format(image, ichunk1))
          rm_fr(casaimage1)
          print("writing CASA image for slice {0} {1} to {2}".format(blc, trc, casaimage1))
          img.subimage(blc, trc, dropdegenerate=False).saveas(casaimage1)
          kw0.update(model=casaimage1)
      else:
          img.unlock()
      # setup imager options
      kw0.update(chanstart=mschanstart, chanstep=chanstep, nchan=msnumchans)
      print("predicting visibilities into MODEL_DATA")

      _run(predict=True, **kw0)
      if len(chunklist) > 1:
        rm_fr(casaimage1)   
  rm_fr(casaimage)
  
  if column != "MODEL_DATA":
    print('Data was predicted to MODEL_DATA column. Will now copy it to the {} column as requested'.format(column))
    msutils.copycol(msname=msname, fromcol="MODEL_DATA", tocol=column)
def predict_vis ( msname , image , column = "MODEL_DATA" , chanchunk = None , chanstart = 0 , chanstep = 1 ) : casaimage = '{0}/{1}.img' . format ( OUTPUT , os . path . basename ( image ) ) img = pyrap . images . image ( image ) img . saveas ( casaimage ) imgshp = img . shape ( ) ftab = table ( msname + '/SPECTRAL_WINDOW' ) numchans = ftab . getcol ( 'NUM_CHAN' ) [ 0 ] chunklist = [ ( 0 , numchans , None , None ) ] if len ( imgshp ) == 4 and imgshp [ 0 ] > 1 : nimgchan = imgshp [ 0 ] print ( "image cube has {0} channels, MS has {1} channels" . format ( nimgchan , numchans ) ) imgchansize = imgshp [ 1 ] * imgshp [ 2 ] * imgshp [ 3 ] * 4 if chanchunk is None : mem_bytes = os . sysconf ( 'SC_PAGE_SIZE' ) * os . sysconf ( 'SC_PHYS_PAGES' ) chanchunk = max ( ( mem_bytes / 20 ) / imgchansize , 1 ) print ( "based on available memory ({0}), max image chunk is {1} channels" . format ( mem_bytes , chanchunk ) ) if chanchunk < nimgchan : mschanstep = numchans * chanstep / nimgchan if numchans % nimgchan : warn ( "MS channels not evenly divisible into $nimgchan image channels, chunking may be incorrect" ) chunklist = [ ] for chan0 in range ( 0 , nimgchan , chanchunk ) : imch0 , imch1 = chan0 , ( min ( chan0 + chanchunk , nimgchan ) - 1 ) msch0 = chanstart + imch0 * mschanstep msnch = ( imch1 - imch0 + 1 ) * mschanstep / chanstep if imch0 : imch0 -= 1 msch0 -= mschanstep / 2 msnch += mschanstep / 2 print ( "image chunk {0}~{1} corresponds to MS chunk {2}~{3}" . format ( imch0 , imch1 , msch0 , msch0 + msnch - 1 ) ) chunklist . append ( ( msch0 , msnch , imch0 , imch1 ) ) kw0 = { } kw0 . update ( ms = msname , model = casaimage , niter = 0 , fixed = 1 , mode = "channel" , operation = "csclean" , img_nchan = 1 , img_chanstart = chanstart , img_chanstep = numchans * chanstep ) kw0 [ 'fillmodel' ] = 1 blc = [ 0 ] * len ( imgshp ) trc = [ x - 1 for x in imgshp ] for ichunk , ( mschanstart , msnumchans , imgch0 , imgch1 ) in enumerate ( chunklist ) : if len ( chunklist ) > 1 : blc [ 0 ] , trc [ 0 ] = imgch0 , imgch1 print ( "writing CASA image for slice {0} {1}" . format ( blc , trc ) ) casaimage1 = II ( "{0}.{1}.img" . format ( image , ichunk1 ) ) rm_fr ( casaimage1 ) print ( "writing CASA image for slice {0} {1} to {2}" . format ( blc , trc , casaimage1 ) ) img . subimage ( blc , trc , dropdegenerate = False ) . saveas ( casaimage1 ) kw0 . update ( model = casaimage1 ) else : img . unlock ( ) kw0 . update ( chanstart = mschanstart , chanstep = chanstep , nchan = msnumchans ) print ( "predicting visibilities into MODEL_DATA" ) _run ( predict = True , ** kw0 ) if len ( chunklist ) > 1 : rm_fr ( casaimage1 ) rm_fr ( casaimage ) if column != "MODEL_DATA" : print ( 'Data was predicted to MODEL_DATA column. Will now copy it to the {} column as requested' . format ( column ) ) msutils . copycol ( msname = msname , fromcol = "MODEL_DATA" , tocol = column )
Converts image into predicted visibilities