C++ ActionSetup::TopAddress方法代码示例

// Action_Unwrap::Setup()
Action::RetType Action_Unwrap::Setup(ActionSetup& setup) {
  // Ensure same number of atoms in current parm and ref parm
  if ( RefParm_!=0 ) {
    if ( setup.Top().Natom() != RefParm_->Natom() ) {
      mprinterr("Error: unwrap: # atoms in reference parm %s is not\n", RefParm_->c_str());
      mprinterr("Error:         equal to # atoms in parm %s\n", setup.Top().c_str());
      return Action::ERR;
    }
  }
  // Check box type
  if (setup.CoordInfo().TrajBox().Type()==Box::NOBOX) {
    mprintf("Error: unwrap: Parm %s does not contain box information.\n",
            setup.Top().c_str());
    return Action::ERR;
  }
  orthogonal_ = false;
  if (setup.CoordInfo().TrajBox().Type()==Box::ORTHO)
    orthogonal_ = true;

  // Setup atom pairs to be unwrapped.
  imageList_ = Image::CreatePairList(setup.Top(), imageMode_, maskExpression_);
  if (imageList_.empty()) {
    mprintf("Warning: Mask selects no atoms for topology '%s'.\n", setup.Top().c_str());
    return Action::SKIP;
  }
  mprintf("\tNumber of %ss to be unwrapped is %zu\n",
          Image::ModeString(imageMode_), imageList_.size()/2);

  // Use current parm as reference if not already set
  if (RefParm_ == 0)
    RefParm_ = setup.TopAddress();
  return Action::OK;
}

/** Determine what atoms each mask pertains to for the current parm file.
  */
Action::RetType Action_LIE::Setup(ActionSetup& setup) {
  if (setup.Top().SetupIntegerMask( Mask1_ )) return Action::ERR;
  if (setup.Top().SetupIntegerMask( Mask2_ )) return Action::ERR;

  mprintf("\tLIE: %i Ligand Atoms, %i Surrounding Atoms\n",
          Mask1_.Nselected(), Mask2_.Nselected());

  if (setup.CoordInfo().TrajBox().Type() == Box::NOBOX) {
    mprinterr("Error: LIE: Must have explicit solvent system with box info\n");
    return Action::ERR;
  }

  if (Mask1_.None() || Mask2_.None()) {
    mprintf("Warning: LIE: One or both masks have no atoms.\n");
    return Action::SKIP;
  }

  if (SetupParms(setup.Top()))
    return Action::ERR;

  // Back up the parm
  CurrentParm_ = setup.TopAddress();

  return Action::OK;
}

/** Set angle up for this parmtop. Get masks etc.
  */
Action::RetType Action_Esander::Setup(ActionSetup& setup) {
  if (currentParm_ != 0 && currentParm_->Pindex() != setup.Top().Pindex())
  {
    mprintf("Warning: Current topology is %i:%s but reference is %i:%s. Skipping.\n",
            setup.Top().Pindex(), setup.Top().c_str(),
            currentParm_->Pindex(), currentParm_->c_str());
    return Action::SKIP;
  }
  // Check for LJ terms
  if (!setup.Top().Nonbond().HasNonbond())
  {
    mprinterr("Error: Topology '%s' does not have non-bonded parameters.\n", setup.Top().c_str());
    return Action::ERR;
  }
  // If reference specified, init now. Otherwise using first frame.
  if (currentParm_ != 0 ) {
    if ( InitForRef() ) return Action::ERR;
  } else
    currentParm_ = setup.TopAddress();
  // If saving of forces is requested, make sure CoordinateInfo has force.
  if (save_forces_) {
    cInfo_ = setup.CoordInfo();
    cInfo_.SetForce( true );
    newFrame_.SetupFrameV( setup.Top().Atoms(), cInfo_ );
    setup.SetCoordInfo( &cInfo_ );
    ret_ = Action::MODIFY_COORDS;
    return Action::MODIFY_TOPOLOGY;
  }
  ret_ = Action::OK;
  return Action::OK;
}

/** Set up solute and solvent masks. If no solvent mask was specified use 
  * solvent information in the current topology.
  */
