LCOV - code coverage report
Current view: top level - synthesis/TransformMachines2 - VisModelData.cc (source / functions) Hit Total Coverage
Test: casacpp_coverage.info Lines: 16 945 1.7 %
Date: 2024-10-10 15:00:01 Functions: 5 36 13.9 %

          Line data    Source code
       1             : //# VisImagingWeight.cc: imaging weight calculation for a give buffer
       2             : //# Copyright (C) 2011
       3             : //# Associated Universities, Inc. Washington DC, USA.
       4             : //#
       5             : //# This library is free software; you can redistribute it and/or modify it
       6             : //# under the terms of the GNU Library General Public License as published by
       7             : //# the Free Software Foundation; either version 2 of the License, or (at your
       8             : //# option) any later version.
       9             : //#
      10             : //# This library is distributed in the hope that it will be useful, but WITHOUT
      11             : //# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
      12             : //# FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Library General Public
      13             : //# License for more details.
      14             : //#
      15             : //# You should have received a copy of the GNU Library General Public License
      16             : //# along with this library; if not, write to the Free Software Foundation,
      17             : //# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
      18             : //#
      19             : //# Correspondence concerning AIPS++ should be adressed as follows:
      20             : //#        Internet email: casa-feedback@nrao.edu.
      21             : //#        Postal address: AIPS++ Project Office
      22             : //#                        National Radio Astronomy Observatory
      23             : //#                        520 Edgemont Road
      24             : //#                        Charlottesville, VA 22903-2475 USA
      25             : //#
      26             : //#
      27             : //# $Id$
      28             : 
      29             : 
      30             : #include <casacore/casa/Utilities/CountedPtr.h>
      31             : #include <casacore/casa/Arrays/ArrayMath.h>
      32             : #include <casacore/casa/Arrays/Vector.h>
      33             : #include <casacore/casa/OS/Timer.h>
      34             : #include <casacore/casa/Containers/Record.h>
      35             : #include <casacore/casa/Logging/LogIO.h>
      36             : #include <casacore/tables/Tables/ScaRecordColDesc.h>
      37             : #include <casacore/tables/Tables/TableUtil.h>
      38             : #include <components/ComponentModels/ComponentList.h>
      39             : #include <casacore/ms/MSSel/MSSelection.h>
      40             : #include <casacore/ms/MSSel/MSSelectionTools.h>
      41             : #include <casacore/ms/MeasurementSets/MSSource.h>
      42             : #include <casacore/ms/MSSel/MSSourceIndex.h>
      43             : #include <casacore/ms/MeasurementSets/MSSourceColumns.h>
      44             : 
      45             : #include <msvis/MSVis/VisBuffer2.h>
      46             : #include <synthesis/TransformMachines2/VisModelData.h>
      47             : #include <synthesis/TransformMachines2/FTMachine.h>
      48             : #include <synthesis/TransformMachines2/SimpleComponentFTMachine.h>
      49             : #include <synthesis/TransformMachines2/GridFT.h>
      50             : //#include <synthesis/TransformMachines/rGridFT.h>
      51             : #include <synthesis/TransformMachines2/MosaicFTNew.h>
      52             : #include <synthesis/TransformMachines2/WProjectFT.h>
      53             : //#include <synthesis/TransformMachines/MultiTermFT.h>
      54             : #include <synthesis/TransformMachines2/MultiTermFTNew.h>
      55             : #include <synthesis/TransformMachines2/SetJyGridFT.h>
      56             : 
      57             : namespace {
      58             : 
      59           1 :   casa::VisModelDataI * createRefImVisModelData (){
      60           1 :     return new casa::refim::VisModelData ();
      61             :   }
      62             : 
      63             :   // FIXME compiler warning re: unused variable initializeVisModelDataFactory 
      64             :   bool initializeVisModelDataFactory = casa::VisModelDataI::setFactory (createRefImVisModelData, 1);
      65             : 
      66             : }
      67             : 
      68             : 
      69             : 
      70             : using namespace casacore;
      71             : namespace casa { //# NAMESPACE CASA - BEGIN
      72             : 
      73             : namespace refim {//namespace for refactor
      74             : 
      75             : using namespace casacore;
      76             : using namespace casa;
      77             : using namespace casacore;
      78             : using namespace casa::refim;
      79             : using namespace casacore;
      80             : using namespace casa::vi;
      81           1 :   VisModelData::VisModelData(): clholder_p(0), ftholder_p(0), flatholder_p(0), version_p(-1){
      82             :   
      83           1 :   cft_p=new SimpleComponentFTMachine();
      84           1 :   }
      85             : 
      86           2 :   VisModelData::~VisModelData(){
      87             : 
      88             : 
      89           2 :   }
      90             : 
      91             : VisModelDataI *
      92           0 : VisModelData::clone ()
      93             : {
      94           0 :     return new VisModelData (* this);
      95             : }
      96             : 
      97             : 
      98           0 :   Bool VisModelData::hasAnyModel(const MeasurementSet& thems, Vector<Int>& fieldids){
      99           0 :     Bool retval=False;
     100           0 :     fieldids.resize();
     101           0 :     MSColumns msc(thems);
     102           0 :     Vector<Int> fields=msc.fieldId().getColumn();
     103           0 :     const Sort::Order order=Sort::Ascending;
     104           0 :     const Int option=Sort::HeapSort | Sort::NoDuplicates;
     105           0 :     Int nfields=GenSort<Int>::sort (fields, order, option);\
     106           0 :     if (nfields>0) {
     107           0 :       for (Int j=0; j< 2; ++j){
     108           0 :         const Table* thetab=&thems;
     109           0 :         if (j==1)
     110           0 :           thetab=&(thems.source());
     111           0 :         for (Int k=0; k< nfields; ++k){
     112           0 :                         if(thetab->keywordSet().isDefined("definedmodel_field_"+String::toString(fields[k])))
     113             :                 {
     114           0 :                   String elkey=thetab->keywordSet().asString("definedmodel_field_"+String::toString(fields[k]));
     115           0 :                   if(thetab->keywordSet().isDefined(elkey)){
     116           0 :                     fieldids.resize(fieldids.nelements()+1, True);
     117           0 :                     fieldids[fieldids.nelements()-1]=fields[k];
     118             :                   }
     119           0 :                 }
     120             :         }
     121             :       }
     122             :     }
     123           0 :     if(fieldids.nelements() >0)
     124           0 :       retval=True;
     125             : 
     126           0 :     return retval;
     127           0 :   }
     128             : 
     129           0 : void VisModelData::listModel(const MeasurementSet& thems){
     130             :  
     131             :   //Table newTab(thems);
     132             : 
     133           0 :   MSColumns msc(thems);
     134           0 :   Vector<String> fldnames=msc.field().name().getColumn();
     135           0 :   Vector<Int> fieldids;
     136           0 :   LogIO logio;
     137           0 :   if(hasAnyModel(thems, fieldids)){
     138           0 :       for (uInt k=0; k< fieldids.nelements(); ++k){
     139             :                         
     140           0 :         logio << " " << fldnames[fieldids[k]] << " (id = " << fieldids[k] << ")" << LogIO::POST;
     141             :                           
     142             :       }
     143             :     }
     144             :     else{
     145           0 :       logio <<  " None." << LogIO::POST;
     146             :     }
     147             :   
     148             :     
     149           0 : }
     150             : 
     151           0 : void VisModelData::deleteDiskImage(MeasurementSet& theMS, const String& theKey){
     152           0 :   TableRecord theRec;
     153           0 :   if(!VisModelData::getModelRecord(theKey, theRec, theMS))
     154           0 :     return;
     155           0 :   if(theRec.isDefined("numft")){
     156           0 :     Int numft=0;
     157           0 :     numft=theRec.asInt("numft");
     158           0 :     for (Int k=0; k < numft; ++k){
     159           0 :       String ftname=String("ft_")+String::toString(k);
     160           0 :       const Record& ftrec=theRec.asRecord(ftname);
     161           0 :       if(ftrec.asRecord("container").isDefined("diskimage")){
     162           0 :         String diskim=ftrec.asRecord("container").asString("diskimage");
     163           0 :         if(TableUtil::canDeleteTable(diskim))
     164           0 :           TableUtil::deleteTable(diskim);
     165           0 :       }
     166             : 
     167           0 :     }
     168             :   }
     169             : 
     170             : 
     171           0 : }
     172           0 : void VisModelData::removeRecordByKey(MeasurementSet& theMS, const String& theKey){
     173             : 
     174             : 
     175           0 :   deleteDiskImage(theMS, theKey);
     176           0 :   if(Table::isReadable(theMS.sourceTableName()) &&theMS.source().nrow() > 0 ){
     177           0 :     if(theMS.source().keywordSet().isDefined(theKey)){
     178           0 :       Int rowid=theMS.source().keywordSet().asInt(theKey);
     179           0 :       TableRecord elrec;
     180             :       //Replace the model with an empty record
     181           0 :       MSSourceColumns srcCol(theMS.source());
     182           0 :       srcCol.sourceModel().put(rowid, elrec);
     183           0 :       theMS.source().rwKeywordSet().removeField(theKey);
     184           0 :     }
     185             :   }
     186             :   //Remove from Main table if it is there
     187           0 :   if(theMS.rwKeywordSet().isDefined(theKey))
     188           0 :     theMS.rwKeywordSet().removeField(theKey);
     189           0 : }
     190             : 
     191           0 : void VisModelData::clearModel(const MeasurementSet& thems){
     192             :   //  Table newTab(thems);
     193           0 :   MeasurementSet& newTab=const_cast<MeasurementSet& >(thems);
     194           0 :   if(!newTab.isWritable())
     195           0 :     return;
     196           0 :   auto part_block = newTab.getPartNames(true);
     197           0 :   Vector<String> theParts(part_block.begin( ),part_block.end( ));
     198           0 :   if(theParts.nelements() > 1){
     199           0 :     for (uInt k=0; k < theParts.nelements(); ++k){
     200           0 :       MeasurementSet subms(theParts[k], newTab.lockOptions(), Table::Update);
     201           0 :       clearModel(subms);
     202           0 :     }
     203           0 :     return;
     204             :   }
     205           0 :   Bool alreadyLocked=newTab.hasLock(True);
     206           0 :   if(!alreadyLocked)
     207           0 :     newTab.lock(True);
     208           0 :   if(Table::isReadable(newTab.sourceTableName())){
     209           0 :     if(!newTab.source().isWritable())
     210           0 :       return;
     211           0 :     newTab.source().lock(True);   
     212             :   }
     213           0 :   LogIO logio;
     214             :   logio << "Clearing all model records in MS header."
     215           0 :           << LogIO::POST;
     216             : 
     217           0 :   MSColumns msc(thems);
     218           0 :   Vector<Int> fields=msc.fieldId().getColumn();
     219           0 :   Vector<String> fldnames=msc.field().name().getColumn();
     220           0 :   const Sort::Order order=Sort::Ascending;
     221           0 :   const Int option=Sort::HeapSort | Sort::NoDuplicates;
     222           0 :   Int nfields=GenSort<Int>::sort (fields, order, option);
     223           0 :   for (Int k=0; k< nfields; ++k){
     224           0 :     String elkey;
     225             :     Int srow;
     226             :     //if(newTab.rwKeywordSet().isDefined("definedmodel_field_"+String::toString(fields[k])))
     227           0 :       if(isModelDefined(fields[k], thems, elkey, srow))
     228             :       {
     229           0 :         logio << " " << fldnames[fields[k]] << " (id = " << fields[k] << ") deleted." << LogIO::POST;
     230             :         //Remove from Source table
     231           0 :         removeRecordByKey(newTab, elkey);
     232           0 :         if(srow > -1){
     233           0 :           if(thems.source().keywordSet().isDefined("definedmodel_field_"+String::toString(fields[k]))){
     234           0 :             newTab.source().rwKeywordSet().removeField("definedmodel_field_"+String::toString(fields[k]));
     235             :           }
     236             :         }
     237             :         //Remove from Main table if it is there
     238           0 :         if(newTab.rwKeywordSet().isDefined("definedmodel_field_"+String::toString(fields[k])))
     239           0 :           newTab.rwKeywordSet().removeField("definedmodel_field_"+String::toString(fields[k]));
     240             :       }
     241             :       else{
     242             :         ///remove keys that point to nowehere
     243             :         //cerr << "The key " << elkey << endl;
     244           0 :         if(!newTab.source().isNull() && newTab.source().rwKeywordSet().isDefined(elkey))
     245           0 :           newTab.source().rwKeywordSet().removeField(elkey);
     246           0 :         if(newTab.rwKeywordSet().isDefined(elkey))
     247           0 :           newTab.rwKeywordSet().removeField(elkey);
     248           0 :         if( !newTab.source().isNull() && newTab.source().rwKeywordSet().isDefined("definedmodel_field_"+String::toString(fields[k])))
     249           0 :           newTab.source().rwKeywordSet().removeField("definedmodel_field_"+String::toString(fields[k]));
     250           0 :         if(newTab.rwKeywordSet().isDefined("definedmodel_field_"+String::toString(fields[k])))
     251           0 :            newTab.rwKeywordSet().removeField("definedmodel_field_"+String::toString(fields[k]));
     252             :       }
     253           0 :   }
     254           0 :   if(!alreadyLocked)
     255           0 :     newTab.unlock();
     256           0 :   if(Table::isReadable(newTab.sourceTableName())){
     257           0 :     newTab.source().unlock();   
     258             :   }
     259             :   ////Cleaning out orphaned image disk models
     260           0 :   String srctable=thems.source().isNull() ? "" : thems.source().tableName();
     261           0 :   Vector<String> possibleFT(2); 
     262           0 :   possibleFT(0)=srctable+"/FT_MODEL*";
     263           0 :   possibleFT(1)=theParts[0]+"/FT_MODEL*";
     264             :   
     265           0 :   Vector<String> ftmods=Directory::shellExpand(possibleFT);
     266           0 :   for (uInt kk=0; kk < ftmods.nelements(); ++kk){
     267           0 :     if(TableUtil::canDeleteTable(ftmods[kk])){
     268           0 :       TableUtil::deleteTable(ftmods[kk]);
     269             :     }
     270             :   }
     271             : 
     272           0 : }
     273           0 :   void VisModelData::clearModel(const MeasurementSet& thems, const String field, const String specwindows){
     274           0 :     MeasurementSet& newTab=const_cast<MeasurementSet& >(thems);
     275           0 :     if(!newTab.isWritable())
     276           0 :       return;
     277           0 :   auto part_block = newTab.getPartNames(true);
     278           0 :   Vector<String> theParts(part_block.begin( ),part_block.end( ));
     279           0 :   if(theParts.nelements() > 1){
     280           0 :     for (uInt k =0; k < theParts.nelements(); ++k){
     281           0 :       MeasurementSet subms(theParts[k], newTab.lockOptions(), Table::Update);
     282           0 :       clearModel(subms, field, specwindows);
     283           0 :     }
     284           0 :     return;
     285             :   }
     286           0 :   if(!newTab.isWritable())
     287           0 :     return;
     288           0 :   Bool alreadyLocked=newTab.hasLock(True);
     289           0 :   if(!alreadyLocked)
     290           0 :     newTab.lock(True);
     291           0 :   if(Table::isReadable(newTab.sourceTableName())){
     292           0 :     if(!Table::isWritable(newTab.sourceTableName()))
     293           0 :       return;
     294           0 :     newTab.source().lock(True);   
     295             :   }
     296             : 
     297           0 :   MSColumns msc(thems);
     298           0 :   Vector<String> fldnames=msc.field().name().getColumn();
     299           0 :   Int nfields=0;
     300           0 :   Vector<Int> fields(0);
     301             :   {
     302             :     // Parse field specification
     303           0 :     MSSelection mss;
     304           0 :     mss.setFieldExpr(field);
     305           0 :     TableExprNode ten=mss.toTableExprNode(&thems);
     306           0 :     fields=mss.getFieldList();
     307           0 :     nfields=fields.nelements();
     308           0 :   }
     309             : 
     310             : 
     311           0 :   if (nfields==0)
     312             :     // Call the method that deletes them all
     313           0 :     VisModelData::clearModel(thems);
     314             :   else {
     315             : 
     316             :     //Now we have the two cases the whole field or specific spws
     317             :     // only delete the specified ones
     318           0 :     Vector<Int> spws(0);
     319           0 :     Int nspws=0;
     320             :     {
     321             :       // Parse field specification
     322           0 :       MSSelection mss;
     323           0 :       mss.setFieldExpr(field);
     324           0 :       mss.setSpwExpr(specwindows);
     325           0 :       TableExprNode ten=mss.toTableExprNode(&thems);
     326           0 :       spws=mss.getSpwList();
     327           0 :       nspws=spws.nelements();
     328           0 :     }
     329             : 
     330             :    
     331           0 :     LogIO logio;
     332             :     logio << "Clearing model records in MS header for selected fields." 
     333           0 :           << LogIO::POST;
     334             :     
     335           0 :     for (Int k=0; k< nfields; ++k){
     336           0 :       String elkey;
     337             :       Int srow;
     338           0 :       if(isModelDefined(fields[k], newTab, elkey, srow))
     339             :         
     340             :         {
     341             :         
     342           0 :           if(nspws==0){
     343           0 :             logio << " " << fldnames[fields[k]] << " (id = " << fields[k] << ") deleted." << LogIO::POST;
     344           0 :             removeRecordByKey(newTab, elkey);
     345           0 :             if(srow > -1 ){
     346           0 :               if(newTab.source().keywordSet().isDefined("definedmodel_field_"+String::toString(fields[k])))
     347           0 :                 newTab.source().rwKeywordSet().removeField("definedmodel_field_"+String::toString(fields[k]));              
     348             :             }
     349           0 :             if(newTab.rwKeywordSet().isDefined("definedmodel_field_"+String::toString(fields[k])))
     350           0 :               newTab.rwKeywordSet().removeField("definedmodel_field_"+String::toString(fields[k]));
     351             :           }
     352             :           else{
     353             :             //if(newTab.rwKeywordSet().isDefined(elkey)){
     354             :             
     355           0 :             TableRecord conteneur;  
     356           0 :             getModelRecord(elkey, conteneur, newTab); 
     357             :            
     358             :             //=newTab.rwKeywordSet().asRecord(elkey);
     359           0 :             if(removeSpw(conteneur, spws, fields)){
     360           0 :                 putRecordByKey(newTab, elkey, conteneur, srow);
     361             :             }
     362             :             else{
     363           0 :                 removeRecordByKey(newTab, elkey);
     364           0 :                 if(newTab.source().keywordSet().isDefined("definedmodel_field_"+String::toString(fields[k])))
     365           0 :                         newTab.source().rwKeywordSet().removeField("definedmodel_field_"+String::toString(fields[k]));
     366             :             }
     367             :               /*
     368             :               Vector<Int> defspws=conteneur.asArrayInt("spws");
     369             :               Vector<Int> newdefspw(defspws.nelements(), -1);
     370             :               Int counter=0;
     371             :               for(uInt k=0; k < defspws.nelements(); ++k){
     372             :                 for (Int j=0; j < nspws; ++j){
     373             :                   if(defspws[k] != spws[j]){
     374             :                     newdefspw[counter]=defspws[k];
     375             :                     ++counter;
     376             :                   }
     377             :                 }
     378             :               }
     379             :               if(counter==0){
     380             :                 //Now we have to remove this ftm or cft
     381             :                 newTab.rwKeywordSet().removeField(elkey);
     382             :               }
     383             :               else{
     384             :                 conteneur.define("spws", newdefspw);
     385             :                 updatespw(conteneur, 
     386             :               }
     387             :               */
     388           0 :           }
     389             :               
     390             :         
     391             :             
     392             :         }
     393             :       else
     394           0 :         logio << " " << fldnames[fields[k]] << " (id = " << fields[k] << ") not found." << LogIO::POST;
     395           0 :     }
     396             :     
     397           0 :   }
     398           0 :   if(!alreadyLocked)
     399           0 :     newTab.unlock();
     400           0 :   if(Table::isReadable(newTab.sourceTableName())){
     401           0 :     newTab.source().unlock();   
     402             :   }
     403             :   
     404             : 
     405           0 :   }
     406             : 
     407           0 :   Bool VisModelData::removeSpw(TableRecord& therec, const Vector<Int>& spws, const Vector<Int>& fields){
     408           0 :      Int numft=0;
     409           0 :      Vector<String> numtype(2);
     410           0 :      Vector<String> keyval(2);
     411           0 :      numtype[0]=String("numft");
     412           0 :      numtype[1]=String("numcl");
     413           0 :      keyval[0]="ft_";
     414           0 :      keyval[1]="cl_";
     415             :      
     416           0 :      for (Int j=0; j< 2; ++j){
     417           0 :          if(therec.isDefined(numtype[j])){
     418           0 :                  numft=therec.asInt(numtype[j]);
     419           0 :                  Vector<String> toberemoved(numft);
     420           0 :                  Int numrem=0;
     421           0 :                  for (Int k=0; k < numft; ++k){
     422           0 :                          RecordInterface& ftrec=therec.asrwRecord(keyval[j]+String::toString(k));
     423           0 :                          if(!ftrec.isDefined("version")){
     424           0 :                            Bool hasField=true;
     425           0 :                            if(fields.nelements() >0){
     426           0 :                              hasField=false;
     427           0 :                              Vector<Int> fieldsinrec;
     428           0 :                              ftrec.get("fields", fieldsinrec);
     429           0 :                              if(anyEQ(fieldsinrec, fields))
     430           0 :                                hasField=true;
     431           0 :                            }
     432           0 :                            if(hasField && !removeSpwFromMachineRec(ftrec, spws)){
     433           0 :                              toberemoved[numrem]=String(keyval[j]+String::toString(k));
     434           0 :                              ++numrem;
     435             :                            }
     436             :                          }
     437             :                          else{
     438             :                            //Version 2
     439           0 :                            const Matrix<Int>& ftcombind=ftrec.asArrayInt("indexcombination");
     440           0 :                            Matrix<Int> newComb;
     441           0 :                            for(uInt row=0; row < ftcombind.nrow(); ++row){
     442           0 :                              Vector<Int> rowcomb=ftcombind.row(row);      
     443           0 :                              Bool rowTobeRemoved=False;
     444           0 :                              for (uInt ff=0; ff < fields.nelements() ; ++ff){
     445           0 :                                for (uInt ss=0; ss < spws.nelements(); ++ss){
     446           0 :                                  if(rowcomb[0]==fields[ff] && rowcomb[1] == spws[ss] )
     447           0 :                                    rowTobeRemoved=True;
     448             :                                  
     449             :                                  
     450             :                                }
     451             :                              }
     452           0 :                              if(!rowTobeRemoved){
     453           0 :                                newComb.resize(newComb.nrow()+1,4, True);
     454           0 :                                newComb.row(newComb.nrow()-1)=rowcomb;
     455             :                              }
     456           0 :                            }
     457           0 :                            if(newComb.nrow()==0){
     458           0 :                              toberemoved[numrem]=String(keyval[j]+String::toString(k));
     459           0 :                              ++numrem;  
     460             :                            }
     461             :                            else{
     462           0 :                              ftrec.define("indexcombination", newComb);
     463             :                            }
     464             :                            
     465             :                            
     466           0 :                          }
     467             :                  }
     468           0 :                  if(numrem > 0){
     469           0 :                          for (Int k=0; k < numrem; ++k)
     470           0 :                                  removeFTFromRec(therec, toberemoved[k], k==(numrem-1));
     471             :                  }
     472           0 :          }
     473             :      }
     474           0 :      numft=0; Int numcl=0;
     475           0 :      if(therec.isDefined("numft")) numft=therec.asInt("numft");
     476           0 :      if(therec.isDefined("numcl")) numcl=therec.asInt("numcl");
     477           0 :      if (numft==0 && numcl==0)
     478           0 :          return false;
     479           0 :      return true;
     480           0 :   }
     481             : 
     482           0 :   Bool VisModelData::removeFTFromRec(TableRecord& therec, const String& keyval, const Bool relabel){
     483             : 
     484             : 
     485           0 :     String *splitkey=new String[2];
     486           0 :     Int nsep=split(keyval, splitkey, 2, String("_"));
     487           0 :     if (nsep <1 || !therec.isDefined(keyval)) 
     488           0 :       return false;
     489           0 :     String eltype=splitkey[0];
     490             :     //Int modInd=String::toInt(splitkey[1]);
     491           0 :     Int numcomp= (eltype==String("ft")) ? therec.asInt("numft"): therec.asInt("numcl");
     492           0 :     therec.removeField(keyval);
     493             :     
     494           0 :     numcomp=numcomp-1;
     495           0 :     if(eltype=="ft")
     496           0 :       therec.define("numft", numcomp);
     497             :     else
     498           0 :       therec.define("numcl", numcomp);
     499           0 :     if(relabel){
     500             : 
     501             :     
     502           0 :         eltype=eltype+String("_");
     503           0 :         Int id=0;
     504           0 :         for(uInt k=0; k < therec.nfields(); ++k){
     505           0 :                 if(therec.name(k).contains(eltype)){
     506           0 :                         therec.renameField(eltype+String::toString(id), k);
     507           0 :                         ++id;
     508             :                 }
     509             :         }
     510             :     }
     511             : 
     512           0 :     delete [] splitkey;
     513           0 :     return true;
     514           0 :   } 
     515             : 
     516           0 :   Bool VisModelData::removeSpwFromMachineRec(RecordInterface& ftclrec, const Vector<Int>& spws){
     517           0 :     Vector<Int> defspws=ftclrec.asArrayInt("spws");
     518           0 :     Vector<Int> newdefspw(defspws.nelements(), -1);
     519           0 :     Int counter=0;
     520           0 :     for(uInt k=0; k < defspws.nelements(); ++k){
     521           0 :       for (uInt j=0; j < spws.nelements(); ++j){
     522           0 :         if(defspws[k] == spws[j]){
     523           0 :           defspws[k]=-1;
     524           0 :           ++counter;
     525             :         }
     526             :       }
     527             :     }
     528           0 :     if(defspws.nelements() == uInt(counter)){
     529             :       //Now we have to remove this ft or cl model
     530             : 
     531           0 :       return false;
     532             : 
     533             :     }    
     534           0 :     newdefspw.resize(defspws.nelements()-counter);
     535           0 :     counter=0;
     536           0 :     for (uInt k=0; k < defspws.nelements(); ++k){
     537           0 :       if(defspws[k] != -1){
     538           0 :         newdefspw[counter]=defspws[k];
     539           0 :         ++counter;
     540             :       }
     541             :     }
     542             :     
     543           0 :     ftclrec.define("spws", newdefspw);
     544           0 :     return true;
     545           0 :   }
     546             : 
     547           0 :   Bool VisModelData::addToRec(TableRecord& therec, const Vector<Int>& spws){
     548             : 
     549           0 :     Int numft=0;
     550           0 :     Int numcl=0;
     551           0 :     Vector<Bool> hasSpw(spws.nelements(), false);
     552           0 :     if(therec.isDefined("numft")){
     553           0 :       numft=therec.asInt("numft");
     554           0 :       Vector<Int> ft_toremove(numft, 0);
     555           0 :       for(Int k=0; k < numft; ++k){
     556           0 :         const Record& ftrec=therec.asRecord("ft_"+String::toString(k));
     557           0 :         const Vector<Int>& ftspws=ftrec.asArrayInt("spws");
     558           0 :         for (uInt i=0; i<spws.nelements(); ++i){
     559           0 :           for (uInt j=0; j<ftspws.nelements(); ++j){
     560           0 :             if(spws[i]==ftspws[j]){
     561           0 :               hasSpw[i]=true;      
     562           0 :               ft_toremove[k]=1;
     563             :             }
     564             :           }     
     565             :         }
     566           0 :       }
     567           0 :       if(sum(ft_toremove) >0){
     568           0 :         for(Int k=0; k < numft; ++k){
     569           0 :           if(ft_toremove[k]==1)
     570           0 :             therec.removeField("ft_"+String::toString(k));
     571             :         }
     572           0 :         numft=numft-sum(ft_toremove);
     573           0 :         therec.define("numft", numft);
     574           0 :         Int id=0;
     575           0 :         for(uInt k=0; k < therec.nfields(); ++k){
     576           0 :           if(therec.name(k).contains("ft_")){
     577           0 :             therec.renameField("ft_"+String::toString(id), k);
     578           0 :             ++id;
     579             :           }
     580             :         }
     581             :       }
     582           0 :     }
     583           0 :     if(therec.isDefined("numcl")){
     584           0 :       numcl=therec.asInt("numcl");
     585           0 :       Vector<Int> cl_toremove(numcl, 0);
     586           0 :       for(Int k=0; k < numcl; ++k){
     587           0 :         const Record& clrec=therec.asRecord("cl_"+String::toString(k));
     588           0 :         const Vector<Int>& clspws=clrec.asArrayInt("spws");
     589           0 :         for (uInt i=0; i<spws.nelements(); ++i){
     590           0 :           for (uInt j=0; j<clspws.nelements(); ++j){
     591           0 :             if(spws[i]==clspws[j]){
     592           0 :               hasSpw[i]=true;       
     593           0 :               cl_toremove[k]=1;
     594             :             }
     595             :           }     
     596             :         }
     597           0 :       }
     598           0 :       if(sum(cl_toremove) >0){
     599           0 :         for(Int k=0; k < numcl; ++k){
     600           0 :           if(cl_toremove[k]==1)
     601           0 :             therec.removeField("cl_"+String::toString(k));
     602             :         }
     603           0 :         numcl=numcl-sum(cl_toremove);
     604           0 :         therec.define("numcl", numcl);
     605           0 :         Int id=0;
     606           0 :         for(uInt k=0; k < therec.nfields(); ++k){
     607           0 :           if(therec.name(k).contains("cl_")){
     608           0 :             therec.renameField("cl_"+String::toString(id), k);
     609           0 :             ++id;
     610             :           }
     611             :         }
     612             :       }
     613           0 :     }
     614           0 :     return (!allTrue(hasSpw) || ((numft+numcl)>0));
     615           0 :   }
     616             : 
     617           0 : Bool VisModelData::addToRec(TableRecord& therec, const Matrix<Int>& combind){
     618           0 :   if(therec.nfields() >0 && ((!therec.isDefined("version")) || (therec.asInt("version") !=2)))
     619           0 :     throw(AipsError("cannot combine with older version virtual model column;\n please clear it before proceeding"));
     620           0 :   Int numft=0;
     621           0 :   Int numcl=0;
     622           0 :   Vector<Bool> hasComb(combind.nrow(), false);
     623           0 :   if(therec.isDefined("numft")){
     624           0 :     numft=therec.asInt("numft");
     625           0 :     Vector<Int> ft_toremove(numft, 0);
     626           0 :     for(Int k=0; k < numft; ++k){
     627           0 :       const Record& ftrec=therec.asRecord("ft_"+String::toString(k));
     628           0 :       const Matrix<Int>& ftcombind=ftrec.asArrayInt("indexcombination");
     629           0 :       for (uInt i=0; i<combind.nrow(); ++i){
     630           0 :         for (uInt j=0; j<ftcombind.nrow(); ++j){
     631           0 :           if(allEQ(combind.row(i), ftcombind.row(j))){
     632             :             //We could upgrade this to just delete the intent from the ft
     633           0 :             hasComb[i]=true;       
     634           0 :             ft_toremove[k]=1;
     635             :           }
     636             :         }       
     637             :       }
     638           0 :     }
     639           0 :     if(sum(ft_toremove) >0){
     640           0 :       for(Int k=0; k < numft; ++k){
     641           0 :         if(ft_toremove[k]==1)
     642           0 :           therec.removeField("ft_"+String::toString(k));
     643             :       }
     644           0 :       numft=numft-sum(ft_toremove);
     645           0 :       therec.define("numft", numft);
     646           0 :       Int id=0;
     647           0 :       for(uInt k=0; k < therec.nfields(); ++k){
     648           0 :         if(therec.name(k).contains("ft_")){
     649           0 :           therec.renameField("ft_"+String::toString(id), k);
     650           0 :           ++id;
     651             :         }
     652             :       }
     653             :     }
     654           0 :   }
     655           0 :   if(therec.isDefined("numcl")){
     656           0 :     numcl=therec.asInt("numcl");
     657           0 :     Vector<Int> cl_toremove(numcl, 0);
     658           0 :     for(Int k=0; k < numcl; ++k){
     659           0 :       const Record& clrec=therec.asRecord("cl_"+String::toString(k));
     660           0 :       const Matrix<Int>& clcombind=clrec.asArrayInt("indexcombination");
     661           0 :       for (uInt i=0; i<combind.nrow(); ++i){
     662           0 :         for (uInt j=0; j<clcombind.nrow(); ++j){
     663           0 :           if(allEQ(combind.row(i),clcombind.row(j))){
     664           0 :             hasComb[i]=true;        
     665           0 :             cl_toremove[k]=1;
     666             :           }
     667             :         }       
     668             :       }
     669           0 :       }
     670           0 :     if(sum(cl_toremove) >0){
     671           0 :       for(Int k=0; k < numcl; ++k){
     672           0 :         if(cl_toremove[k]==1)
     673           0 :           therec.removeField("cl_"+String::toString(k));
     674             :       }
     675           0 :       numcl=numcl-sum(cl_toremove);
     676           0 :       therec.define("numcl", numcl);
     677           0 :       Int id=0;
     678           0 :       for(uInt k=0; k < therec.nfields(); ++k){
     679           0 :         if(therec.name(k).contains("cl_")){
     680           0 :           therec.renameField("cl_"+String::toString(id), k);
     681           0 :           ++id;
     682             :         }
     683             :       }
     684             :     }
     685           0 :   }
     686           0 :   return (!allTrue(hasComb) || ((numft+numcl)>0));
     687           0 : }
     688             :  
     689             :   
     690        1200 : Bool VisModelData::isModelDefined(const Int fieldId, const MeasurementSet& thems, String& thekey, Int& sourceRow){
     691        1200 :   sourceRow=-1;
     692        2400 :   String modelkey=String("definedmodel_field_")+String::toString(fieldId);
     693        1200 :   thekey="";
     694        1200 :   if(Table::isReadable(thems.sourceTableName()) && thems.source().nrow() > 0 && (thems.source().keywordSet().isDefined(modelkey))){
     695             :       {
     696           0 :       thekey=thems.source().keywordSet().asString(modelkey);
     697           0 :       if(thems.source().keywordSet().isDefined(thekey))
     698           0 :         sourceRow=thems.source().keywordSet().asInt(thekey);
     699             :     }
     700             :   }else{
     701        1200 :     if(thems.keywordSet().isDefined(modelkey)){
     702           0 :       thekey=thems.keywordSet().asString(modelkey);  
     703             :     }
     704             :   }
     705        1200 :   if(thekey != "" )
     706           0 :     return VisModelData::isModelDefined(thekey, thems);
     707        1200 :   return false;
     708             :   
     709        1200 : }
     710             : 
     711           0 :   Bool VisModelData::isModelDefined(const String& elkey, const MeasurementSet& thems){
     712             :     //Let's try the Source table
     713           0 :     if(Table::isReadable(thems.sourceTableName()) &&thems.source().nrow() > 0 ){
     714           0 :       if(thems.source().keywordSet().isDefined(elkey))
     715           0 :         return true;      
     716             :     }
     717             :     //Let's try the Main table 
     718           0 :     if(thems.keywordSet().isDefined(elkey))
     719           0 :       return true;
     720           0 :     return false;
     721             :   }
     722             : 
     723           0 :   Bool VisModelData::getModelRecord(const String& theKey, TableRecord& theRec, const MeasurementSet& theMs){
     724             :     //Let's try the Source table
     725           0 :     if(Table::isReadable(theMs.sourceTableName()) &&theMs.source().nrow() > 0 && (theMs.source().keywordSet().isDefined(theKey))){
     726             :         {
     727             :         //Get the row for the model 
     728           0 :         Int row=theMs.source().keywordSet().asInt(theKey);
     729             :         //MSSourceColumns srcCol(theMs.source());
     730             :      
     731           0 :         ScalarColumn<TableRecord> scol(theMs.source(), "SOURCE_MODEL");
     732           0 :         scol.get(row, theRec);
     733           0 :       }
     734           0 :       return true;
     735             :     }
     736             :     //Let's try the Main table 
     737           0 :     if(theMs.keywordSet().isDefined(theKey)){
     738           0 :       theRec=theMs.keywordSet().asRecord(theKey);
     739           0 :       return true;
     740             :     }
     741           0 :     return false;
     742             : 
     743             : 
     744             :   }
     745             : 
     746           0 :   void VisModelData::putRecordByKey(MeasurementSet& theMS, const String& theKey, const TableRecord& theRec, const Int sourceRowNum){
     747             :     //Prefer the Source table first    
     748           0 :     if( (sourceRowNum> -1) && Table::isReadable(theMS.sourceTableName()) && Int(theMS.source().nrow()) > sourceRowNum ){
     749           0 :       MSSource& mss=theMS.source();
     750           0 :       if(!mss.isColumn(MSSource::SOURCE_MODEL) ){
     751           0 :         mss.addColumn(ScalarRecordColumnDesc("SOURCE_MODEL"), true);
     752             :       }
     753           0 :       if(mss.rwKeywordSet().isDefined(theKey))
     754           0 :         mss.rwKeywordSet().removeField(theKey);
     755           0 :       mss.rwKeywordSet().define(theKey, sourceRowNum);
     756           0 :       MSSourceColumns srcCol(mss);
     757           0 :       srcCol.sourceModel().put(sourceRowNum, theRec);
     758           0 :       return;
     759           0 :     }
     760             :     //Oh well no source table so will add it to the main table
     761           0 :     theMS.rwKeywordSet().defineRecord(theKey, theRec);
     762             :     
     763             : 
     764             :   }
     765             : 
     766           0 :   Bool VisModelData::putModelRecord(const Vector<Int>& fieldIds, TableRecord& theRec, MeasurementSet& theMS){
     767           0 :     auto part_block = theMS.getPartNames(true);
     768           0 :     Vector<String> theParts(part_block.begin( ),part_block.end( ));
     769           0 :     if(theParts.nelements() > 1){
     770           0 :       Bool retval=true;
     771           0 :       for (uInt k =0; k < theParts.nelements(); ++k){
     772           0 :           MeasurementSet subms(theParts[k], theMS.lockOptions(), Table::Update);
     773           0 :           retval= retval && putModelRecord(fieldIds, theRec, subms);
     774           0 :       }
     775           0 :       return retval;
     776             :     }
     777             :  
     778           0 :     String elkey="model";
     779           0 :     for (uInt k=0; k < fieldIds.nelements();  ++k){
     780           0 :       elkey=elkey+"_"+String::toString(fieldIds[k]);
     781             :     }
     782           0 :     if(theMS.rwKeywordSet().isDefined(elkey))
     783           0 :       theMS.rwKeywordSet().removeField(elkey);
     784           0 :     Int row=-1;
     785             :     //Prefer the Source table first    
     786           0 :     MSFieldColumns fCol(theMS.field());
     787           0 :     if(Table::isReadable(theMS.sourceTableName()) && (theMS.source().nrow() > 0) &&  (!fCol.sourceId().isNull()) && (fCol.sourceId().get(fieldIds[0]) != -1) ){
     788             :       //
     789           0 :       row=0;
     790           0 :       MSSource& mss=theMS.source();
     791             :      
     792           0 :       Int sid=fCol.sourceId().get(fieldIds[0]);
     793           0 :       auto rows=MSSourceIndex(mss).getRowNumbersOfSourceID(sid);
     794           0 :       if(rows.nelements() > 0) 
     795           0 :         row=rows[0];
     796             :       else{
     797           0 :         LogIO logio;
     798           0 :         logio << "Invalid Source_id "+String::toString(sid)+" found in FIELD table\n" 
     799             :               <<"Model is being written at Source ID 0 position which will be erased\n" 
     800             :               << "Fix the FIELD table before proceeding " 
     801           0 :               <<  LogIO::WARN << LogIO::POST;
     802             :         
     803           0 :       }
     804             :           
     805             :         
     806             :       
     807           0 :       putRecordByKey(theMS, elkey, theRec, row);
     808           0 :       for (uInt k=0; k < fieldIds.nelements();  ++k){
     809           0 :         mss.rwKeywordSet().define("definedmodel_field_"+String::toString(fieldIds[k]), elkey);
     810             :       }
     811           0 :       return true;
     812             :       
     813           0 :     }
     814             :     //Oh well no source table so will add it to the main table
     815           0 :     putRecordByKey(theMS, elkey, theRec, -1);
     816           0 :     for (uInt k=0; k < fieldIds.nelements();  ++k){
     817           0 :       theMS.rwKeywordSet().define("definedmodel_field_"+String::toString(fieldIds[k]), elkey);        
     818             :     }
     819             : 
     820           0 :     return true;
     821           0 :   }
     822             : 
     823             : 
     824           0 : void VisModelData::putModel(const MeasurementSet& thems, const RecordInterface& rec, const Vector<Int>& validfieldids, const Vector<Int>& spws, const Vector<Int>& starts, const Vector<Int>& nchan,  const Vector<Int>& incr, Bool iscomponentlist, Bool incremental){
     825             : 
     826           0 :   LogIO logio;
     827             : 
     828             :   try{
     829             :     //A field can have multiple FTmachines and ComponentList associated with it 
     830             :     //For example having many flanking images for the model
     831             :     //For componentlist it may have multiple componentlist ...for different spw
     832             :   //Timer tim;
     833             :   //tim.mark();
     834           0 :     Int counter=0;
     835           0 :     Record modrec;
     836           0 :     modrec.define("fields", validfieldids);
     837           0 :     modrec.define("spws", spws);
     838           0 :     modrec.define("start", starts);
     839           0 :     modrec.define("nchan", nchan);
     840           0 :     modrec.define("incr", incr);
     841           0 :     modrec.defineRecord("container", rec);
     842           0 :     String elkey="model";
     843           0 :     for (uInt k=0; k < validfieldids.nelements();  ++k){
     844           0 :       elkey=elkey+"_"+String::toString(validfieldids[k]);
     845             :     }
     846           0 :     TableRecord outRec; 
     847           0 :     Bool addtorec=false;
     848           0 :     MeasurementSet& newTab=const_cast<MeasurementSet& >(thems);
     849           0 :     newTab.lock(True);
     850           0 :     if(Table::isReadable(newTab.sourceTableName())){
     851           0 :       newTab.source().lock(True);   
     852             :     }
     853           0 :     if(isModelDefined(elkey, newTab)){ 
     854           0 :       getModelRecord(elkey, outRec, thems);
     855             :       //if incremental no need to check & remove what is in the record
     856           0 :       if(!incremental)
     857           0 :           addtorec=addToRec(outRec, spws);
     858             :     }
     859             : 
     860             :     ///////even if it is not defined some other field model might be sitting on that
     861             :     //////model key
     862           0 :     Int hasSourceRecord=firstSourceRowRecord(validfieldids[0], thems, outRec);
     863           0 :     if(hasSourceRecord > -1 && outRec.nfields() > 0)
     864           0 :       addtorec=true;
     865             :     //cerr << "has Source " << hasSourceRecord << endl;
     866             :     ////
     867             : 
     868           0 :     incremental=incremental || addtorec;
     869           0 :     if(iscomponentlist){
     870           0 :       modrec.define("type", "componentlist");
     871           0 :       if(outRec.isDefined("numcl"))
     872           0 :         counter=incremental ? outRec.asInt("numcl") : 0;
     873             :             
     874             :     }
     875             :     else{
     876           0 :       modrec.define("type", "ftmachine");
     877           0 :       if(outRec.isDefined("numft"))
     878           0 :         counter=incremental ? outRec.asInt("numft") : 0;
     879             :     }
     880           0 :     iscomponentlist ? outRec.define("numcl", counter+1) : outRec.define("numft", counter+1); 
     881             :   
     882             :     //for (uInt k=0; k < validfieldids.nelements();  ++k){
     883             :     //  newTab.rwKeywordSet().define("definedmodel_field_"+String::toString(validfieldids[k]), elkey);
     884             :       
     885             :     // }
     886           0 :     iscomponentlist ? outRec.defineRecord("cl_"+String::toString(counter), modrec):
     887           0 :       outRec.defineRecord("ft_"+String::toString(counter), modrec);
     888             :     //////////////////;
     889             :     //for (uInt k=0; k < newTab.rwKeywordSet().nfields() ; ++k){
     890             :     //  cerr << "keys " << k << "  is  " << newTab.rwKeywordSet().name(k) << " type " << newTab.rwKeywordSet().dataType(k) << endl;
     891             :     //}
     892             :     ////////////////////////
     893             :     //if image for a given key is on disk and not incrementing ...lets remove it
     894           0 :     if(!incremental) 
     895           0 :       deleteDiskImage(newTab, elkey);
     896           0 :     putModelRecord(validfieldids, outRec, newTab);  
     897             :     //if(newTab.rwKeywordSet().isDefined(elkey))
     898             :     //  newTab.rwKeywordSet().removeField(elkey);
     899             :     //newTab.rwKeywordSet().defineRecord(elkey, outRec);
     900             :     
     901             :     // tim.show("Time taken to save record ");
     902             :  
     903             :       /*
     904             :       String subName=newTab.tableName()+"/"+elkey;
     905             :       if(Table::isReadable(subName) && !Table::canDeleteTable(subName, true))
     906             :         throw(AipsError("Cannot save model into MS"));
     907             :       else if ((Table::isReadable(subName) && Table::canDeleteTable(subName, true)))
     908             :         Table::deleteTable(subName, true);
     909             :       TableDesc td1 ("td1", TableDesc::New);
     910             :       //td1.addColumn (ArrayColumnDesc<Int> ("SPECTRAL_WINDOW_ID"));
     911             :       //td1.addColumn(ScalarColumnDesc<TableRecord>("MODEL"));
     912             :       SetupNewTable newtab1 (subName, td1, Table::New);
     913             :       Table tab1 (newtab1);
     914             :       newTab.rwKeywordSet().defineTable(elkey, tab1);
     915             :       tab1.rwKeywordSet().defineRecord(elkey, outRec);
     916             :       */
     917             :       //ArrayColumn<Int> spwCol(tab1, "SPECTRAL_WINDOW_ID");
     918             :       //ScalarColumn<TableRecord> modCol(tab1, "MODEL");
     919             :       //newTab.addRow(1,false);
     920             :       //spwCol.put(0,spws);
     921             :       //modCol.put(0,outRec);
     922             :       //newTab.flush();
     923             :     // MSSource& mss=newTab.source();
     924             :     //  cerr << "has model_source " << mss.isColumn(MSSource::SOURCE_MODEL) << endl;
     925             :     //  if(!mss.isColumn(MSSource::SOURCE_MODEL) ){
     926             :     //mss.addColumn(ScalarRecordColumnDesc("SOURCE_MODEL"), true);
     927             :     //  }
     928             :     //  MSSourceColumns srcCol(mss);
     929             :     //  srcCol.sourceModel().put(0, outRec);
     930           0 :     newTab.unlock();
     931           0 :     if(Table::isReadable(newTab.sourceTableName())){
     932           0 :       newTab.source().unlock();   
     933             :     }
     934           0 :   }
     935           0 :   catch(...){
     936           0 :     logio << "Could not save virtual model data column due to an artificial virtual model size limit. \nYou may need to use the scratch column if you need model visibilities" << LogIO::WARN << LogIO::POST ;
     937           0 :    const_cast<MeasurementSet& >(thems).unlock(); 
     938           0 :   }
     939             : 
     940           0 : }
     941           0 :   void VisModelData::putModel(const MeasurementSet& thems,const RecordInterface& rec, const Matrix<Int>& indexComb, const Matrix<Int>& chanSel,Bool iscomponentlist, Bool incremental){
     942           0 :     LogIO logio;
     943             :     
     944             :     try{
     945             :       //A field can have multiple FTmachines and ComponentList associated with it 
     946             :       //For example having many flanking images for the model
     947             :       //For componentlist it may have multiple componentlist ...for different spw
     948             :       //Timer tim;
     949             :       //tim.mark();
     950           0 :       Int counter=0;
     951           0 :       Record modrec;
     952           0 :       modrec.define("version", Int(2));
     953           0 :       Vector<Int> validfieldids;
     954           0 :       validfieldids=indexComb.column(0);
     955           0 :       Int nfields=GenSort<Int>::sort (validfieldids, Sort::Ascending, 
     956           0 :                                       Sort::QuickSort|Sort::NoDuplicates);
     957           0 :       validfieldids.resize(nfields, True);
     958           0 :       modrec.define("indexcombination", indexComb);
     959             :       //cerr << "modrec after indexcombination " << modrec << endl;
     960           0 :       modrec.define("channelselection", chanSel);
     961           0 :       modrec.defineRecord("container", rec);
     962           0 :       String elkey="model";
     963           0 :       for (uInt k=0; k < validfieldids.nelements();  ++k){
     964           0 :         elkey=elkey+"_"+String::toString(validfieldids[k]);
     965             :       }
     966           0 :       TableRecord outRec; 
     967           0 :       Bool addtorec=false;
     968           0 :       MeasurementSet& newTab=const_cast<MeasurementSet& >(thems);
     969           0 :       newTab.lock(True);
     970           0 :       if(Table::isReadable(newTab.sourceTableName())){
     971           0 :         newTab.source().lock(True);   
     972             :       }
     973             :       //cerr << elkey << " incr " << incremental << endl;
     974           0 :       if(isModelDefined(elkey, newTab)){ 
     975           0 :         getModelRecord(elkey, outRec, thems);
     976             :         //if incremental no need to check & remove what is in the record
     977           0 :         if(!incremental)
     978           0 :           addtorec=addToRec(outRec, indexComb);
     979             :         //cerr << "addToRec " << addtorec << endl;
     980             :       }
     981             :       else{
     982             :         ///////even if it is not defined some other field model might be sitting on that
     983             :         //////model key
     984           0 :         Int hasSourceRecord=firstSourceRowRecord(validfieldids[0], thems, outRec);
     985           0 :         if(hasSourceRecord > -1 && outRec.nfields() > 0)
     986           0 :           addtorec=true;
     987             :         //cerr << "has Source " << hasSourceRecord << " addToRec " << addtorec << endl;
     988             :         ////
     989             :       }
     990           0 :       if(outRec.nfields() >0 && ((!outRec.isDefined("version")) || (outRec.asInt("version") !=2)))
     991           0 :         throw(AipsError("cannot combine with older version virtual model column;\n please clear it before proceeding"));
     992           0 :       outRec.define("version", Int(2));
     993           0 :       incremental=incremental || addtorec;
     994           0 :       if(iscomponentlist){
     995           0 :         modrec.define("type", "componentlist");
     996           0 :         if(outRec.isDefined("numcl"))
     997           0 :           counter=incremental ? outRec.asInt("numcl") : 0;
     998             :         
     999             :       }
    1000             :       else{
    1001           0 :         modrec.define("type", "ftmachine");
    1002           0 :         if(outRec.isDefined("numft"))
    1003           0 :           counter=incremental ? outRec.asInt("numft") : 0;
    1004             :       }
    1005           0 :       iscomponentlist ? outRec.define("numcl", counter+1) : outRec.define("numft", counter+1); 
    1006             :       
    1007             :       //for (uInt k=0; k < validfieldids.nelements();  ++k){
    1008             :       //  newTab.rwKeywordSet().define("definedmodel_field_"+String::toString(validfieldids[k]), elkey);
    1009             :       
    1010             :       // }
    1011           0 :       iscomponentlist ? outRec.defineRecord("cl_"+String::toString(counter), modrec):
    1012           0 :         outRec.defineRecord("ft_"+String::toString(counter), modrec);
    1013             :       //////////////////;
    1014             :       //for (uInt k=0; k < newTab.rwKeywordSet().nfields() ; ++k){
    1015             :       //  cerr << "keys " << k << "  is  " << newTab.rwKeywordSet().name(k) << " type " << newTab.rwKeywordSet().dataType(k) << endl;
    1016             :       //}
    1017             :       ////////////////////////
    1018             :       //if image for a given key is on disk and not incrementing ...lets remove it
    1019           0 :       if(!incremental) 
    1020           0 :         deleteDiskImage(newTab, elkey);
    1021           0 :       putModelRecord(validfieldids, outRec, newTab);  
    1022           0 :       newTab.unlock();
    1023           0 :       if(Table::isReadable(newTab.sourceTableName())){
    1024           0 :         newTab.source().unlock();   
    1025             :     }
    1026           0 :     }
    1027           0 :     catch(...){
    1028           0 :       logio << "Could not save virtual model data for some reason \nYou may need clear the model and redo or  use the scratch column if you need model visibilities" << LogIO::WARN << LogIO::POST ;
    1029           0 :       const_cast<MeasurementSet& >(thems).unlock(); 
    1030           0 :     }
    1031             :     
    1032           0 :   }
    1033             : 
    1034           0 :   void VisModelData::modifyDiskImagePath(Record& rec, const VisBuffer2& vb){
    1035           0 :   Record& ftmac= rec.rwSubRecord("container");
    1036           0 :   if(ftmac.isDefined("diskimage")){
    1037           0 :     String theDiskImage=ftmac.asString("diskimage");
    1038           0 :     if(File(theDiskImage).exists())
    1039           0 :       return;
    1040           0 :     String subPathname[30];
    1041           0 :     String sep = "/";
    1042           0 :     uInt nw = split(theDiskImage, subPathname, 20, sep);
    1043           0 :     String theposs=(subPathname[nw-1]);
    1044             :     ///String msname=vb.msName(false);
    1045             :     ///As the above does not work ..returns some bogus reference table
    1046           0 :     String msname=(vb.ms()).antenna().tableName();
    1047           0 :     msname.erase(msname.length()-8);
    1048           0 :     for (uInt i=nw-2 ; i>0; --i){ 
    1049           0 :       if(File(msname+"/"+theposs).exists()){
    1050           0 :                 theDiskImage=msname+"/"+theposs;
    1051           0 :                 ftmac.define("diskmage", theDiskImage);
    1052           0 :                 return;
    1053             :       }
    1054           0 :       theposs=subPathname[i]+"/"+theposs;
    1055             :     }
    1056           0 :   }
    1057             : }
    1058             : 
    1059             : 
    1060             : 
    1061           0 :   void VisModelData::addModel(const RecordInterface& rec,  const Vector<Int>& /*msids*/, const VisBuffer2& vb){
    1062             :     
    1063             : 
    1064             : 
    1065           0 :     Int indexft=-1;
    1066           0 :     if(rec.isDefined("numft")){
    1067           0 :       Int numft=rec.asInt("numft");
    1068           0 :       if(numft >0){
    1069           0 :         for(Int ftk=0; ftk < numft; ++ftk){
    1070           0 :           Record ftrec(rec.asRecord("ft_"+String::toString(ftk)));
    1071           0 :           modifyDiskImagePath(ftrec, vb);
    1072           0 :           Vector<Int>fields;
    1073           0 :           Vector<Int> spws;
    1074           0 :           ftrec.get("fields", fields);
    1075           0 :           ftrec.get("spws", spws);
    1076           0 :           if(anyEQ(spws, vb.spectralWindows()(0))){
    1077           0 :             indexft=ftholder_p.nelements();
    1078           0 :             ftholder_p.resize(indexft+1, false, true);
    1079           0 :             ftholder_p[indexft].resize(1);
    1080           0 :             ftholder_p[indexft][0]=NEW_FT(ftrec.asRecord("container"));
    1081           0 :             if(!(ftholder_p[indexft][0]))
    1082           0 :               throw(AipsError("Unknown FTMachine saved in virtual model"));
    1083           0 :             ftholder_p[indexft][0]->initMaps(vb);
    1084             :             
    1085           0 :             for( uInt fi=0; fi < fields.nelements(); ++fi){
    1086           0 :               for(uInt spi=0; spi < spws.nelements(); ++spi){
    1087           0 :                 Int indx=-1;
    1088           0 :                 Int ftindx=-1;
    1089           0 :                 if(hasModel(vb.msId(), fields[fi], spws[spi]) && (ftindex_p(spws[spi], fields[fi], vb.msId()) >= 0 )){
    1090             :                  
    1091           0 :                   indx=ftindex_p(spws[spi], fields[fi], vb.msId());
    1092           0 :                   ftindx=ftholder_p[indx].nelements();
    1093           0 :                   Bool alreadyAdded=false;
    1094           0 :                   for (Int kk=1; kk < ftindx; ++kk){
    1095           0 :                     alreadyAdded= alreadyAdded || (ftholder_p[indexft][0]==ftholder_p[indx][kk]);
    1096             :                   }
    1097           0 :                   if(!alreadyAdded){
    1098           0 :                     ftholder_p[indx].resize(ftindx+1, true);
    1099           0 :                     ftholder_p[indx][ftindx]=ftholder_p[indexft][0];
    1100             :                   }
    1101             :                 }
    1102             :                 else{
    1103           0 :                   ftindex_p(spws[spi], fields[fi], vb.msId())=indexft;
    1104             :                 }
    1105             :               }
    1106             :             }
    1107             :           }
    1108             :           else{
    1109           0 :             if(hasModel(vb.msId(), vb.fieldId()(0), vb.spectralWindows()(0)) < 0)
    1110           0 :               ftindex_p(vb.spectralWindows()(0), vb.fieldId()(0), vb.msId())=-2;
    1111             :           }
    1112             : 
    1113             :           
    1114           0 :         }
    1115             :       }       
    1116             :     }
    1117           0 :     Int indexcl=-1;
    1118           0 :     if(rec.isDefined("numcl")){
    1119           0 :       Int numcl=rec.asInt("numcl");
    1120           0 :       if(numcl >0){
    1121           0 :         for(Int clk=0; clk < numcl; ++clk){
    1122           0 :           Vector<Int>fields;
    1123           0 :           Vector<Int> spws;
    1124           0 :           Record clrec(rec.asRecord("cl_"+String::toString(clk)));
    1125           0 :           clrec.get("fields", fields);
    1126           0 :           clrec.get("spws", spws);
    1127           0 :           if(anyEQ(spws, vb.spectralWindows()(0))){
    1128           0 :             indexcl=clholder_p.nelements();
    1129           0 :             clholder_p.resize(indexcl+1, false, true);
    1130           0 :             clholder_p[indexcl].resize(1);
    1131           0 :             clholder_p[indexcl][0]=new ComponentList();
    1132           0 :             String err;
    1133           0 :             if(!((clholder_p[indexcl][0])->fromRecord(err, clrec.asRecord("container"))))
    1134           0 :               throw(AipsError("Component model failed to load for field "+String::toString(fields)));
    1135           0 :             for( uInt fi=0; fi < fields.nelements(); ++fi){
    1136           0 :               for(uInt spi=0; spi < spws.nelements(); ++spi){
    1137           0 :                 Int indx=-1;
    1138           0 :                 Int clindx=-1;
    1139           0 :                 if(hasModel(vb.msId(), fields[fi], spws[spi]) && (clindex_p(spws[spi], fields[fi], vb.msId()) >= 0 )){
    1140           0 :                   indx=clindex_p(spws[spi], fields[fi], vb.msId());
    1141           0 :                   clindx=clholder_p[indx].nelements();
    1142           0 :                   Bool alreadyAdded=false;
    1143           0 :                   for (Int kk=1; kk < clindx; ++kk){
    1144           0 :                     alreadyAdded= alreadyAdded || (clholder_p[indexcl][0]==clholder_p[indx][kk]);
    1145             :                   } 
    1146           0 :                   if(!alreadyAdded){
    1147           0 :                     clholder_p[indx].resize(clindx+1, true);
    1148           0 :                     clholder_p[indx][clindx]=clholder_p[indexcl][0];
    1149             :                   }
    1150             :                 }
    1151             :                 else{
    1152             :                   //              cerr << "spws, field, msid " << spws[spi] << " " << fields[fi] << " " <<  vb.msId() << " index " << indexcl << endl;
    1153           0 :                   clindex_p(spws[spi], fields[fi], vb.msId())=indexcl;
    1154             :                 }
    1155             :               }
    1156             :             }
    1157           0 :           }
    1158             :           else{
    1159           0 :             if(hasModel(vb.msId(), vb.fieldId()(0), vb.spectralWindows()(0)) < 0)
    1160           0 :               clindex_p(vb.spectralWindows()(0), vb.fieldId()(0), vb.msId())=-2;
    1161             :           }
    1162             : 
    1163           0 :         }
    1164             :       }
    1165             :     }
    1166             : 
    1167             : 
    1168           0 :   }
    1169             : 
    1170           0 :   FTMachine* VisModelData::NEW_FT(const Record& ftrec){
    1171           0 :     String name=ftrec.asString("name");
    1172           0 :     if(name=="GridFT")
    1173           0 :       return new GridFT(ftrec);
    1174             :     //if(name=="rGridFT")
    1175             :     //  return new rGridFT(ftrec);
    1176           0 :     if(name=="WProjectFT")
    1177           0 :       return new WProjectFT(ftrec);
    1178             :     //if(name=="MultiTermFT")
    1179             :     //  return new MultiTermFT(ftrec);
    1180           0 :     if(name=="MosaicFT")
    1181           0 :       return new MosaicFT(ftrec);
    1182           0 :     if(name=="MosaicFTNew")
    1183           0 :       return new MosaicFTNew(ftrec);
    1184           0 :     if(name=="SetJyGridFT")
    1185           0 :       return new SetJyGridFT(ftrec);
    1186           0 :     if(name=="MultiTermFTNew")
    1187           0 :       return new MultiTermFTNew(ftrec);
    1188             :     //When the following have constructors from Record they should be uncommented
    1189             :     //   if(name=="AWProjectFT")
    1190             :     //  return new AWProjectFT(ftrec);
    1191             :     //if(name=="AWProjectWBFT")
    1192             :     //  return new  AWProjectWBFT(ftrec);
    1193             :     //if(name=="MultiTermAWProjectWBFT")
    1194             :     //  return new MultiTermAWProjectWBFT(ftrec);
    1195           0 :     return NULL;
    1196           0 :   }
    1197             : 
    1198           0 :   Int VisModelData::hasModel(Int msid, Int field, Int spw){
    1199             : 
    1200           0 :     IPosition oldcubeShape=ftindex_p.shape();
    1201           0 :     if(oldcubeShape(0) <(spw+1) || oldcubeShape(1) < (field+1) || oldcubeShape(2) < (msid+1)){
    1202           0 :       Cube<Int> newind(max((spw+1), oldcubeShape(0)), max((field+1),oldcubeShape(1)) , max((msid+1), oldcubeShape(2)));
    1203           0 :       newind.set(-1);
    1204           0 :       newind(IPosition(3, 0,0,0), (oldcubeShape-1))=ftindex_p;
    1205           0 :       ftindex_p.assign(newind);
    1206           0 :       newind.set(-1);
    1207           0 :       newind(IPosition(3, 0,0,0), (oldcubeShape-1))=clindex_p;
    1208           0 :       clindex_p.assign(newind);
    1209           0 :     }
    1210             : 
    1211           0 :     if( (clindex_p(spw, field, msid) + ftindex_p(spw, field, msid)) < -2)
    1212           0 :       return -2;
    1213           0 :     else if( (clindex_p(spw, field, msid) ==-1)  &&  (ftindex_p(spw, field, msid) ==-1))
    1214           0 :       return -1;
    1215           0 :     return 1;
    1216             : 
    1217             : 
    1218           0 :   }
    1219             : 
    1220           0 :    void VisModelData::initializeToVis(){
    1221             :     
    1222             : 
    1223           0 :   }
    1224             : 
    1225             : 
    1226             : 
    1227           0 :   void VisModelData::getMatchingMachines(Vector<CountedPtr<FTMachine> >& ft, Vector<CountedPtr<ComponentList> >& cl, const vi::VisBuffer2& vb){
    1228           0 :     if(isVersion2()){
    1229           0 :       Vector<Vector<Int> > combindx;
    1230           0 :       getUniqueIndicesComb(vb, combindx);
    1231           0 :       if(combindx.nelements() != 1)
    1232           0 :         throw(AipsError("Cannot deal with multiple intent per visbuffer "));
    1233           0 :       std::vector<Int> indexInBuf=combindx[0].tovector();
    1234           0 :       indexInBuf.push_back(vb.msId());
    1235           0 :       auto itFTMap = ftindex2_p.find(indexInBuf);
    1236           0 :       if(itFTMap != ftindex2_p.end() ) {
    1237           0 :         if((itFTMap->second) < 0){
    1238           0 :           ft.resize(0);
    1239             :         }
    1240             :         else{
    1241           0 :           ft=ftholder_p[itFTMap->second];
    1242             :                   }
    1243             :       }
    1244           0 :       auto itCLMap = clindex2_p.find(indexInBuf);
    1245           0 :       if(itCLMap != clindex2_p.end() ) {
    1246           0 :         if((itCLMap->second)< 0){
    1247           0 :           cl.resize(0);
    1248             :         }
    1249             :         else{
    1250             :           //cerr << "matching " << Vector<Int>(itCLMap->first) << "   " <<  itCLMap->second << endl;
    1251           0 :                         cl=clholder_p[itCLMap->second];
    1252             :                         //              cerr << "returned length  " << cl.nelements() << endl;
    1253             :         }
    1254             :       }
    1255             :       ///Now let's deal with a key that has not been visited before
    1256           0 :       if((itCLMap == clindex2_p.end()) && (itFTMap == ftindex2_p.end() )){
    1257             :         //cerr << "no matching holder " <<  Vector<Int>(indexInBuf) << " num of cl " << clholder_p.nelements() << endl;
    1258           0 :                   updateHolders(vb, indexInBuf);
    1259           0 :                   getMatchingMachines(ft, cl, vb);
    1260             :                   
    1261             :       }
    1262           0 :     }
    1263             :     else{
    1264           0 :       cl=getCL(vb.msId(), vb.fieldId()(0), vb.spectralWindows()(0));
    1265           0 :       ft=getFT(vb.msId(), vb.fieldId()(0), vb.spectralWindows()(0));
    1266             :     }
    1267           0 :     return;
    1268             :           
    1269             :   }
    1270           0 :   void VisModelData::updateHolders(const vi::VisBuffer2& vb, const std::vector<Int>& indexInBuf){
    1271           0 :           Int fieldId=vb.fieldId()[0];
    1272           0 :           const MeasurementSet& thems= vb.ms();
    1273           0 :           Int snum=-1;
    1274           0 :           Bool hasmodkey=False;
    1275           0 :           String modelkey;
    1276             :           //cerr <<"IndexInBuf " << Vector<Int>(indexInBuf) << " nelem " << ftindex2_p.size() << "   " << clindex2_p.size() << endl; 
    1277           0 :           hasmodkey=isModelDefined(fieldId, thems, modelkey, snum);
    1278           0 :           if(!hasmodkey){
    1279             :             //cerr << "NOHAS KEY fieldId " << fieldId << endl;
    1280           0 :                   clindex2_p[indexInBuf]=-1;
    1281           0 :                   ftindex2_p[indexInBuf]=-1;
    1282           0 :                   return;
    1283             :           }
    1284             :           //if we have already filled for this field
    1285           0 :          for (auto it=ftindex2_p.begin(); it != ftindex2_p.end(); ++it){
    1286             :             //  cerr << Vector<Int>(it->first) << "   val " << it->second << endl;
    1287           0 :            if((it->first)[0]==fieldId){
    1288           0 :              clindex2_p[indexInBuf]=-2;
    1289           0 :              ftindex2_p[indexInBuf]=-2;
    1290           0 :              return;
    1291             :            }
    1292             :              
    1293             :          } 
    1294             :           //We do have this key
    1295           0 :           TableRecord therec;
    1296           0 :           getModelRecord(modelkey, therec, thems);
    1297           0 :           if(!therec.isDefined("version") || (therec.asInt("version")!=2)){
    1298             :                 ///Set something versioning so as it can go to do version 1 way 
    1299           0 :                 version_p=1;
    1300           0 :                 clindex2_p[indexInBuf]=-1;
    1301           0 :                 ftindex2_p[indexInBuf]=-1;
    1302           0 :                 return;
    1303             :           }
    1304             :           else{
    1305           0 :             version_p=2;
    1306           0 :             Int indexft=-1;
    1307           0 :             if(therec.isDefined("numft")){
    1308           0 :               Int numft=therec.asInt("numft");
    1309           0 :               if(numft >0){
    1310           0 :                 for(Int ftk=0; ftk < numft; ++ftk){
    1311           0 :                   Record ftrec(therec.asRecord("ft_"+String::toString(ftk)));
    1312           0 :                   modifyDiskImagePath(ftrec, vb);
    1313           0 :                   const Matrix<Int>& ftcombind=ftrec.asArrayInt("indexcombination");
    1314             :                   //cerr << " ftcombind " << ftcombind << endl;
    1315           0 :                   indexft=ftholder_p.nelements();
    1316           0 :                   ftholder_p.resize(indexft+1, false, true);
    1317           0 :                   ftholder_p[indexft].resize(1);
    1318           0 :                   ftholder_p[indexft][0]=NEW_FT(ftrec.asRecord("container"));
    1319           0 :                   if(!( ftholder_p[indexft][0]))
    1320           0 :                     throw(AipsError("Unsupported FTMachine found in virtual MODEL_DATA column")); 
    1321           0 :                   ftholder_p[indexft][0]->initMaps(vb);
    1322           0 :                   for(uInt row=0; row < ftcombind.nrow(); ++row){
    1323           0 :                     std::vector<int> key=ftcombind.row(row).tovector();
    1324           0 :                     key.push_back(int(vb.msId()));
    1325             :                     
    1326           0 :                     if(ftindex2_p.count(key) >0){
    1327           0 :                       Int numftforkey=ftholder_p[ftindex2_p[key]].nelements();
    1328           0 :                       Int indx=ftindex2_p[key];
    1329           0 :                       Bool alreadyAdded=false;
    1330           0 :                       for (Int kk=1; kk < numftforkey; ++kk){
    1331           0 :                         alreadyAdded= alreadyAdded || (ftholder_p[indexft][0]==ftholder_p[indx][kk]);
    1332             :                       }
    1333           0 :                       if(!alreadyAdded){
    1334           0 :                         ftholder_p[indx].resize(numftforkey+1, true);
    1335           0 :                         ftholder_p[indx][numftforkey]=ftholder_p[indexft][0];
    1336             :                       }
    1337             :                     }
    1338             :                     else{
    1339           0 :                       ftindex2_p[key]=indexft;
    1340             :                     }
    1341           0 :                   }
    1342             :                   
    1343           0 :                 }
    1344             :               }
    1345             :             }
    1346           0 :             if(ftindex2_p.count(indexInBuf)==0)
    1347           0 :               ftindex2_p[indexInBuf]=-2;
    1348           0 :             Int indexcl=-1;
    1349             :             ////////////
    1350             :             //cerr <<"map " << endl;
    1351             :             //for (auto it=clindex2_p.begin(); it != clindex2_p.end(); ++it){
    1352             :             // cerr << Vector<Int>(it->first) << "   val " << it->second << endl;
    1353             : 
    1354             :             //}
    1355             :             ///////////
    1356           0 :             if(therec.isDefined("numcl")){
    1357           0 :                   Int numcl=therec.asInt("numcl");
    1358             :                   //cerr << "numcl " << numcl << endl;
    1359           0 :                   if(numcl >0){
    1360           0 :                     for(Int clk=0; clk < numcl; ++clk){
    1361             :           
    1362           0 :                       Record clrec(therec.asRecord("cl_"+String::toString(clk)));
    1363           0 :                       const Matrix<Int>& clcombind=clrec.asArrayInt("indexcombination");
    1364             :                       
    1365             :                       
    1366           0 :                       indexcl=clholder_p.nelements();
    1367           0 :                       clholder_p.resize(indexcl+1, false, true);
    1368           0 :                       clholder_p[indexcl].resize(1);
    1369           0 :                       clholder_p[indexcl][0]=new ComponentList();
    1370           0 :                       String err;
    1371           0 :                       if(!((clholder_p[indexcl][0])->fromRecord(err, clrec.asRecord("container"))))
    1372           0 :                         throw(AipsError("Component model failed to load for field "+String::toString(clcombind.column(0))));
    1373           0 :                       for(uInt row=0; row < clcombind.nrow(); ++row){
    1374           0 :                         std::vector<int> key=clcombind.row(row).tovector();
    1375             :                         
    1376           0 :                         key.push_back(int(vb.msId()));
    1377             :                         //cerr << "key " <<  Vector<Int>(key) << endl;
    1378           0 :                         if(clindex2_p.count(key) >0){
    1379           0 :                           Int numclforkey=clholder_p[clindex2_p[key]].nelements();
    1380           0 :                           Int indx=clindex2_p[key];
    1381           0 :                           Bool alreadyAdded=false;
    1382           0 :                           for (Int kk=1; kk < numclforkey; ++kk){
    1383           0 :                             alreadyAdded= alreadyAdded || (clholder_p[indexcl][0]==clholder_p[indx][kk]);
    1384             :                           }
    1385           0 :                           if(!alreadyAdded){
    1386           0 :                             clholder_p[indx].resize(numclforkey+1, true);
    1387           0 :                             clholder_p[indx][numclforkey]=clholder_p[indexcl][0];
    1388             :                           }
    1389             :                         }
    1390             :                         else{
    1391           0 :                           clindex2_p[key]=indexcl;
    1392             :                         }
    1393           0 :                       }
    1394             :                 
    1395           0 :                     }
    1396             :                   }
    1397             :                 }
    1398           0 :                 if(clindex2_p.count(indexInBuf)==0)
    1399           0 :                   clindex2_p[indexInBuf]=-2;
    1400             :           }
    1401             :           //cerr <<"update holder " << clholder_p.nelements() << endl; 
    1402             :          
    1403             :           
    1404           0 :   }
    1405           0 :   void VisModelData::getUniqueIndicesComb(const vi::VisBuffer2& vb, Vector< Vector<Int> >& retval){
    1406           0 :            const Vector<Int>& state = vb.stateId();
    1407           0 :            const Vector<Int>& scan=vb.scan();
    1408           0 :            const Vector<Double>& t=vb.time();
    1409           0 :            const Vector<Int>& fldid=vb.fieldId();
    1410           0 :            const Vector<Int>& spwid=vb.spectralWindows();
    1411           0 :            Vector<uInt>  uniqIndx;
    1412           0 :            uInt nTimes=GenSortIndirect<Double>::sort (uniqIndx, t, Sort::Ascending, Sort::QuickSort|Sort::NoDuplicates);
    1413           0 :            Vector<Int> comb(4);
    1414             :            
    1415           0 :            for (uInt k=0; k < nTimes; ++k){
    1416           0 :              comb(0)=fldid(uniqIndx[k]);
    1417           0 :              comb(1)=spwid(uniqIndx[k]);
    1418           0 :              comb(2)=scan(uniqIndx[k]);
    1419           0 :              comb(3)=state(uniqIndx[k]);
    1420           0 :              Bool hasComb=False;
    1421           0 :              if(retval.nelements() >0){
    1422           0 :                for (uInt j=0; j < retval.nelements(); ++j){
    1423           0 :                  if(allEQ(retval[j], comb))
    1424           0 :                    hasComb=True;
    1425             :                }
    1426             :                                    
    1427             :              }
    1428           0 :              if(!hasComb){
    1429           0 :                retval.resize(retval.nelements()+1, True);
    1430           0 :                retval[retval.nelements()-1]=comb;
    1431             :              }
    1432             :            }
    1433             :            
    1434           0 :   }
    1435             : 
    1436           0 :   void VisModelData::init(const VisBuffer2& vb){
    1437           0 :     if(version_p < 1){
    1438           0 :       updateHolders(vb, std::vector<Int>({vb.fieldId()(0), vb.spectralWindows()(0), vb.scan()(0), vb.stateId()(0), vb.msId()}));
    1439             : 
    1440             :     }
    1441             : 
    1442           0 :   }
    1443             :   
    1444             :   
    1445           0 :   Bool VisModelData::getModelVis(VisBuffer2& vb){
    1446             : 
    1447             :     //Vector<CountedPtr<ComponentList> >cl=getCL(vb.msId(), vb.fieldId()(0), vb.spectralWindows()(0));
    1448             :     //Vector<CountedPtr<FTMachine> > ft=getFT(vb.msId(), vb.fieldId()(0), vb.spectralWindows()(0));
    1449           0 :     Vector<CountedPtr<ComponentList> >cl(0);
    1450           0 :     Vector<CountedPtr<FTMachine> > ft(0);
    1451           0 :     getMatchingMachines(ft, cl, vb);
    1452             :     //Fill the buffer with 0.0; also prevents reading from disk if MODEL_DATA exists
    1453             :     ///Oh boy this is really dangerous...
    1454             :     //nCorr etc are public..who know who changed these values before reaching here.
    1455           0 :     Cube<Complex> mod(vb.nCorrelations(), vb.nChannels(), vb.nRows(), Complex(0.0));
    1456           0 :     vb.setVisCubeModel(mod);
    1457           0 :     Bool incremental=false;
    1458             :     //cerr << "CL nelements visbuff " << cl.nelements() << endl;
    1459           0 :     if( cl.nelements()>0){
    1460             :       //     cerr << "In cft " << cl.nelements() << endl;
    1461           0 :       for (uInt k=0; k < cl.nelements(); ++k)
    1462           0 :         if(!cl[k].null()){
    1463           0 :           cft_p->get(vb, *(cl[k]), -1); 
    1464             :           //cerr << "max " << max(vb.visCubeModel()) << endl;
    1465             :           //cout << "max " << max(vb.visCubeModel()) << endl;
    1466           0 :           incremental=true;
    1467             :         }
    1468             :     }
    1469           0 :     if(ft.nelements()>0){
    1470           0 :       Cube<Complex> tmpModel;
    1471           0 :       if(incremental || ft.nelements() >1)
    1472           0 :         tmpModel.assign(vb.visCubeModel());
    1473           0 :       Bool allnull=true;
    1474           0 :       for (uInt k=0; k < ft.nelements(); ++k){
    1475           0 :         if(!ft[k].null()){
    1476           0 :           if(k >0) vb.setVisCubeModel(Cube<Complex> (vb.nCorrelations(), vb.nChannels(), vb.nRows(), Complex(0.0)));
    1477             :     
    1478           0 :           ft[k]->get(vb, -1);
    1479           0 :           if(ft.nelements()>1 || incremental){
    1480           0 :             tmpModel+=vb.visCubeModel();
    1481             :           }
    1482           0 :           allnull=false;
    1483             :         }
    1484             :       }
    1485             :       //cerr << "min max after ft " << min(vb.visCubeModel()) << max(vb.visCubeModel()) << endl;
    1486           0 :       if(!allnull){
    1487           0 :         if(ft.nelements()>1 || incremental)
    1488           0 :           vb.setVisCubeModel(tmpModel);
    1489           0 :         incremental=true;
    1490             :       }      
    1491           0 :     }
    1492           0 :     if(!incremental){
    1493             :       //No model was set so....
    1494             :       ///Set the Model to 1.0 for parallel hand and 0.0 for x-hand
    1495             :       
    1496           0 :       vb.setVisCubeModel(Complex(1.0));
    1497           0 :       Cube<Complex> modelCube(vb.visCubeModel());
    1498           0 :       Vector<Stokes::StokesTypes> corrType=vb.getCorrelationTypesSelected();
    1499           0 :       uInt nCorr = corrType.nelements();
    1500           0 :       for (uInt i=0; i<nCorr; i++) {
    1501           0 :           if(corrType[i]==Stokes::RL || corrType[i]==Stokes::LR ||
    1502           0 :              corrType[i]==Stokes::XY || corrType[i]==Stokes::YX){
    1503           0 :                   modelCube.yzPlane(i)=0.0;
    1504             :           }
    1505             :       }
    1506           0 :     }
    1507             :     
    1508           0 :     return true;
    1509             :     
    1510           0 :   }
    1511             : 
    1512             : 
    1513           0 :   Vector<CountedPtr<ComponentList> > VisModelData::getCL(const Int msId, const Int fieldId, const Int spwId){
    1514           0 :     if(!hasModel(msId, fieldId, spwId))
    1515           0 :       return Vector<CountedPtr<ComponentList> >(0);
    1516           0 :     Int indx=clindex_p(spwId, fieldId, msId);
    1517             :     //  cerr << "indx " << indx << "   " << clholder_p[indx].nelements() <<  " spw " << spwId << " field " << fieldId << endl;
    1518           0 :     if(indx <0)
    1519           0 :       return Vector<CountedPtr<ComponentList> >(0);
    1520             :     // cerr << "CL " << clholder_p[indx][0]->summarize(0) << endl;
    1521           0 :     return clholder_p[indx];
    1522             :         
    1523             : 
    1524             :   }
    1525             : 
    1526           0 :   Vector<CountedPtr<FTMachine> >VisModelData::getFT(const Int msId, const Int fieldId, Int spwId){
    1527             : 
    1528           0 :     if(!hasModel(msId, fieldId, spwId))
    1529           0 :       return Vector<CountedPtr<FTMachine> >(0);
    1530           0 :     Int indx=ftindex_p(spwId, fieldId, msId);
    1531             :     //cerr << "indx " << indx << endl;
    1532           0 :     if(indx <0)
    1533           0 :       return Vector<CountedPtr<FTMachine> >(0);
    1534           0 :     return ftholder_p[indx];
    1535             :   }
    1536             : 
    1537           0 :   Int VisModelData::firstSourceRowRecord(const Int field, const MeasurementSet& theMS, TableRecord& rec){
    1538           0 :     Int row=-1;
    1539             :     
    1540             :     //Prefer the Source table first    
    1541           0 :     MSFieldColumns fCol(theMS.field());
    1542           0 :     if(Table::isReadable(theMS.sourceTableName()) && (theMS.source().nrow() > 0) &&  (!fCol.sourceId().isNull()) && (fCol.sourceId().get(field) != -1) ){
    1543             :     
    1544           0 :       const MSSource& mss=theMS.source();
    1545             :       
    1546           0 :       Int sid=fCol.sourceId().get(field);
    1547           0 :       auto rows=MSSourceIndex(mss).getRowNumbersOfSourceID(sid);
    1548           0 :       if(rows.nelements() > 0) 
    1549           0 :         row=rows[0];
    1550           0 :       const TableRecord& keywords=mss.keywordSet();
    1551           0 :       Bool rowIsUsed=false;
    1552           0 :       for (uInt n=0; n< keywords.nfields(); ++n){
    1553           0 :         if(keywords.dataType(n) == TpInt){
    1554           0 :           if(row==keywords.asInt(n)){
    1555           0 :             rowIsUsed=true;
    1556             :           }
    1557             :         }
    1558             :       }
    1559             :       //nobody is using that row ..any record there must be dangling
    1560           0 :       if(!rowIsUsed)
    1561           0 :         return -1;
    1562           0 :       if(row >-1 && mss.isColumn(MSSource::SOURCE_MODEL)){
    1563           0 :         ScalarColumn<TableRecord> scol(theMS.source(), "SOURCE_MODEL");
    1564           0 :         scol.get(row, rec);
    1565           0 :       }
    1566           0 :     }
    1567           0 :     return row;
    1568           0 :   }
    1569             : 
    1570             : 
    1571             : }// end namespace refim
    1572             : }//# NAMESPACE CASA - END
    1573             : 

Generated by: LCOV version 1.16