From 2bb50b0c47e1db9720663755e88c6a1d21d8e43f Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 12:25:23 -0400 Subject: [PATCH 01/20] DRR - Cpptraj: Allow states to have multiple criteria --- src/Analysis_State.cpp | 88 +++++++++++++++++++++++++++++------------- src/Analysis_State.h | 37 +++++++++++++----- 2 files changed, 90 insertions(+), 35 deletions(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index 0cf98d9cfc..c6c2f9a15d 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -27,10 +27,13 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set while (!state_arg.empty()) { // Expect form <#>, ArgList argtmp(state_arg, ","); - if (argtmp.Nargs() != 4) { + if (argtmp.Nargs() < 4) { mprinterr("Error: Malformed state argument '%s': expect ,,,\n", state_arg.c_str()); return Analysis::ERR; + } else if (((argtmp.Nargs()-1) % 3) != 0) { + mprinterr("Error: Each state arg must have a and value.\n"); + return Analysis::ERR; } std::string state_id = argtmp.GetStringNext(); // TODO: Check duplicate names @@ -38,19 +41,24 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set mprinterr("Error: In state argument, could not get ID.\n"); return Analysis::ERR; } - DataSet* ds = setup.DSL().GetDataSet( argtmp.GetStringNext() ); - if (ds == 0) return Analysis::ERR; - if (ds->Ndim() != 1) { - mprinterr("Error: Only 1D data sets allowed.\n"); - return Analysis::ERR; + StateType currentState( state_id ); + int nSets = (argtmp.Nargs()-1) / 3; + for (int setArg = 0; setArg < nSets; setArg++) { + DataSet* ds = setup.DSL().GetDataSet( argtmp.GetStringNext() ); + if (ds == 0) return Analysis::ERR; + if (ds->Ndim() != 1) { + mprinterr("Error: Only 1D data sets allowed.\n"); + return Analysis::ERR; + } + double min = argtmp.getNextDouble(0.0); + double max = argtmp.getNextDouble(0.0); + if (max < min) { + mprinterr("Error: max value cannot be less than min.\n"); + return Analysis::ERR; + } + currentState.AddCriterion( (DataSet_1D*)ds, min, max ); } - double min = argtmp.getNextDouble(0.0); - double max = argtmp.getNextDouble(0.0); - if (max < min) { - mprinterr("Error: max value cannot be less than min.\n"); - return Analysis::ERR; - } - States_.push_back( StateType(state_id, (DataSet_1D*)ds, min, max) ); + States_.push_back( currentState ); state_arg = analyzeArgs.GetStringKey("state"); } } @@ -64,8 +72,10 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set mprintf(" STATE: The following states have been set up:\n"); for (StateArray::const_iterator state = States_.begin(); state != States_.end(); ++state) - mprintf("\t%u: %20s %12.4f <= %-20s < %12.4f\n", state - States_.begin(), state->DS().legend(), - state->Min(), state->id(), state->Max()); + { + mprintf("\t%u: ", state - States_.begin()); + state->PrintState(); + } mprintf("\tState data set: %s\n", state_data_->legend()); if (outfile != 0) mprintf("\tStates vs time output to file '%s'\n", outfile->DataFilename().full()); @@ -90,14 +100,17 @@ std::string const& Analysis_State::StateName(int idx) const { Analysis::RetType Analysis_State::Analyze() { // Only process as much data as is in the smallest data set. - size_t nframes = States_.front().DS().Size(); + size_t nframes = States_.front().Nframes(); for (StateArray::const_iterator state = States_.begin(); state != States_.end(); ++state) { - if ( state != States_.begin() && nframes != state->DS().Size() ) - mprintf("Warning: Set '%s' for state '%s' has a different # of frames (%zu)" - " than previous sets (%zu).\n", state->DS().legend(), state->id(), - state->DS().Size(), nframes); - nframes = std::min( nframes, state->DS().Size() ); + if ( state != States_.begin() ) { + size_t nf = state->Nframes(); + if (nframes != nf) { + mprintf("Warning: State '%s' has a different # of frames (%zu) " + "than previous states (%zu).\n", state->id(), nf, nframes); + nframes = std::min( nframes, nf ); + } + } } mprintf("\tProcessing %zu frames.\n", nframes); if (nframes < 1) return Analysis::ERR; @@ -121,12 +134,13 @@ Analysis::RetType Analysis_State::Analyze() { int state_num = -1; for (StateArray::const_iterator state = States_.begin(); state != States_.end(); ++state) { - double dVal = state->DS().Dval( frm ); - if (dVal >= state->Min() && dVal < state->Max()) { // TODO: Periodic - if (state_num != -1) + if (state->InState( frm )) { + if (state_num != -1) { mprintf("Warning: Frame %zu already defined as state '%s', also matches state '%s'.\n", frm, States_[state_num].id(), state->id()); - else + mprintf("Warning: Overriding previous state.\n"); + state_num = (int)(state - States_.begin()); + } else state_num = (int)(state - States_.begin()); } } @@ -200,4 +214,26 @@ Analysis::RetType Analysis_State::Analyze() { trans->second.NormalizeCurve(); } return Analysis::OK; -} +} + +// ----- Analysis_State::StateType --------------------------------------------- +size_t Analysis_State::StateType::Nframes() const +{ + if (Sets_.empty()) return 0; + size_t nframes = Sets_[0]->Size(); + for (unsigned int idx = 1; idx < Sets_.size(); idx++) + if (nframes != Sets_[idx]->Size()) { + mprintf("Warning: State '%s', set '%s' has differing # frames (%zu).", + id(), Sets_[idx]->legend(), Sets_[idx]->Size()); + nframes = std::min( nframes, Sets_[idx]->Size() ); + mprintf(" Only using %zu frames.\n", nframes); + } + return nframes; +} + +void Analysis_State::StateType::PrintState() const { + mprintf("%s", id()); + for (unsigned int idx = 0; idx != Sets_.size(); idx++) + mprintf(" {%.4f <= %s < %.4f}", Min_[idx], Sets_[idx]->legend(), Max_[idx]); + mprintf("\n"); +} diff --git a/src/Analysis_State.h b/src/Analysis_State.h index e8c3f51f1a..461fd671cf 100644 --- a/src/Analysis_State.h +++ b/src/Analysis_State.h @@ -2,6 +2,7 @@ #define INC_ANALYSIS_STATE_H #include #include "Analysis.h" +#include "Array1D.h" #include "DataSet_double.h" /// Analyze transitions between states class Analysis_State : public Analysis { @@ -17,19 +18,37 @@ class Analysis_State : public Analysis { /// Hold the definition of a state and associated data class StateType { public: - StateType() : set_(0), min_(0.0), max_(0.0) {} - StateType(std::string const& i, DataSet_1D* d, double m, double x) : - id_(i), set_(d), min_(m), max_(x) {} + StateType() {} + StateType(std::string const& i) : id_(i) {} + /// Add criterion for state + void AddCriterion(DataSet_1D* d, double m, double x) { + Sets_.push_back( d ); + Min_.push_back( m ); + Max_.push_back( x ); + } + /// \return smallest number of frames among all criterion DataSets + size_t Nframes() const; + /// \return State ID std::string const& ID() const { return id_; } + /// \return State ID as const char* const char* id() const { return id_.c_str(); } - DataSet_1D const& DS() const { return *set_; } - double Min() const { return min_; } - double Max() const { return max_; } + /// \return true if given frame satifies all criteria + bool InState(int n) const { + for (unsigned int idx = 0; idx != Sets_.size(); idx++) { + double dval = Sets_[idx]->Dval(n); + if (dval < Min_[idx] || dval >= Max_[idx]) // TODO periodic + return false; + } + return true; + } + /// Print state info to stdout + void PrintState() const; private: + typedef std::vector Darray; std::string id_; - DataSet_1D* set_; - double min_; - double max_; + Array1D Sets_; ///< DataSets used to determine if we are in State + Darray Min_; ///< Above this value we are in state. + Darray Max_; ///< Below this value we are in state. }; /// Hold information about a transition from state0 to state1 class Transition { From 0975c80459986f702cb039844d6c74fb810ec997 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 12:56:21 -0400 Subject: [PATCH 02/20] DRR - Cpptraj: Calculate state totals. --- src/Analysis_State.cpp | 28 +++++++++++++++++++++++----- src/Analysis_State.h | 3 ++- 2 files changed, 25 insertions(+), 6 deletions(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index c6c2f9a15d..3445497fc6 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -2,6 +2,15 @@ #include "Analysis_State.h" #include "CpptrajStdio.h" +Analysis_State::Analysis_State() : + state_data_(0), + curveOut_(0), + stateOut_(0), + transOut_(0), + countOut_(0), + debug_(0) +{} + void Analysis_State::Help() const { mprintf("\t{state ,,,} [out ] [name ]\n" "\t[curveout ] [stateout ] [transout ]\n" @@ -15,10 +24,12 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set masterDSL_ = setup.DslPtr(); DataFile* outfile = setup.DFL().AddDataFile( analyzeArgs.GetStringKey("out"), analyzeArgs ); curveOut_ = setup.DFL().AddDataFile( analyzeArgs.GetStringKey("curveout"), analyzeArgs ); - stateOut_ = setup.DFL().AddCpptrajFile( analyzeArgs.GetStringKey("stateout"), "State Output", - DataFileList::TEXT, true); - transOut_ = setup.DFL().AddCpptrajFile( analyzeArgs.GetStringKey("transout"), "Transitions Output", - DataFileList::TEXT, true); + stateOut_ = setup.DFL().AddCpptrajFile(analyzeArgs.GetStringKey("stateout"), + "State Lifetime Output", DataFileList::TEXT, true); + transOut_ = setup.DFL().AddCpptrajFile(analyzeArgs.GetStringKey("transout"), + "Transitions Output", DataFileList::TEXT, true); + countOut_ = setup.DFL().AddCpptrajFile(analyzeArgs.GetStringKey("countout"), + "State Count Output", DataFileList::TEXT, true); normalize_ = analyzeArgs.hasKey("norm"); // Get definitions of states if present. // Define states as 'state <#>,,,' @@ -81,8 +92,9 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set mprintf("\tStates vs time output to file '%s'\n", outfile->DataFilename().full()); if (curveOut_ != 0) mprintf("\tCurves output to file '%s'\n", curveOut_->DataFilename().full()); - mprintf("\tState output to file '%s'\n", stateOut_->Filename().full()); + mprintf("\tState lifetime output to file '%s'\n", stateOut_->Filename().full()); mprintf("\tTransitions output to file '%s'\n", transOut_->Filename().full()); + mprintf("\tState counts output to file '%s'\n", countOut_->Filename().full()); if (normalize_) mprintf("\tCurves will be normalized to 1.0\n"); @@ -129,6 +141,7 @@ Analysis::RetType Analysis_State::Analyze() { int last_state = -2; size_t last_State_Start = 0; size_t final_frame = nframes - 1; + std::vector stateFrames(States_.size()+1, 0); for (size_t frm = 0; frm != nframes; frm++) { // Determine which state we are in. int state_num = -1; @@ -145,6 +158,7 @@ Analysis::RetType Analysis_State::Analyze() { } } state_data_->Add( frm, &state_num ); + stateFrames[state_num+1]++; // Determine if there has been a transition. if (last_state == -2) @@ -195,6 +209,10 @@ Analysis::RetType Analysis_State::Analyze() { trans->first.second, trans->second.Max(), trans->second.Sum(), trans->second.Nlifetimes(), trans->second.Avg()); } + countOut_->Printf("%-8s %12s %12s %s\n", "#Index", "Count", "Frac", "State"); + for (int idx = 0; idx != (int)stateFrames.size(); idx++) + countOut_->Printf("%-8i %12i %12.4f %s\n", idx-1, stateFrames[idx], + (double)stateFrames[idx]/(double)nframes, stateName(idx-1)); stateOut_->Printf("%-8s %12s %12s %12s %s\n", "#Index", "N", "Average", "Max", "State"); for (int idx = 0; idx != (int)Status.size(); idx++) stateOut_->Printf("%-8i %12i %12.4f %12i %s\n", idx-1, Status[idx].Nlifetimes(), diff --git a/src/Analysis_State.h b/src/Analysis_State.h index 461fd671cf..8216998740 100644 --- a/src/Analysis_State.h +++ b/src/Analysis_State.h @@ -7,7 +7,7 @@ /// Analyze transitions between states class Analysis_State : public Analysis { public: - Analysis_State() : state_data_(0), curveOut_(0), stateOut_(0), transOut_(0), debug_(0) {} + Analysis_State(); DispatchObject* Alloc() const { return (DispatchObject*)new Analysis_State(); } void Help() const; @@ -113,6 +113,7 @@ class Analysis_State : public Analysis { DataFile* curveOut_; CpptrajFile* stateOut_; CpptrajFile* transOut_; + CpptrajFile* countOut_; int debug_; bool normalize_; }; From 8354a4b94750098268fabcf2467804262d060591 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 13:22:17 -0400 Subject: [PATCH 03/20] DRR - Cpptraj: Allow states to be defined multiple times; this effectively acts as an 'OR' operator in addition to the implicit 'AND' --- src/Analysis_State.cpp | 29 ++++++++++++++++++++--------- src/Analysis_State.h | 19 +++++++++++++++---- 2 files changed, 35 insertions(+), 13 deletions(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index 3445497fc6..50e18343c3 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -52,7 +52,17 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set mprinterr("Error: In state argument, could not get ID.\n"); return Analysis::ERR; } - StateType currentState( state_id ); + // Determine state number + int state_num = -1; + IdxMapType::const_iterator it = NameMap_.find( state_id ); + if (it == NameMap_.end()) { + state_num = NameMap_.size(); + NameMap_.insert( IdxPair( state_id, state_num ) ); + // For mapping state number back to name + StateNames_.push_back( state_id ); + } else + state_num = it->second; + StateType currentState( state_id, state_num ); int nSets = (argtmp.Nargs()-1) / 3; for (int setArg = 0; setArg < nSets; setArg++) { DataSet* ds = setup.DSL().GetDataSet( argtmp.GetStringNext() ); @@ -81,7 +91,7 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set if (state_data_ == 0) return Analysis::ERR; if (outfile != 0) outfile->AddDataSet( state_data_ ); - mprintf(" STATE: The following states have been set up:\n"); + mprintf(" STATE: The following state definitions have been set up:\n"); for (StateArray::const_iterator state = States_.begin(); state != States_.end(); ++state) { mprintf("\t%u: ", state - States_.begin()); @@ -105,7 +115,7 @@ std::string Analysis_State::UNDEFINED_ = "Undefined"; std::string const& Analysis_State::StateName(int idx) const { if (idx > -1) - return States_[idx].ID(); + return StateNames_[idx]; else // Should never be called with any other negative number but -1 return UNDEFINED_; } @@ -128,8 +138,9 @@ Analysis::RetType Analysis_State::Analyze() { if (nframes < 1) return Analysis::ERR; std::vector Status; - Status.reserve( States_.size() + 1 ); // +1 for state -1, undefined - for (int i = 0; i != (int)States_.size() + 1; i++) { + int numStates = (int)NameMap_.size() + 1; // + 1 for state -1 (undefined) + Status.reserve( numStates ); + for (int i = 0; i != numStates; i++) { DataSet* ds = masterDSL_->AddSet(DataSet::DOUBLE, MetaData(state_data_->Meta().Name(), "sCurve", i)); if (ds == 0) return Analysis::ERR; @@ -141,7 +152,7 @@ Analysis::RetType Analysis_State::Analyze() { int last_state = -2; size_t last_State_Start = 0; size_t final_frame = nframes - 1; - std::vector stateFrames(States_.size()+1, 0); + std::vector stateFrames(numStates, 0); for (size_t frm = 0; frm != nframes; frm++) { // Determine which state we are in. int state_num = -1; @@ -152,9 +163,9 @@ Analysis::RetType Analysis_State::Analyze() { mprintf("Warning: Frame %zu already defined as state '%s', also matches state '%s'.\n", frm, States_[state_num].id(), state->id()); mprintf("Warning: Overriding previous state.\n"); - state_num = (int)(state - States_.begin()); + state_num = state->Num(); } else - state_num = (int)(state - States_.begin()); + state_num = state->Num(); } } state_data_->Add( frm, &state_num ); @@ -250,7 +261,7 @@ size_t Analysis_State::StateType::Nframes() const } void Analysis_State::StateType::PrintState() const { - mprintf("%s", id()); + mprintf("%s (%i)", id(), num_); for (unsigned int idx = 0; idx != Sets_.size(); idx++) mprintf(" {%.4f <= %s < %.4f}", Min_[idx], Sets_[idx]->legend(), Max_[idx]); mprintf("\n"); diff --git a/src/Analysis_State.h b/src/Analysis_State.h index 8216998740..c194391792 100644 --- a/src/Analysis_State.h +++ b/src/Analysis_State.h @@ -15,11 +15,12 @@ class Analysis_State : public Analysis { Analysis::RetType Analyze(); private: // ----- PRIVATE CLASS DEFINITIONS --------------------------------------------- + // ------------------------------------------- /// Hold the definition of a state and associated data class StateType { public: - StateType() {} - StateType(std::string const& i) : id_(i) {} + StateType() : num_(-1) {} + StateType(std::string const& i, int n) : id_(i), num_(n) {} /// Add criterion for state void AddCriterion(DataSet_1D* d, double m, double x) { Sets_.push_back( d ); @@ -32,6 +33,8 @@ class Analysis_State : public Analysis { std::string const& ID() const { return id_; } /// \return State ID as const char* const char* id() const { return id_.c_str(); } + /// \return Unique state number + int Num() const { return num_; } /// \return true if given frame satifies all criteria bool InState(int n) const { for (unsigned int idx = 0; idx != Sets_.size(); idx++) { @@ -47,9 +50,11 @@ class Analysis_State : public Analysis { typedef std::vector Darray; std::string id_; Array1D Sets_; ///< DataSets used to determine if we are in State - Darray Min_; ///< Above this value we are in state. - Darray Max_; ///< Below this value we are in state. + Darray Min_; ///< Above this value we are in state. + Darray Max_; ///< Below this value we are in state. + int num_; ///< Unique state identifier }; + // ------------------------------------------- /// Hold information about a transition from state0 to state1 class Transition { public: @@ -95,18 +100,24 @@ class Analysis_State : public Analysis { int Nlifetimes_; ///< Number of state0 lifetimes before going to state1 DataSet_double* curve_; ///< Lifetime curve for state0->state1 }; + // ------------------------------------------- typedef std::vector StateArray; typedef std::pair StatePair; typedef std::map TransMapType; typedef std::pair TransPair; typedef std::vector Iarray; + typedef std::pair IdxPair; + typedef std::map IdxMapType; + typedef std::vector Sarray; std::string const& StateName(int) const; const char* stateName(int i) const { return StateName(i).c_str(); } static std::string UNDEFINED_; StateArray States_; + IdxMapType NameMap_; + Sarray StateNames_; TransMapType TransitionMap_; DataSet* state_data_; DataSetList* masterDSL_; From ae5124de36cab8b4a16d6c22d778f392ea686c2f Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 13:54:39 -0400 Subject: [PATCH 04/20] DRR - Cpptraj: Start adding state-specific data sets --- src/Analysis_State.cpp | 15 ++++++++++++++- src/Analysis_State.h | 1 + 2 files changed, 15 insertions(+), 1 deletion(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index 50e18343c3..fd58a70ee6 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -18,6 +18,7 @@ void Analysis_State::Help() const { " will be assigned state <#>.\n"); } +// Analysis_State::Setup() Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& setup, int debugIn) { debug_ = debugIn; @@ -90,6 +91,13 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set state_data_ = setup.DSL().AddSet(DataSet::INTEGER, analyzeArgs.GetStringKey("name"), "State"); if (state_data_ == 0) return Analysis::ERR; if (outfile != 0) outfile->AddDataSet( state_data_ ); + // Set up state data sets + Dimension Xdim(-1, 1, "Index"); + DataSet* ds = setup.DSL().AddSet(DataSet::INTEGER, + MetaData(state_data_->Meta().Name(), "Count")); + if (ds == 0) return Analysis::ERR; + ds->SetDim(0, Xdim); + state_counts_ = ds; mprintf(" STATE: The following state definitions have been set up:\n"); for (StateArray::const_iterator state = States_.begin(); state != States_.end(); ++state) @@ -111,8 +119,10 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set return Analysis::OK; } +/** Default state name for state -1 */ std::string Analysis_State::UNDEFINED_ = "Undefined"; +// Analysis_State::StateName() std::string const& Analysis_State::StateName(int idx) const { if (idx > -1) return StateNames_[idx]; @@ -120,6 +130,7 @@ std::string const& Analysis_State::StateName(int idx) const { return UNDEFINED_; } +// Analysis_State::Analyze() Analysis::RetType Analysis_State::Analyze() { // Only process as much data as is in the smallest data set. size_t nframes = States_.front().Nframes(); @@ -221,9 +232,11 @@ Analysis::RetType Analysis_State::Analyze() { trans->second.Nlifetimes(), trans->second.Avg()); } countOut_->Printf("%-8s %12s %12s %s\n", "#Index", "Count", "Frac", "State"); - for (int idx = 0; idx != (int)stateFrames.size(); idx++) + for (int idx = 0; idx != (int)stateFrames.size(); idx++) { countOut_->Printf("%-8i %12i %12.4f %s\n", idx-1, stateFrames[idx], (double)stateFrames[idx]/(double)nframes, stateName(idx-1)); + state_counts_->Add(idx, &(stateFrames[idx])); + } stateOut_->Printf("%-8s %12s %12s %12s %s\n", "#Index", "N", "Average", "Max", "State"); for (int idx = 0; idx != (int)Status.size(); idx++) stateOut_->Printf("%-8i %12i %12.4f %12i %s\n", idx-1, Status[idx].Nlifetimes(), diff --git a/src/Analysis_State.h b/src/Analysis_State.h index c194391792..772493a9a1 100644 --- a/src/Analysis_State.h +++ b/src/Analysis_State.h @@ -120,6 +120,7 @@ class Analysis_State : public Analysis { Sarray StateNames_; TransMapType TransitionMap_; DataSet* state_data_; + DataSet* state_counts_; DataSetList* masterDSL_; DataFile* curveOut_; CpptrajFile* stateOut_; From f4075687308b2489ce11da45f89a121e628d6275 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 14:19:05 -0400 Subject: [PATCH 05/20] DRR - Cpptraj: Add more data sets for state-related data --- src/Analysis_State.cpp | 31 ++++++++++++++++++++++++++++++- src/Analysis_State.h | 5 +++++ 2 files changed, 35 insertions(+), 1 deletion(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index fd58a70ee6..4f61fffeb5 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -98,6 +98,26 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set if (ds == 0) return Analysis::ERR; ds->SetDim(0, Xdim); state_counts_ = ds; + ds = setup.DSL().AddSet(DataSet::DOUBLE, MetaData(state_data_->Meta().Name(), "Frac")); + if (ds == 0) return Analysis::ERR; + ds->SetDim(0, Xdim); + state_fracs_ = ds; + ds = setup.DSL().AddSet(DataSet::INTEGER, MetaData(state_data_->Meta().Name(), "Nlifetimes")); + if (ds == 0) return Analysis::ERR; + ds->SetDim(0, Xdim); + state_lifetimes_ = ds; + ds = setup.DSL().AddSet(DataSet::DOUBLE, MetaData(state_data_->Meta().Name(), "Avglife")); + if (ds == 0) return Analysis::ERR; + ds->SetDim(0, Xdim); + state_avglife_ = ds; + ds = setup.DSL().AddSet(DataSet::INTEGER, MetaData(state_data_->Meta().Name(), "Maxlife")); + if (ds == 0) return Analysis::ERR; + ds->SetDim(0, Xdim); + state_maxlife_ = ds; + ds = setup.DSL().AddSet(DataSet::STRING, MetaData(state_data_->Meta().Name(), "Name")); + if (ds == 0) return Analysis::ERR; + ds->SetDim(0, Xdim); + state_names_ = ds; mprintf(" STATE: The following state definitions have been set up:\n"); for (StateArray::const_iterator state = States_.begin(); state != States_.end(); ++state) @@ -232,10 +252,19 @@ Analysis::RetType Analysis_State::Analyze() { trans->second.Nlifetimes(), trans->second.Avg()); } countOut_->Printf("%-8s %12s %12s %s\n", "#Index", "Count", "Frac", "State"); - for (int idx = 0; idx != (int)stateFrames.size(); idx++) { + for (int idx = 0; idx != numStates; idx++) { countOut_->Printf("%-8i %12i %12.4f %s\n", idx-1, stateFrames[idx], (double)stateFrames[idx]/(double)nframes, stateName(idx-1)); state_counts_->Add(idx, &(stateFrames[idx])); + double dval = (double)stateFrames[idx]/(double)nframes; + state_fracs_->Add(idx, &dval); + state_names_->Add(idx, stateName(idx-1)); + int ival = Status[idx].Nlifetimes(); + state_lifetimes_->Add(idx, &ival); + dval = Status[idx].Avg(); + state_avglife_->Add(idx, &dval); + ival = Status[idx].Max(); + state_maxlife_->Add(idx, &ival); } stateOut_->Printf("%-8s %12s %12s %12s %s\n", "#Index", "N", "Average", "Max", "State"); for (int idx = 0; idx != (int)Status.size(); idx++) diff --git a/src/Analysis_State.h b/src/Analysis_State.h index 772493a9a1..49145ba09c 100644 --- a/src/Analysis_State.h +++ b/src/Analysis_State.h @@ -121,6 +121,11 @@ class Analysis_State : public Analysis { TransMapType TransitionMap_; DataSet* state_data_; DataSet* state_counts_; + DataSet* state_fracs_; + DataSet* state_lifetimes_; + DataSet* state_avglife_; + DataSet* state_maxlife_; + DataSet* state_names_; DataSetList* masterDSL_; DataFile* curveOut_; CpptrajFile* stateOut_; From e6e796ddb83a98a3ceae965dff31b13d573dc93b Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 14:32:38 -0400 Subject: [PATCH 06/20] DRR - Cpptraj: Finish transitioning state-related data to data sets. --- src/Analysis_State.cpp | 45 ++++++++++++++++++++++++++++++------------ src/Analysis_State.h | 4 ++-- 2 files changed, 34 insertions(+), 15 deletions(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index 4f61fffeb5..fd3e806122 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -4,10 +4,16 @@ Analysis_State::Analysis_State() : state_data_(0), + state_counts_(0), + state_fracs_(0), + state_lifetimes_(0), + state_avglife_(0), + state_maxlife_(0), + state_names_(0), curveOut_(0), - stateOut_(0), - transOut_(0), + lifeOut_(0), countOut_(0), + transOut_(0), debug_(0) {} @@ -25,12 +31,10 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set masterDSL_ = setup.DslPtr(); DataFile* outfile = setup.DFL().AddDataFile( analyzeArgs.GetStringKey("out"), analyzeArgs ); curveOut_ = setup.DFL().AddDataFile( analyzeArgs.GetStringKey("curveout"), analyzeArgs ); - stateOut_ = setup.DFL().AddCpptrajFile(analyzeArgs.GetStringKey("stateout"), - "State Lifetime Output", DataFileList::TEXT, true); + lifeOut_ = setup.DFL().AddDataFile( analyzeArgs.GetStringKey("stateout"), analyzeArgs ); + countOut_ = setup.DFL().AddDataFile( analyzeArgs.GetStringKey("countout"), analyzeArgs ); transOut_ = setup.DFL().AddCpptrajFile(analyzeArgs.GetStringKey("transout"), "Transitions Output", DataFileList::TEXT, true); - countOut_ = setup.DFL().AddCpptrajFile(analyzeArgs.GetStringKey("countout"), - "State Count Output", DataFileList::TEXT, true); normalize_ = analyzeArgs.hasKey("norm"); // Get definitions of states if present. // Define states as 'state <#>,,,' @@ -119,6 +123,18 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set ds->SetDim(0, Xdim); state_names_ = ds; + if (countOut_ != 0) { + countOut_->AddDataSet( state_counts_ ); + countOut_->AddDataSet( state_fracs_ ); + if (lifeOut_ != countOut_) countOut_->AddDataSet( state_names_ ); + } + if (lifeOut_ != 0) { + lifeOut_->AddDataSet( state_lifetimes_ ); + lifeOut_->AddDataSet( state_avglife_ ); + lifeOut_->AddDataSet( state_maxlife_ ); + lifeOut_->AddDataSet( state_names_ ); + } + mprintf(" STATE: The following state definitions have been set up:\n"); for (StateArray::const_iterator state = States_.begin(); state != States_.end(); ++state) { @@ -130,9 +146,11 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set mprintf("\tStates vs time output to file '%s'\n", outfile->DataFilename().full()); if (curveOut_ != 0) mprintf("\tCurves output to file '%s'\n", curveOut_->DataFilename().full()); - mprintf("\tState lifetime output to file '%s'\n", stateOut_->Filename().full()); + if (lifeOut_ != 0) + mprintf("\tState lifetime output to file '%s'\n", lifeOut_->DataFilename().full()); + if (countOut_ != 0) + mprintf("\tState counts output to file '%s'\n", countOut_->DataFilename().full()); mprintf("\tTransitions output to file '%s'\n", transOut_->Filename().full()); - mprintf("\tState counts output to file '%s'\n", countOut_->Filename().full()); if (normalize_) mprintf("\tCurves will be normalized to 1.0\n"); @@ -172,6 +190,7 @@ Analysis::RetType Analysis_State::Analyze() { int numStates = (int)NameMap_.size() + 1; // + 1 for state -1 (undefined) Status.reserve( numStates ); for (int i = 0; i != numStates; i++) { + // TODO can this be done in Setup() DataSet* ds = masterDSL_->AddSet(DataSet::DOUBLE, MetaData(state_data_->Meta().Name(), "sCurve", i)); if (ds == 0) return Analysis::ERR; @@ -251,10 +270,10 @@ Analysis::RetType Analysis_State::Analyze() { trans->first.second, trans->second.Max(), trans->second.Sum(), trans->second.Nlifetimes(), trans->second.Avg()); } - countOut_->Printf("%-8s %12s %12s %s\n", "#Index", "Count", "Frac", "State"); + //countOut_->Printf("%-8s %12s %12s %s\n", "#Index", "Count", "Frac", "State"); for (int idx = 0; idx != numStates; idx++) { - countOut_->Printf("%-8i %12i %12.4f %s\n", idx-1, stateFrames[idx], - (double)stateFrames[idx]/(double)nframes, stateName(idx-1)); + /*countOut_->Printf("%-8i %12i %12.4f %s\n", idx-1, stateFrames[idx], + (double)stateFrames[idx]/(double)nframes, stateName(idx-1));*/ state_counts_->Add(idx, &(stateFrames[idx])); double dval = (double)stateFrames[idx]/(double)nframes; state_fracs_->Add(idx, &dval); @@ -266,10 +285,10 @@ Analysis::RetType Analysis_State::Analyze() { ival = Status[idx].Max(); state_maxlife_->Add(idx, &ival); } - stateOut_->Printf("%-8s %12s %12s %12s %s\n", "#Index", "N", "Average", "Max", "State"); + /*stateOut_->Printf("%-8s %12s %12s %12s %s\n", "#Index", "N", "Average", "Max", "State"); for (int idx = 0; idx != (int)Status.size(); idx++) stateOut_->Printf("%-8i %12i %12.4f %12i %s\n", idx-1, Status[idx].Nlifetimes(), - Status[idx].Avg(), Status[idx].Max(), stateName(idx-1)); + Status[idx].Avg(), Status[idx].Max(), stateName(idx-1));*/ transOut_->Printf("%-12s %12s %12s %s\n", "#N", "Average", "Max", "Transition"); for (TransMapType::const_iterator trans = TransitionMap_.begin(); trans != TransitionMap_.end(); ++trans) diff --git a/src/Analysis_State.h b/src/Analysis_State.h index 49145ba09c..b32f593926 100644 --- a/src/Analysis_State.h +++ b/src/Analysis_State.h @@ -128,9 +128,9 @@ class Analysis_State : public Analysis { DataSet* state_names_; DataSetList* masterDSL_; DataFile* curveOut_; - CpptrajFile* stateOut_; + DataFile* lifeOut_; + DataFile* countOut_; CpptrajFile* transOut_; - CpptrajFile* countOut_; int debug_; bool normalize_; }; From 8fae74a2fb30895a7beda34bb4ee8c1ff8b53fb9 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 14:56:37 -0400 Subject: [PATCH 07/20] DRR - Cpptraj: Correctly determine length of final lifetime when transition has not occurred. --- src/Analysis_State.cpp | 14 ++++++++++++-- 1 file changed, 12 insertions(+), 2 deletions(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index fd3e806122..5932809273 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -243,13 +243,23 @@ Analysis::RetType Analysis_State::Analyze() { } else // Update previous transition entry->second.Update(length); - } + // If final frame, this new state length has to be 1 and may not be + // representative of a "real" lifetime - print a warning and ignore. + if (frm == final_frame) + mprintf("Warning: Transition occurred on final frame. Final frame will" + "Warning: not count as a lifetime for state %i '%s' since it" + "Warning: likely does not represent a \"real\" lifetime.\n", + state_num, stateName(state_num)); + } else if (frm == final_frame) + // No transition - include last frame in current state length + length++; // Update single state information. + //mprintf("DEBUG: Calling update at %zu for state %i with length %i\n", frm+1, last_state, length); Status[last_state + 1].Update(length); last_State_Start = frm; last_state = state_num; } - } + } // END loop over frames // DEBUG: Print single state info. if (debug_ > 0) { From fc120b9818abc18e9929e72885bb803a92b3e83f Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 14:57:34 -0400 Subject: [PATCH 08/20] DRR - Cpptraj: Update for fixed final frame lifetime calc and new format --- test/Test_CalcState/RunTest.sh | 1 + test/Test_CalcState/curve.dat.save | 2 +- test/Test_CalcState/states.dat.save | 8 ++++---- 3 files changed, 6 insertions(+), 5 deletions(-) diff --git a/test/Test_CalcState/RunTest.sh b/test/Test_CalcState/RunTest.sh index 91974125ee..4b8f66eccc 100755 --- a/test/Test_CalcState/RunTest.sh +++ b/test/Test_CalcState/RunTest.sh @@ -14,6 +14,7 @@ calcstate name State \ out state.dat \ curveout curve.dat \ stateout states.dat \ + countout states.dat \ transout trans.dat EOF RunCpptraj "Calcstate test" diff --git a/test/Test_CalcState/curve.dat.save b/test/Test_CalcState/curve.dat.save index 4d4beec940..ff3bb6a6ed 100644 --- a/test/Test_CalcState/curve.dat.save +++ b/test/Test_CalcState/curve.dat.save @@ -1,4 +1,4 @@ #Frame Undefined C0 C1 C0->Undefined Undefined->C0 C0->C1 1 2.0000 3.0000 1.0000 2.0000 2.0000 1.0000 2 0.0000 1.0000 1.0000 0.0000 0.0000 1.0000 - 3 0.0000 1.0000 0.0000 0.0000 0.0000 1.0000 + 3 0.0000 1.0000 1.0000 0.0000 0.0000 1.0000 diff --git a/test/Test_CalcState/states.dat.save b/test/Test_CalcState/states.dat.save index ca67265101..fa464dfaac 100644 --- a/test/Test_CalcState/states.dat.save +++ b/test/Test_CalcState/states.dat.save @@ -1,4 +1,4 @@ -#Index N Average Max State --1 2 1.0000 1 Undefined -0 3 1.6667 3 C0 -1 1 2.0000 2 C1 +#Index State[Count] State[Frac] State[Nlifetimes] State[Avglife] State[Maxlife] State[Name] + -1 2 0.2000 2 1.0000 1 Undefined + 0 5 0.5000 3 1.6667 3 C0 + 1 3 0.3000 1 3.0000 3 C1 From 3af6e5e8ed3035d64e06d18cb3a77d9e07529963 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 14:58:12 -0400 Subject: [PATCH 09/20] DRR - Cpptraj: Minor version bump since 'calcstate' output format has changed slightly and certain parts are not automatically written to stdout anymore. --- src/Version.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/Version.h b/src/Version.h index f0bb225c87..394cb4fba0 100644 --- a/src/Version.h +++ b/src/Version.h @@ -20,5 +20,5 @@ * Whenever a number that precedes is incremented, all subsequent * numbers should be reset to 0. */ -#define CPPTRAJ_INTERNAL_VERSION "V4.6.2" +#define CPPTRAJ_INTERNAL_VERSION "V4.7.0" #endif From 7b59fa9058c17d108e64d312cd4739c654fc1899 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 15:45:37 -0400 Subject: [PATCH 10/20] DRR - Cpptraj: Remove old code. --- src/Analysis_State.cpp | 9 ++------- 1 file changed, 2 insertions(+), 7 deletions(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index 5932809273..3ef290994e 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -280,10 +280,8 @@ Analysis::RetType Analysis_State::Analyze() { trans->first.second, trans->second.Max(), trans->second.Sum(), trans->second.Nlifetimes(), trans->second.Avg()); } - //countOut_->Printf("%-8s %12s %12s %s\n", "#Index", "Count", "Frac", "State"); + // Calculate state-dependent stuff for (int idx = 0; idx != numStates; idx++) { - /*countOut_->Printf("%-8i %12i %12.4f %s\n", idx-1, stateFrames[idx], - (double)stateFrames[idx]/(double)nframes, stateName(idx-1));*/ state_counts_->Add(idx, &(stateFrames[idx])); double dval = (double)stateFrames[idx]/(double)nframes; state_fracs_->Add(idx, &dval); @@ -295,10 +293,7 @@ Analysis::RetType Analysis_State::Analyze() { ival = Status[idx].Max(); state_maxlife_->Add(idx, &ival); } - /*stateOut_->Printf("%-8s %12s %12s %12s %s\n", "#Index", "N", "Average", "Max", "State"); - for (int idx = 0; idx != (int)Status.size(); idx++) - stateOut_->Printf("%-8i %12i %12.4f %12i %s\n", idx-1, Status[idx].Nlifetimes(), - Status[idx].Avg(), Status[idx].Max(), stateName(idx-1));*/ + // Calculate transitions transOut_->Printf("%-12s %12s %12s %s\n", "#N", "Average", "Max", "Transition"); for (TransMapType::const_iterator trans = TransitionMap_.begin(); trans != TransitionMap_.end(); ++trans) From ae55b9e0d144a49715999c8ba6e40946e28073b5 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 15:46:02 -0400 Subject: [PATCH 11/20] DRR - Cpptraj: New 'dataset' command mode: invert --- src/Exec_DataSetCmd.cpp | 75 ++++++++++++++++++++++++++++++++++++++--- src/Exec_DataSetCmd.h | 1 + 2 files changed, 72 insertions(+), 4 deletions(-) diff --git a/src/Exec_DataSetCmd.cpp b/src/Exec_DataSetCmd.cpp index 55ab0184a6..2a89713518 100644 --- a/src/Exec_DataSetCmd.cpp +++ b/src/Exec_DataSetCmd.cpp @@ -86,6 +86,9 @@ Exec::RetType Exec_DataSetCmd::Execute(CpptrajState& State, ArgList& argIn) { } else if (argIn.hasKey("dim")) { // Modify dimension of set(s) err = ChangeDim(State, argIn); // --------------------------------------------- + } else if (argIn.hasKey("invert")) { // Invert set(s) X/Y, create new sets + err = InvertSets(State, argIn); + // --------------------------------------------- } else { // Default: change mode/type for one or more sets. err = ChangeModeType(State, argIn); } @@ -622,10 +625,7 @@ Exec::RetType Exec_DataSetCmd::Concatenate(CpptrajState& State, ArgList& argIn) double XY[2]; for (DataSetList::const_iterator ds = dsl.begin(); ds != dsl.end(); ++ds) { - if ( (*ds)->Type() != DataSet::INTEGER && - (*ds)->Type() != DataSet::DOUBLE && - (*ds)->Type() != DataSet::FLOAT && - (*ds)->Type() != DataSet::XYMESH ) + if ( (*ds)->Group() != DataSet::SCALAR_1D ) { mprintf("Warning: '%s': Concatenation only supported for 1D scalar data sets.\n", (*ds)->legend()); @@ -784,3 +784,70 @@ Exec::RetType Exec_DataSetCmd::ChangeModeType(CpptrajState const& State, ArgList } return CpptrajState::OK; } + +/** Syntax: dataset invert ... name */ +Exec::RetType Exec_DataSetCmd::InvertSets(CpptrajState& State, ArgList& argIn) { + DataSetList& DSL = State.DSL(); + // Get keywords + std::string dsname = argIn.GetStringKey("name"); + if (dsname.empty()) { + mprinterr("Error: 'invert' requires that 'name ' be specified.\n"); + return CpptrajState::ERR; + } + DataFile* outfile = State.DFL().AddDataFile( argIn.GetStringKey("out"), argIn ); + // TODO determine type some other way + DataSet::DataType outtype = DataSet::DOUBLE; + // Get input DataSets + std::vector input_sets; + std::string dsarg = argIn.GetStringNext(); + while (!dsarg.empty()) { + DataSetList sets = DSL.GetMultipleSets( dsarg ); + for (DataSetList::const_iterator ds = sets.begin(); ds != sets.end(); ++ds) + { + if ( (*ds)->Group() != DataSet::SCALAR_1D ) { + mprintf("Warning: '%s': Inversion only supported for 1D scalar data sets.\n", + (*ds)->legend()); + } else { + if (!input_sets.empty()) { + if ( (*ds)->Size() != input_sets.back()->Size() ) { + mprinterr("Error: Set '%s' has different size (%zu) than previous set (%zu)\n", + (*ds)->legend(), (*ds)->Size(), input_sets.back()->Size()); + return CpptrajState::ERR; + } + } + input_sets.push_back( (DataSet_1D*)*ds ); + } + } + dsarg = argIn.GetStringNext(); + } + if (input_sets.empty()) { + mprinterr("Error: No sets selected.\n"); + return CpptrajState::ERR; + } + // Need an output data set for each point in input sets + std::vector output_sets; + int column = 1; + for (int idx = 0; idx != (int)input_sets[0]->Size(); idx++, column++) { + DataSet* ds = DSL.AddSet(outtype, MetaData(dsname, column)); + if (ds == 0) return CpptrajState::ERR; + output_sets.push_back( ds ); + if (outfile != 0) outfile->AddDataSet( ds ); + } + // Create a data set containing names of each input data set + DataSet* nameset = DSL.AddSet(DataSet::STRING, MetaData(dsname, column)); + if (nameset == 0) return CpptrajState::ERR; + if (outfile != 0) outfile->AddDataSet( nameset ); + // Populate output data sets + for (int jdx = 0; jdx != (int)input_sets.size(); jdx++) + { + DataSet_1D const& INP = static_cast( *(input_sets[jdx]) ); + nameset->Add( jdx, INP.legend() ); + for (unsigned int idx = 0; idx != INP.Size(); idx++) + { + double dval = INP.Dval( idx ); + output_sets[idx]->Add( jdx, &dval ); + } + } + + return CpptrajState::OK; +} diff --git a/src/Exec_DataSetCmd.h b/src/Exec_DataSetCmd.h index d3c5d9867a..291ee29638 100644 --- a/src/Exec_DataSetCmd.h +++ b/src/Exec_DataSetCmd.h @@ -30,5 +30,6 @@ class Exec_DataSetCmd : public Exec { static void Help_ChangeDim(); RetType ChangeDim(CpptrajState const&, ArgList&); RetType ChangeModeType(CpptrajState const&, ArgList&); + RetType InvertSets(CpptrajState&, ArgList&); }; #endif From 91281a1889bd067f56ae13368f3ef1d88ad91c68 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Thu, 23 Aug 2018 16:07:10 -0400 Subject: [PATCH 12/20] DRR - Cpptraj: Add keyword for specifying set containing legends. --- src/Exec_DataSetCmd.cpp | 28 ++++++++++++++++++++++++++-- 1 file changed, 26 insertions(+), 2 deletions(-) diff --git a/src/Exec_DataSetCmd.cpp b/src/Exec_DataSetCmd.cpp index 2a89713518..7bd1999de7 100644 --- a/src/Exec_DataSetCmd.cpp +++ b/src/Exec_DataSetCmd.cpp @@ -4,6 +4,7 @@ #include "DataSet_1D.h" #include "DataSet_MatrixDbl.h" #include "DataSet_Vector.h" +#include "DataSet_string.h" #include "StringRoutines.h" void Exec_DataSetCmd::Help() const { @@ -789,7 +790,20 @@ Exec::RetType Exec_DataSetCmd::ChangeModeType(CpptrajState const& State, ArgList Exec::RetType Exec_DataSetCmd::InvertSets(CpptrajState& State, ArgList& argIn) { DataSetList& DSL = State.DSL(); // Get keywords - std::string dsname = argIn.GetStringKey("name"); + DataSet* inputNames = 0; + std::string dsname = argIn.GetStringKey("legendset"); + if (!dsname.empty()) { + inputNames = DSL.GetDataSet( dsname ); + if (inputNames == 0) { + mprinterr("Error: Name set '%s' not found.\n", dsname.c_str()); + return CpptrajState::ERR; + } + if (inputNames->Type() != DataSet::STRING) { + mprinterr("Error: Set '%s' does not contain strings.\n", inputNames->legend()); + return CpptrajState::ERR; + } + } + dsname = argIn.GetStringKey("name"); if (dsname.empty()) { mprinterr("Error: 'invert' requires that 'name ' be specified.\n"); return CpptrajState::ERR; @@ -824,18 +838,28 @@ Exec::RetType Exec_DataSetCmd::InvertSets(CpptrajState& State, ArgList& argIn) { mprinterr("Error: No sets selected.\n"); return CpptrajState::ERR; } + if (inputNames != 0 && inputNames->Size() != input_sets.front()->Size()) { + mprinterr("Error: Name set '%s' size (%zu) differs from # data points (%zu).\n", + inputNames->legend(), inputNames->Size(), input_sets.front()->Size()); + return CpptrajState::ERR; + } // Need an output data set for each point in input sets std::vector output_sets; int column = 1; for (int idx = 0; idx != (int)input_sets[0]->Size(); idx++, column++) { - DataSet* ds = DSL.AddSet(outtype, MetaData(dsname, column)); + DataSet* ds = 0; + ds = DSL.AddSet(outtype, MetaData(dsname, column)); if (ds == 0) return CpptrajState::ERR; + if (inputNames != 0) + ds->SetLegend( (*((DataSet_string*)inputNames))[idx] ); output_sets.push_back( ds ); if (outfile != 0) outfile->AddDataSet( ds ); } // Create a data set containing names of each input data set DataSet* nameset = DSL.AddSet(DataSet::STRING, MetaData(dsname, column)); if (nameset == 0) return CpptrajState::ERR; + if (inputNames != 0) + nameset->SetLegend("Names"); if (outfile != 0) outfile->AddDataSet( nameset ); // Populate output data sets for (int jdx = 0; jdx != (int)input_sets.size(); jdx++) From 434284f50c53b4f051a6e57ec7ed5cc307f7e716 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Fri, 24 Aug 2018 08:42:10 -0400 Subject: [PATCH 13/20] DRR - Cpptraj: When writing in grace format, if a single string set is present use it as data set labels. --- src/DataIO_Grace.cpp | 24 +++++++++++++++++++++++- src/cpptrajdepend | 6 +++--- 2 files changed, 26 insertions(+), 4 deletions(-) diff --git a/src/DataIO_Grace.cpp b/src/DataIO_Grace.cpp index c724af0fc1..db4b3aa4b8 100644 --- a/src/DataIO_Grace.cpp +++ b/src/DataIO_Grace.cpp @@ -4,6 +4,7 @@ #include "CpptrajStdio.h" #include "BufferedLine.h" #include "DataSet_double.h" +#include "DataSet_string.h" // TODO: Set dimension labels // Dont assume anything about set ordering @@ -130,10 +131,29 @@ int DataIO_Grace::WriteDataNormal(CpptrajFile& file, DataSetList const& Sets) { file.Printf("@with g0\n@ xaxis label \"%s\"\n@ yaxis label \"%s\"\n" "@ legend 0.2, 0.995\n@ legend char size 0.60\n", Sets[0]->Dim(0).Label().c_str(), ""); + // Determine if we have a single STRING DataSet. Assume these contain labels. + DataSet_string* labelSet = 0; + for (DataSetList::const_iterator set = Sets.begin(); set != Sets.end(); ++set) + { + if (labelSet == 0) { + if ( (*set)->Type() == DataSet::STRING ) + labelSet = (DataSet_string*)(*set); + } else { + if ( (*set)->Type() == DataSet::STRING ) { + labelSet = 0; + break; + } + } + } + if (labelSet != 0) + mprintf("\tUsing string set '%s' for data point labels\n", labelSet->legend()); // Loop over DataSets unsigned int setnum = 0; DataSet::SizeArray frame(1); - for (DataSetList::const_iterator set = Sets.begin(); set != Sets.end(); ++set, ++setnum) { + for (DataSetList::const_iterator set = Sets.begin(); set != Sets.end(); ++set, ++setnum) + { + // Skip the label set if defined. + if (*set == labelSet) continue; size_t maxFrames = (*set)->Size(); // Set information file.Printf("@ s%u legend \"%s\"\n@target G0.S%u\n@type xy\n", @@ -148,6 +168,8 @@ int DataIO_Grace::WriteDataNormal(CpptrajFile& file, DataSetList const& Sets) { for (frame[0] = 0; frame[0] < maxFrames; frame[0]++) { file.Printf(xfmt.fmt(), (*set)->Coord(0, frame[0])); (*set)->WriteBuffer(file, frame); + if (labelSet != 0) + file.Printf(" \"%s\"", (*labelSet)[frame[0]].c_str()); file.Printf("\n"); } } diff --git a/src/cpptrajdepend b/src/cpptrajdepend index 9808883637..ebb6c70eaf 100644 --- a/src/cpptrajdepend +++ b/src/cpptrajdepend @@ -117,7 +117,7 @@ Analysis_RmsAvgCorr.o : Analysis_RmsAvgCorr.cpp ActionState.h Analysis.h Analysi Analysis_Rotdif.o : Analysis_Rotdif.cpp ActionState.h Analysis.h AnalysisState.h Analysis_Rotdif.h ArgList.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h ComplexArray.h Constants.h CoordinateInfo.h Corr.h CpptrajFile.h CpptrajStdio.h CurveFit.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_Mat3x3.h DataSet_Mesh.h DataSet_Vector.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h ProgressBar.h PubFFT.h Random.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SimplexMin.h Spline.h StringRoutines.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h Analysis_RunningAvg.o : Analysis_RunningAvg.cpp ActionState.h Analysis.h AnalysisState.h Analysis_RunningAvg.h ArgList.h Array1D.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_Mesh.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h Spline.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h Analysis_Spline.o : Analysis_Spline.cpp ActionState.h Analysis.h AnalysisState.h Analysis_Spline.h ArgList.h Array1D.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h Constants.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_Mesh.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h Spline.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h -Analysis_State.o : Analysis_State.cpp ActionState.h Analysis.h AnalysisState.h Analysis_State.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_double.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h +Analysis_State.o : Analysis_State.cpp ActionState.h Analysis.h AnalysisState.h Analysis_State.h ArgList.h Array1D.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_double.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h Analysis_Statistics.o : Analysis_Statistics.cpp ActionState.h Analysis.h AnalysisState.h Analysis_Statistics.h ArgList.h Array1D.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_float.h DataSet_integer.h DataSet_string.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h Analysis_TI.o : Analysis_TI.cpp ActionState.h Analysis.h AnalysisState.h Analysis_TI.h ArgList.h Array1D.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_Mesh.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Random.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h Spline.h StringRoutines.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h Analysis_Timecorr.o : Analysis_Timecorr.cpp ActionState.h Analysis.h AnalysisState.h Analysis_Timecorr.h ArgList.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h ComplexArray.h CoordinateInfo.h Corr.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_Vector.h DataSet_double.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h PubFFT.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h @@ -177,7 +177,7 @@ DataIO_Cmatrix.o : DataIO_Cmatrix.cpp ArgList.h ArrayIterator.h AssociatedData.h DataIO_Cpout.o : DataIO_Cpout.cpp ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h BufferedLine.h CharMask.h CoordinateInfo.h Cph.h CpptrajFile.h CpptrajStdio.h DataIO.h DataIO_Cpout.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_PHREMD.h DataSet_PHREMD_Explicit.h DataSet_PHREMD_Implicit.h DataSet_pH.h Dimension.h FileIO.h FileName.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h DataIO_Evecs.o : DataIO_Evecs.cpp ArgList.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h BufferedFrame.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataIO.h DataIO_Evecs.h DataSet.h DataSetList.h DataSet_2D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_MatrixDbl.h DataSet_Modes.h Dimension.h FileIO.h FileName.h Frame.h MaskToken.h Matrix.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h DataIO_Gnuplot.o : DataIO_Gnuplot.cpp ArgList.h Array1D.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h BufferedLine.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataIO.h DataIO_Gnuplot.h DataSet.h DataSetList.h DataSet_1D.h DataSet_2D.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h FileIO.h FileName.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h -DataIO_Grace.o : DataIO_Grace.cpp ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h BufferedLine.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataIO.h DataIO_Grace.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_double.h Dimension.h FileIO.h FileName.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h +DataIO_Grace.o : DataIO_Grace.cpp ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h BufferedLine.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataIO.h DataIO_Grace.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_double.h DataSet_string.h Dimension.h FileIO.h FileName.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h DataIO_Mdout.o : DataIO_Mdout.cpp ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h BufferedLine.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataIO.h DataIO_Mdout.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_double.h Dimension.h FileIO.h FileName.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h StringRoutines.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h DataIO_NC_Cmatrix.o : DataIO_NC_Cmatrix.cpp ArgList.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMap.h AtomMask.h BaseIOtype.h Box.h CharMask.h ClusterDist.h ClusterSieve.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataIO.h DataIO_NC_Cmatrix.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Cmatrix.h DataSet_Cmatrix_MEM.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h FileIO.h FileName.h Frame.h Hungarian.h MapAtom.h MaskToken.h Matrix.h Matrix_3x3.h MetaData.h Molecule.h NC_Cmatrix.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h SymmetricRmsdCalc.h TextFormat.h Timer.h Topology.h Vec3.h DataIO_OpenDx.o : DataIO_OpenDx.cpp ArgList.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h BufferedLine.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataIO.h DataIO_OpenDx.h DataSet.h DataSetList.h DataSet_3D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_GridFlt.h Dimension.h FileIO.h FileName.h Frame.h Grid.h GridBin.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h ProgressBar.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h @@ -245,7 +245,7 @@ Exec_CrdOut.o : Exec_CrdOut.cpp Action.h ActionFrameCounter.h ActionList.h Actio Exec_CreateSet.o : Exec_CreateSet.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_CreateSet.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h RPNcalc.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h Vec3.h Exec_DataFile.o : Exec_DataFile.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_DataFile.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h StringRoutines.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h Vec3.h Exec_DataFilter.o : Exec_DataFilter.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Action_FilterByData.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h Array1D.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_DataFilter.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h ProgressBar.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h Vec3.h -Exec_DataSetCmd.o : Exec_DataSetCmd.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h ComplexArray.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_2D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_MatrixDbl.h DataSet_Vector.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_DataSetCmd.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h StringRoutines.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h Vec3.h +Exec_DataSetCmd.o : Exec_DataSetCmd.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h ComplexArray.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_2D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_MatrixDbl.h DataSet_Vector.h DataSet_string.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_DataSetCmd.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h StringRoutines.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h Vec3.h Exec_GenerateAmberRst.o : Exec_GenerateAmberRst.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h Constants.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h DistRoutines.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_GenerateAmberRst.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TorsionRoutines.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h Vec3.h Exec_Help.o : Exec_Help.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h Cmd.h CmdList.h Command.h Control.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_Help.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h VariableArray.h Vec3.h Exec_LoadCrd.o : Exec_LoadCrd.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_LoadCrd.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h Trajin_Single.h TrajoutList.h Trajout_Single.h Vec3.h From ecc8daeff09c651d71719e97fb7bd47d077271a5 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Fri, 24 Aug 2018 10:14:53 -0400 Subject: [PATCH 14/20] DRR - Cpptraj: Condense the state superseded warning --- src/Analysis_State.cpp | 14 +++++++++++--- 1 file changed, 11 insertions(+), 3 deletions(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index 3ef290994e..0db275fbfe 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -203,6 +203,7 @@ Analysis::RetType Analysis_State::Analyze() { size_t last_State_Start = 0; size_t final_frame = nframes - 1; std::vector stateFrames(numStates, 0); + std::vector overwriteCount(NameMap_.size(), 0); for (size_t frm = 0; frm != nframes; frm++) { // Determine which state we are in. int state_num = -1; @@ -210,9 +211,12 @@ Analysis::RetType Analysis_State::Analyze() { { if (state->InState( frm )) { if (state_num != -1) { - mprintf("Warning: Frame %zu already defined as state '%s', also matches state '%s'.\n", - frm, States_[state_num].id(), state->id()); - mprintf("Warning: Overriding previous state.\n"); + if (debug_ > 0) { + mprintf("Warning: Frame %zu already defined as state '%s', also matches state '%s'.\n", + frm, States_[state_num].id(), state->id()); + mprintf("Warning: Overriding previous state.\n"); + } + overwriteCount[state_num]++; state_num = state->Num(); } else state_num = state->Num(); @@ -260,6 +264,10 @@ Analysis::RetType Analysis_State::Analyze() { last_state = state_num; } } // END loop over frames + for (unsigned int idx = 0; idx != overwriteCount.size(); idx++) + if (overwriteCount[idx] > 0) + mprintf("\tState %s (%u) was superseded by a subsequent state %u times\n", + stateName(idx), idx, overwriteCount[idx]); // DEBUG: Print single state info. if (debug_ > 0) { From 10dc48bf766f1ae7a209f32b173dc9e0ad50e9cd Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Fri, 24 Aug 2018 11:06:45 -0400 Subject: [PATCH 15/20] DRR - Cpptraj: Add std data io groupbyname option --- src/DataIO_Std.cpp | 57 +++++++++++++++++++++++++++++++++++++--------- src/DataIO_Std.h | 1 + src/DataSetList.h | 4 ++-- 3 files changed, 49 insertions(+), 13 deletions(-) diff --git a/src/DataIO_Std.cpp b/src/DataIO_Std.cpp index 3cd6e39a72..3d7cf14c2a 100644 --- a/src/DataIO_Std.cpp +++ b/src/DataIO_Std.cpp @@ -24,6 +24,7 @@ DataIO_Std::DataIO_Std() : prec_(UNSPEC), indexcol_(-1), isInverted_(false), + groupByName_(false), hasXcolumn_(true), writeHeader_(true), square2d_(false), @@ -600,13 +601,14 @@ int DataIO_Std::Read_Mat3x3(std::string const& fname, // ----------------------------------------------------------------------------- void DataIO_Std::WriteHelp() { - mprintf("\tnoheader : Do not print header line.\n" - "\tinvert : Flip X/Y axes (1D).\n" - "\tnoxcol : Do not print X (index) column (1D).\n" - "\tsquare2d : Write 2D data sets in matrix-like format.\n" - "\tnosquare2d : Write 2D data sets as ' '.\n" - "\tnosparse : Write all 3D grid voxels (default).\n" - "\tsparse : Only write 3D grid voxels with value > cutoff (default 0).\n" + mprintf("\tnoheader : Do not print header line.\n" + "\tinvert : Flip X/Y axes (1D).\n" + "\tgroupbyname : (1D) group data sets by name,\n" + "\tnoxcol : Do not print X (index) column (1D).\n" + "\tsquare2d : Write 2D data sets in matrix-like format.\n" + "\tnosquare2d : Write 2D data sets as ' '.\n" + "\tnosparse : Write all 3D grid voxels (default).\n" + "\tsparse : Only write 3D grid voxels with value > cutoff (default 0).\n" "\t\tcut : Cutoff for 'sparse'; default 0.\n"); } @@ -614,6 +616,8 @@ void DataIO_Std::WriteHelp() { int DataIO_Std::processWriteArgs(ArgList &argIn) { if (!isInverted_ && argIn.hasKey("invert")) isInverted_ = true; + if (!groupByName_ && argIn.hasKey("groupbyname")) + groupByName_ = true; if (hasXcolumn_ && argIn.hasKey("noxcol")) hasXcolumn_ = false; if (writeHeader_ && argIn.hasKey("noheader")) @@ -677,10 +681,41 @@ int DataIO_Std::WriteData(FileName const& fname, DataSetList const& SetList) // Special case of 2D - may have sieved frames. err = WriteCmatrix(file, SetList); } else if (SetList[0]->Ndim() == 1) { - if (isInverted_) - err = WriteDataInverted(file, SetList); - else - err = WriteDataNormal(file, SetList); + if (groupByName_) { + DataSetList tmpdsl; + std::vector setIsWritten(SetList.size(), false); + unsigned int startIdx = 0; + unsigned int nWritten = 0; + while (nWritten < SetList.size()) { + std::string currentName = SetList[startIdx]->Meta().Name(); + int firstNonMatch = -1; + for (unsigned int idx = startIdx; idx != SetList.size(); idx++) + { + if (!setIsWritten[idx]) + { + if (currentName == SetList[idx]->Meta().Name()) + { + tmpdsl.AddCopyOfSet( SetList[idx] ); + setIsWritten[idx] = true; + nWritten++; + } else if (firstNonMatch == -1) + firstNonMatch = (int)idx; + } + } + if (firstNonMatch > -1) + startIdx = (unsigned int)firstNonMatch; + if (isInverted_) + err += WriteDataInverted(file, tmpdsl); + else + err += WriteDataNormal(file, tmpdsl); + tmpdsl.ClearAll(); + } + } else { + if (isInverted_) + err = WriteDataInverted(file, SetList); + else + err = WriteDataNormal(file, SetList); + } } else if (SetList[0]->Ndim() == 2) err = WriteData2D(file, SetList); else if (SetList[0]->Ndim() == 3) diff --git a/src/DataIO_Std.h b/src/DataIO_Std.h index 8c4e4165d7..b6653b3acb 100644 --- a/src/DataIO_Std.h +++ b/src/DataIO_Std.h @@ -35,6 +35,7 @@ class DataIO_Std : public DataIO { precType prec_; ///< 3d reads, data set precision int indexcol_; ///< Read: column containing index (X) values bool isInverted_; ///< For 1D writes invert X/Y. + bool groupByName_; ///< For 1D writes, group sets with same name together bool hasXcolumn_; bool writeHeader_; bool square2d_; diff --git a/src/DataSetList.h b/src/DataSetList.h index 8c81d05efe..a2997fab23 100644 --- a/src/DataSetList.h +++ b/src/DataSetList.h @@ -38,6 +38,8 @@ class DataSetList { const_iterator end() const { return DataList_.end(); } /// Clear all non-Topology and non-Reference DataSets void Clear(); + /// Clear entire DataSetList + void ClearAll(); /// Sort sets in the DataSetList void Sort(); /// True if no DataSets in list. @@ -144,8 +146,6 @@ class DataSetList { void Timing() const; # endif private: - /// Clear entire DataSetList - void ClearAll(); /// Search for and remove specified data set if found, optionally free memory. DataSet* EraseSet( DataSet*, bool ); /// Warn if DataSet not found but may be pending. From e77bab422ba4b3c6ad197bfbacfdc458471ca0df Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Fri, 24 Aug 2018 12:52:36 -0400 Subject: [PATCH 16/20] DRR - Cpptraj: Update help. --- src/Analysis_State.cpp | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index 0db275fbfe..4e77e7c243 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -18,10 +18,12 @@ Analysis_State::Analysis_State() : {} void Analysis_State::Help() const { - mprintf("\t{state ,,,} [out ] [name ]\n" - "\t[curveout ] [stateout ] [transout ]\n" + mprintf("\t{state ,,,[,,,]} ...\n" + "\t[out ] [name ]\n" + "\t[curveout ] [stateout ]\n" + "\t[transout ] [countout ]\n" " Data for the specified data set(s) that matches the given criteria\n" - " will be assigned state <#>.\n"); + " will be assigned state .\n"); } // Analysis_State::Setup() From 4854363118b48d648c525c193c05075b19081c3d Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Mon, 27 Aug 2018 09:28:25 -0400 Subject: [PATCH 17/20] DRR - Cpptraj: Fix up debug output --- src/Analysis_State.cpp | 13 +++++-------- 1 file changed, 5 insertions(+), 8 deletions(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index 4e77e7c243..89f57d9f70 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -274,14 +274,11 @@ Analysis::RetType Analysis_State::Analyze() { // DEBUG: Print single state info. if (debug_ > 0) { mprintf(" States:\n"); - mprintf("\t%i: %s max= %i sum= %i n= %i Avg= %g\n", -1, "Undefined", - Status.front().Max(), Status.front().Sum(), - Status.front().Nlifetimes(), Status.front().Avg()); - StateArray::const_iterator state = States_.begin(); - for (std::vector::const_iterator s = Status.begin() + 1; - s != Status.end(); ++s, ++state) - mprintf("\t%u: %s max= %i sum= %i n= %i Avg= %g\n", state - States_.begin(), - state->id(), s->Max(), s->Sum(), s->Nlifetimes(), s->Avg()); + for (int stateid = -1; stateid < (int)NameMap_.size(); stateid++) + mprintf("\t%i: %s max= %i sum= %i n= %i Avg= %g\n", + stateid, stateName(stateid), Status[stateid+1].Max(), + Status[stateid+1].Sum(), Status[stateid+1].Nlifetimes(), + Status[stateid+1].Avg()); // DEBUG: Print transitions. mprintf(" Transitions:\n"); for (TransMapType::const_iterator trans = TransitionMap_.begin(); From 162ba066831d6663f2d626f69769a6f40fcc3146 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Mon, 27 Aug 2018 13:22:24 -0400 Subject: [PATCH 18/20] DRR - Cpptraj: Do not print definitions index to avoid confusion with state index --- src/Analysis_State.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/Analysis_State.cpp b/src/Analysis_State.cpp index 89f57d9f70..a9c2b6208d 100644 --- a/src/Analysis_State.cpp +++ b/src/Analysis_State.cpp @@ -140,7 +140,7 @@ Analysis::RetType Analysis_State::Setup(ArgList& analyzeArgs, AnalysisSetup& set mprintf(" STATE: The following state definitions have been set up:\n"); for (StateArray::const_iterator state = States_.begin(); state != States_.end(); ++state) { - mprintf("\t%u: ", state - States_.begin()); + mprintf("\t ", state - States_.begin()); state->PrintState(); } mprintf("\tState data set: %s\n", state_data_->legend()); From a599bf25f5ff632bc1b374f4d964eff46a25db11 Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Mon, 27 Aug 2018 13:32:51 -0400 Subject: [PATCH 19/20] DRR - Cpptraj: Add help and more feedback for dataset invert. --- src/Exec_DataSetCmd.cpp | 15 ++++++++++++++- src/Exec_DataSetCmd.h | 1 + 2 files changed, 15 insertions(+), 1 deletion(-) diff --git a/src/Exec_DataSetCmd.cpp b/src/Exec_DataSetCmd.cpp index 7bd1999de7..1169b72ed6 100644 --- a/src/Exec_DataSetCmd.cpp +++ b/src/Exec_DataSetCmd.cpp @@ -9,7 +9,7 @@ void Exec_DataSetCmd::Help() const { mprintf("\t{legend|makexy|vectorcoord|cat|make2d|droppoints|keeppoints|remove|\n" - "\t dim|outformat|mode|type} \n"); + "\t dim|outformat|invert|mode|type} \n"); mprintf(" legend \n" " Set the legend for a single data set.\n"); mprintf(" makexy [name ]\n" @@ -35,6 +35,7 @@ void Exec_DataSetCmd::Help() const { " double - \"Normal\" output, e.g. 0.4032\n" " scientific - Scientific \"E\" notation output, e.g. 4.032E-1\n" " general - Use 'double' or 'scientific', whichever is shortest.\n"); + Help_InvertSets(); mprintf(" [mode ] [type ] [ ...]\n"); mprintf(" : "); for (int i = 0; i != (int)MetaData::UNKNOWN_MODE; i++) @@ -786,6 +787,12 @@ Exec::RetType Exec_DataSetCmd::ChangeModeType(CpptrajState const& State, ArgList return CpptrajState::OK; } +void Exec_DataSetCmd::Help_InvertSets() { + mprintf(" invert ... name [legendset ]\n" + " Given M input sets of length N, create N new sets of length M by\n" + " inverting the input sets.\n"); +} + /** Syntax: dataset invert ... name */ Exec::RetType Exec_DataSetCmd::InvertSets(CpptrajState& State, ArgList& argIn) { DataSetList& DSL = State.DSL(); @@ -802,13 +809,17 @@ Exec::RetType Exec_DataSetCmd::InvertSets(CpptrajState& State, ArgList& argIn) { mprinterr("Error: Set '%s' does not contain strings.\n", inputNames->legend()); return CpptrajState::ERR; } + mprintf("\tUsing names from set '%s' as legends for inverted sets.\n", inputNames->legend()); } dsname = argIn.GetStringKey("name"); if (dsname.empty()) { mprinterr("Error: 'invert' requires that 'name ' be specified.\n"); return CpptrajState::ERR; } + mprintf("\tNew sets will be named '%s'\n", dsname.c_str()); DataFile* outfile = State.DFL().AddDataFile( argIn.GetStringKey("out"), argIn ); + if (outfile != 0) + mprintf("\tNew sets will be output to '%s'\n", outfile->DataFilename().full()); // TODO determine type some other way DataSet::DataType outtype = DataSet::DOUBLE; // Get input DataSets @@ -843,6 +854,8 @@ Exec::RetType Exec_DataSetCmd::InvertSets(CpptrajState& State, ArgList& argIn) { inputNames->legend(), inputNames->Size(), input_sets.front()->Size()); return CpptrajState::ERR; } + mprintf("\t%zu input sets; creating %zu output sets.\n", + input_sets.size(), input_sets.front()->Size()); // Need an output data set for each point in input sets std::vector output_sets; int column = 1; diff --git a/src/Exec_DataSetCmd.h b/src/Exec_DataSetCmd.h index 291ee29638..b5d74b7959 100644 --- a/src/Exec_DataSetCmd.h +++ b/src/Exec_DataSetCmd.h @@ -30,6 +30,7 @@ class Exec_DataSetCmd : public Exec { static void Help_ChangeDim(); RetType ChangeDim(CpptrajState const&, ArgList&); RetType ChangeModeType(CpptrajState const&, ArgList&); + static void Help_InvertSets(); RetType InvertSets(CpptrajState&, ArgList&); }; #endif From ac9397e75cef9e76686545d1e38a2a29483b25ab Mon Sep 17 00:00:00 2001 From: "Daniel R. Roe" Date: Mon, 27 Aug 2018 13:37:27 -0400 Subject: [PATCH 20/20] DRR - Cpptraj: dataset invert test --- test/Test_DataSetCmd/Invert.dat.save | 3 +++ test/Test_DataSetCmd/RunTest.sh | 6 +++++- 2 files changed, 8 insertions(+), 1 deletion(-) create mode 100644 test/Test_DataSetCmd/Invert.dat.save diff --git a/test/Test_DataSetCmd/Invert.dat.save b/test/Test_DataSetCmd/Invert.dat.save new file mode 100644 index 0000000000..9aa4ba360c --- /dev/null +++ b/test/Test_DataSetCmd/Invert.dat.save @@ -0,0 +1,3 @@ +#Frame Invert:1 Invert:2 Invert:3 + 1 1.0000 2.0000 D1:1 + 2 3.0000 4.0000 D2:1 diff --git a/test/Test_DataSetCmd/RunTest.sh b/test/Test_DataSetCmd/RunTest.sh index 9ddec55ccd..2caad41eab 100755 --- a/test/Test_DataSetCmd/RunTest.sh +++ b/test/Test_DataSetCmd/RunTest.sh @@ -2,7 +2,8 @@ . ../MasterTest.sh -CleanFiles data.in avg.dat matrix.dat matrix2.dat VXYZ.dat Keep?.dat Drop?.dat D4.dat +CleanFiles data.in avg.dat matrix.dat matrix2.dat VXYZ.dat Keep?.dat \ + Drop?.dat D4.dat Invert.dat INPUT='-i data.in' cat > data.in <