Action::RetType Action_Watershell::Setup(ActionSetup& setup) {
  // Set up solute mask
  if (setup.Top().SetupIntegerMask( soluteMask_ )) return Action::ERR;
  if ( soluteMask_.None() ) {
    mprintf("Warning: No atoms in solute mask [%s].\n",soluteMask_.MaskString());
    return Action::SKIP;
  }
  // Set up solvent mask
  if (!solventmaskexpr_.empty()) {
    if (setup.Top().SetupIntegerMask( solventMask_ )) return Action::ERR;
  } else {
    solventMask_.ResetMask();
    for (Topology::mol_iterator mol = setup.Top().MolStart();
                                mol != setup.Top().MolEnd(); ++mol)
    {
      if ( mol->IsSolvent() )
        solventMask_.AddAtomRange( mol->BeginAtom(), mol->EndAtom() );
    }
  }
  if ( solventMask_.None() ) {
    if (!solventmaskexpr_.empty())
      mprintf("Warning: No solvent atoms selected by mask [%s]\n",
                solventmaskexpr_.c_str());
    else
      mprintf("Warning: No solvent atoms in topology %s\n",setup.Top().c_str());
    return Action::SKIP;
  }
  SetupImaging( setup.CoordInfo().TrajBox().Type() );
  // Create space for residues
# ifdef _OPENMP
  // Only re-allocate for larger # of residues
  if ( setup.Top().Nres() > NactiveResidues_ ) {
    if (activeResidues_thread_ != 0) {
      // Deallocate each thread
      for (int i = 0; i < NactiveResidues_; ++i)
        delete[] activeResidues_thread_[i];
    } else {
      // Initial thread allocation needed
      activeResidues_thread_ = new int*[ numthreads_ ];
    }
    // Allocate each thread
    for (int i = 0; i < numthreads_; ++i) {
      activeResidues_thread_[i] = new int[ setup.Top().Nres() ];
      std::fill( activeResidues_thread_[i], activeResidues_thread_[i] + setup.Top().Nres(), 0 );
    }
  }
  NactiveResidues_ = setup.Top().Nres();
# else
  activeResidues_.resize( setup.Top().Nres(), 0 );
# endif
  // Store current Parm
  CurrentParm_ = setup.TopAddress();
  return Action::OK;    
}

// Action_Outtraj::Setup()
Action::RetType Action_Outtraj::Setup(ActionSetup& setup) {
  if (associatedParm_->Pindex() != setup.Top().Pindex())
    return Action::SKIP;
  if (!isSetup_) { // TODO: Trajout IsOpen?
    if (outtraj_.SetupTrajWrite(setup.TopAddress(), setup.CoordInfo(), setup.Nframes()))
      return Action::ERR;
    outtraj_.PrintInfo(0);
    isSetup_ = true;
  }
  return Action::OK;
}

// Action_CheckStructure::Setup()
Action::RetType Action_CheckStructure::Setup(ActionSetup& setup) {
  CurrentParm_ = setup.TopAddress();
  if (SeparateSetup( setup.Top(), setup.CoordInfo().TrajBox().Type(),bondcheck_ ))
    return Action::ERR;
  // Print imaging info for this parm
  if (bondcheck_)
    mprintf("\tChecking %u bonds.\n", bondList_.size());
  if (image_.ImagingEnabled())
    mprintf("\tImaging on.\n");
  else
    mprintf("\timaging off.\n");
  return Action::OK;
}

// Action_Outtraj::Setup()
Action::RetType Action_Outtraj::Setup(ActionSetup& setup) {
  if (!isActive_ || associatedParm_->Pindex() != setup.Top().Pindex()) {
    mprintf("\tOutput trajectory not active for topology '%s'\n", setup.Top().c_str());
    return Action::SKIP;
  }
  if (!isSetup_) { // TODO: Trajout IsOpen?
    if (outtraj_.SetupTrajWrite(setup.TopAddress(), setup.CoordInfo(), setup.Nframes()))
      return Action::ERR;
    mprintf("      "); //TODO this is a kludge; PrintInfo should be a string.
    outtraj_.PrintInfo(0);
    mprintf("\tHas %s\n", outtraj_.Traj().CoordInfo().InfoString().c_str());
    isSetup_ = true;
  }
  return Action::OK;
}

// Action_NativeContacts::Setup()
Action::RetType Action_NativeContacts::Setup(ActionSetup& setup) {
  // Setup potential contact lists for this topology
  if (SetupContactLists( setup.Top(), Frame()))
    return Action::SKIP;
  mprintf("\t%zu potential contact sites for '%s'\n", Mask1_.Nselected(), Mask1_.MaskString());
  if (Mask2_.MaskStringSet())
    mprintf("\t%zu potential contact sites for '%s'\n", Mask2_.Nselected(), Mask2_.MaskString());
  // Set up imaging info for this parm
  image_.SetupImaging( setup.CoordInfo().TrajBox().Type() );
  if (image_.ImagingEnabled())
    mprintf("\tImaging enabled.\n");
  else
    mprintf("\tImaging disabled.\n");
  CurrentParm_ = setup.TopAddress();
  return Action::OK;
}

/** Set up mask, allocate memory for exclusion list.
  */
Action::RetType Action_Pairwise::Setup(ActionSetup& setup) {
  // Set up mask
  if ( setup.Top().SetupIntegerMask( Mask0_ ) ) return Action::ERR;
  if (Mask0_.None()) {
    mprintf("Warning: Mask has no atoms.\n");
    return Action::SKIP;
  }

  // Set up exclusion list and determine total # interactions.
  int N_interactions = SetupNonbondParm(Mask0_, setup.Top());
  if (N_interactions == -1) return Action::ERR;
  // Allocate matrix memory
  int previous_size = (int)vdwMat_->Size();
  if (previous_size == 0) {
    vdwMat_->AllocateTriangle( setup.Top().Natom() );
    eleMat_->AllocateTriangle( setup.Top().Natom() );
  } else {
    if (previous_size != N_interactions) {
      mprinterr("Error: Attempting to reallocate matrix with different size.\n"
                "Error:   Original size= %i, new size= %i\n"
                "Error:   This can occur when different #s of atoms are selected in\n"
                "Error:   different topology files.\n", previous_size, N_interactions);
      return Action::ERR;
    }
  }

  // If comparing to a reference frame for atom-by-atom comparison make sure
  // the number of interactions is the same in reference and parm.
  if (nb_calcType_==COMPARE_REF) {
    if (N_interactions != N_ref_interactions_) {
      mprinterr("Error: # reference interactions (%i) != # interactions for this parm (%i)\n",
                N_ref_interactions_, N_interactions);
      return Action::ERR;
    }
  }
  // Set up cumulative energy arrays
  atom_eelec_.clear();
  atom_eelec_.resize(setup.Top().Natom(), 0.0);
  atom_evdw_.clear();
  atom_evdw_.resize(setup.Top().Natom(), 0.0);
  // Print pairwise info for this parm
  Mask0_.MaskInfo();
  CurrentParm_ = setup.TopAddress();
  return Action::OK;  
}

/** Set angle up for this parmtop. Get masks etc.
  */
Action::RetType Action_Energy::Setup(ActionSetup& setup) {
  if (setup.Top().SetupCharMask(Mask1_)) return Action::ERR;
  if (Mask1_.None()) {
    mprintf("Warning: Mask '%s' selects no atoms.\n", Mask1_.MaskString());
    return Action::SKIP;
  }
  Mask1_.MaskInfo();
  Imask_ = AtomMask(Mask1_.ConvertToIntMask(), Mask1_.Natom());
  // Check for LJ terms
  for (calc_it calc = Ecalcs_.begin(); calc != Ecalcs_.end(); ++calc)
    if ((*calc == N14 || *calc == NBD) && !setup.Top().Nonbond().HasNonbond())
    {
      mprinterr("Error: Nonbonded energy calc requested but topology '%s'\n"
                "Error:   does not have non-bonded parameters.\n", setup.Top().c_str());
      return Action::ERR;
    }
  currentParm_ = setup.TopAddress();
  return Action::OK;
}

// Action_Vector::Setup()
Action::RetType Action_Vector::Setup(ActionSetup& setup) {
  if (needBoxInfo_) {
    // Check for box info
    if (setup.CoordInfo().TrajBox().Type() == Box::NOBOX) {
      mprinterr("Error: vector box: No box information.\n",
                setup.Top().c_str());
      return Action::ERR;
    }
  }
  if (mask_.MaskStringSet()) {
    // Setup mask 1
    if (setup.Top().SetupIntegerMask(mask_)) return Action::ERR;
    mask_.MaskInfo();
    if (mask_.None()) {
      mprinterr("Error: First vector mask is empty.\n");
      return Action::ERR;
    }
  }

  // Allocate space for CORRPLANE. 
  if (mode_ == CORRPLANE) {
    if (vcorr_!=0) delete[] vcorr_;
    vcorr_ = new double[ 3 * mask_.Nselected() ];
  }

  // Setup mask 2
  if (mask2_.MaskStringSet()) {
    if (setup.Top().SetupIntegerMask(mask2_)) return Action::ERR;
    mask2_.MaskInfo();
    if (mask2_.None()) {
      mprinterr("Error: Second vector mask is empty.\n");
      return Action::ERR;
    }
  }
  CurrentParm_ = setup.TopAddress();
  return Action::OK;
}

/** Called every time the trajectory changes. Set up FrameMask for the new 
  * parmtop and allocate space for selected atoms from the Frame.
  */
Action::RetType Action_Rmsd::Setup(ActionSetup& setup) {
  // Target setup
  if ( setup.Top().SetupIntegerMask( tgtMask_ ) ) return Action::ERR;
  mprintf("\tTarget mask:");
  tgtMask_.BriefMaskInfo();
  mprintf("\n");
  if ( tgtMask_.None() ) {
    mprintf("Warning: No atoms in mask '%s'.\n", tgtMask_.MaskString());
    return Action::SKIP;
  }
  // Allocate space for selected atoms in the frame. This will also put the
  // correct masses in based on the mask.
  tgtFrame_.SetupFrameFromMask(tgtMask_, setup.Top().Atoms());
  // Reference setup
  if (REF_.SetupRef(setup.Top(), tgtMask_.Nselected(), "rmsd"))
    return Action::SKIP;
 
  // Per residue rmsd setup
  if (perres_) {
    // If RefParm is still NULL probably 'first', set now.
    if (RefParm_ == 0)
      RefParm_ = setup.TopAddress();
    int err = perResSetup(setup.Top(), *RefParm_);
    if      (err == 1) return Action::SKIP;
    else if (err == 2) return Action::ERR;
  }

  // Warn if PBC and rotating
  if (rotate_ && setup.CoordInfo().TrajBox().Type() != Box::NOBOX) {
    mprintf("Warning: Coordinates are being rotated and box coordinates are present.\n"
            "Warning: Unit cell vectors are NOT rotated; imaging will not be possible\n"
            "Warning:  after the RMS-fit is performed.\n");
  }

  return Action::OK;
}

Action::RetType Action_Contacts::Setup(ActionSetup& setup) {
  //if (first_) 
  //  RefParm_ = currentParm;
  // Set up atom mask 
  if (setup.Top().SetupIntegerMask(Mask_)) return Action::ERR;

  // Determine which residues are active based on the mask
  activeResidues_.clear();
  for (AtomMask::const_iterator atom = Mask_.begin();
                                atom != Mask_.end(); ++atom)
  {
    int resnum = setup.Top()[*atom].ResNum();
    activeResidues_.insert( resnum );
  }

  // byresidue header - only on first time through
  if (residueContacts_.empty() && byResidue_) {
    outfile_->Printf("#time");
    outfile2_->Printf("#time");
    for (std::set<int>::iterator res = activeResidues_.begin();
                                 res != activeResidues_.end(); ++res)
    {
      outfile_->Printf("\tresidue %i", *res);
      outfile2_->Printf("\tresidue %i", *res);
    }
    outfile_->Printf("\tTotal\n");
    outfile2_->Printf("\tTotal\n");
  }

  // Reserve space for residue contact counts
  residueContacts_.reserve( setup.Top().Nres() );
  residueNative_.reserve( setup.Top().Nres() );
  
  CurrentParm_ = setup.TopAddress();
  return Action::OK;
}

/** Determine what atoms each mask pertains to for the current parm file.
  */
Action::RetType Action_NMRrst::Setup(ActionSetup& setup) {
  if (!viewrst_.empty() && rsttop_ == 0) rsttop_ = setup.TopAddress();
  // ---------------------------------------------
  // Set up NOEs from file.
  for (noeDataArray::iterator noe = NOEs_.begin(); noe != NOEs_.end(); ++noe) {
    if (setup.Top().SetupIntegerMask( noe->dMask1_ )) return Action::ERR;
    if (setup.Top().SetupIntegerMask( noe->dMask2_ )) return Action::ERR;
    if (noe->dMask1_.None() || noe->dMask2_.None()) {
      mprintf("Warning: One or both masks for NOE '%s' have no atoms (%i and %i).\n",
              noe->dist_->legend(), noe->dMask1_.Nselected(),
              noe->dMask2_.Nselected());
      noe->active_ = false; 
    } else
      noe->active_ = true;
  }
  // ---------------------------------------------
  // Set up potential NOE sites.
  if (findNOEs_) {
    if (setup.Top().SetupCharMask( Mask_ )) return Action::ERR;
    Mask_.MaskInfo();
    if (Mask_.None()) return Action::SKIP;
    SiteArray potentialSites; // .clear();
    AtomMap resMap;
    resMap.SetDebug( debug_ );
    std::vector<bool> selected;
    Range soluteRes = setup.Top().SoluteResidues();
    for (Range::const_iterator res = soluteRes.begin(); res != soluteRes.end(); ++res)
    {
      int res_first_atom = setup.Top().Res(*res).FirstAtom();
      selected.assign( setup.Top().Res(*res).NumAtoms(), false );
      // Find symmetric atom groups.
      AtomMap::AtomIndexArray symmGroups;
      if (resMap.SymmetricAtoms(setup.Top(), symmGroups, *res)) return Action::ERR;
      // DEBUG
      if (debug_ > 0) {
        mprintf("DEBUG: Residue %i: symmetric atom groups:\n", *res + 1);
        for (AtomMap::AtomIndexArray::const_iterator grp = symmGroups.begin();
                                                     grp != symmGroups.end(); ++grp)
        {
          mprintf("\t\t");
          for (AtomMap::Iarray::const_iterator at = grp->begin();
                                               at != grp->end(); ++at)
            mprintf(" %s", setup.Top().TruncAtomNameNum( *at ).c_str());
          mprintf("\n");
        }
      }
      // Each symmetric hydrogen atom group is a site.
      for (AtomMap::AtomIndexArray::const_iterator grp = symmGroups.begin();
                                                   grp != symmGroups.end(); ++grp)
      { // NOTE: If first atom is H all should be H.
        if ( setup.Top()[ grp->front() ].Element() == Atom::HYDROGEN )
        {
          Iarray symmAtomGroup;
          for (Iarray::const_iterator at = grp->begin();
                                      at != grp->end(); ++at)
            if (Mask_.AtomInCharMask( *at ))
              symmAtomGroup.push_back( *at );
          if (!symmAtomGroup.empty()) {
            potentialSites.push_back( Site(*res, symmAtomGroup) );
            // Mark symmetric atoms as selected.
            for (AtomMap::Iarray::const_iterator at = grp->begin();
                                                 at != grp->end(); ++at)
              selected[ *at - res_first_atom ] = true;
          }
        }
      }
      // All other non-selected hydrogens bonded to same heavy atom are sites.
      for (int ratom = res_first_atom; ratom != setup.Top().Res(*res).LastAtom(); ++ratom)
      {
        if ( setup.Top()[ratom].Element() != Atom::HYDROGEN ) {
          Iarray heavyAtomGroup;
          for (Atom::bond_iterator ba = setup.Top()[ratom].bondbegin();
                                   ba != setup.Top()[ratom].bondend(); ++ba)
            if ( Mask_.AtomInCharMask(*ba) && 
                 *ba >= res_first_atom && 
                 *ba < setup.Top().Res(*res).LastAtom() )
            {
              if ( !selected[ *ba - res_first_atom ] &&
                   setup.Top()[ *ba ].Element() == Atom::HYDROGEN )
                heavyAtomGroup.push_back( *ba );
            }
          if (!heavyAtomGroup.empty())
            potentialSites.push_back( Site(*res, heavyAtomGroup) );
        }
      }
    }
    mprintf("\t%zu potential NOE sites:\n", potentialSites.size());
    for (SiteArray::const_iterator site = potentialSites.begin();
                                   site != potentialSites.end(); ++site)
    {
      mprintf("  %u\tRes %i:", site - potentialSites.begin(), site->ResNum()+1);
      for (unsigned int idx = 0; idx != site->Nindices(); ++idx)
        mprintf(" %s", setup.Top().TruncAtomNameNum( site->Idx(idx) ).c_str());
      mprintf("\n");
    }
    if (noeArray_.empty()) {
      size_t siteArraySize = 0;
      // Set up all potential NOE pairs. Keep track of size.
//.........这里部分代码省略.........

/** Determine from selected mask atoms which dihedrals will be rotated. */
Action::RetType Action_DihedralScan::Setup(ActionSetup& setup) {
  DihedralScanType dst;
  // If range is empty (i.e. no resrange arg given) look through all 
  // solute residues.
  Range actualRange;
  if (resRange_.Empty())
    actualRange = setup.Top().SoluteResidues();
  else 
    actualRange = resRange_;
  // Search for dihedrals
  if (dihSearch_.FindDihedrals(setup.Top(), actualRange))
    return Action::ERR;
  // For each found dihedral, set up mask of atoms that will move upon 
  // rotation. Also set up mask of atoms in this residue that will not
  // move, including atom2.
  if (debug_>0)
    mprintf("DEBUG: Dihedrals:\n");
  for (DihedralSearch::mask_it dih = dihSearch_.begin();
                               dih != dihSearch_.end(); ++dih)
  {
    dst.checkAtoms.clear();
    // Set mask of atoms that will move during dihedral rotation.
    dst.Rmask = DihedralSearch::MovingAtoms(setup.Top(), dih->A1(), dih->A2());
    // If randomly rotating angles, check for atoms that are in the same
    // residue as A1 but will not move. They need to be checked for clashes
    // since further rotations will not help them.
    if (mode_ == RANDOM && check_for_clashes_) {
      CharMask cMask( dst.Rmask.ConvertToCharMask(), dst.Rmask.Nselected() );
      int a1res = setup.Top()[dih->A1()].ResNum();
      for (int maskatom = setup.Top().Res(a1res).FirstAtom();
               maskatom < setup.Top().Res(a1res).LastAtom(); ++maskatom)
        if (!cMask.AtomInCharMask(maskatom))
          dst.checkAtoms.push_back( maskatom );
      dst.checkAtoms.push_back(dih->A1()); // TODO: Does this need to be added first?
      // Since only the second atom and atoms it is bonded to move during 
      // rotation, base the check on the residue of the second atom.
      dst.resnum = a1res;
    }
    dst.atom0 = dih->A0(); // FIXME: This duplicates info
    dst.atom1 = dih->A1();
    dst.atom2 = dih->A2();
    dst.atom3 = dih->A3();
    BB_dihedrals_.push_back(dst);
    // DEBUG: List dihedral info.
    if (debug_ > 0) {
      mprintf("\t%s-%s-%s-%s\n", 
              setup.Top().TruncResAtomName(dih->A0()).c_str(),
              setup.Top().TruncResAtomName(dih->A1()).c_str(),
              setup.Top().TruncResAtomName(dih->A2()).c_str(),
              setup.Top().TruncResAtomName(dih->A3()).c_str() );
      if (debug_ > 1 && mode_ == RANDOM && check_for_clashes_) {
        mprintf("\t\tCheckAtoms=");
        for (std::vector<int>::const_iterator ca = dst.checkAtoms.begin();
                                              ca != dst.checkAtoms.end(); ++ca)
          mprintf(" %i", *ca + 1);
        mprintf("\n");
      }
      if (debug_ > 2) {
        mprintf("\t\t");
        dst.Rmask.PrintMaskAtoms("Rmask:");
      }
    }
  }

  // Set up CheckStructure for this parm (false = nobondcheck)
  if (checkStructure_.SeparateSetup(setup.Top(),
                                    setup.CoordInfo().TrajBox().Type(), false) != Action::OK)
    return Action::ERR;

  // Set the overall max number of rotations to try
  max_rotations_ = (int) BB_dihedrals_.size();
  max_rotations_ *= max_factor_;

  // Set up simple structure check. First step is coarse; check distances 
  // between a certain atom in each residue (first, COM, CA, some other atom?)
  // to see if residues are in each others neighborhood. Second step is to 
  // check the atoms in each close residue.
  if (check_for_clashes_) {
    ResidueCheckType rct;
    int res = 0;
    for (Topology::res_iterator residue = setup.Top().ResStart();
                                residue != setup.Top().ResEnd(); ++residue)
    {
      rct.resnum = res++;
      rct.start = residue->FirstAtom();
      rct.stop = residue->LastAtom();
      rct.checkatom = rct.start;
      ResCheck_.push_back(rct);
    }
  }

  if (!outfilename_.empty() && CurrentParm_ == 0) // FIXME: Correct frames for # of rotations
    outtraj_.SetupTrajWrite(setup.TopAddress(), setup.CoordInfo(), setup.Nframes());

  CurrentParm_ = setup.TopAddress();
  return Action::OK;  
}