More comments for trajectoryanalysis.
authorTeemu Murtola <teemu.murtola@gmail.com>
Mon, 2 Apr 2012 18:59:29 +0000 (21:59 +0300)
committerTeemu Murtola <teemu.murtola@gmail.com>
Tue, 10 Apr 2012 04:25:40 +0000 (07:25 +0300)
- Added/updated/made more uniform Doxygen comments for general classes
  in the trajectoryanalysis subdirectory.
- Removed an unnecessary method to avoid copying documentation for it.
- Removed some commented-out code.

Part of issue #640 (child of #867).

Change-Id: I2fb4a03a2054ba0e63581b3815ad6f7176ce57a2

src/gromacs/trajectoryanalysis.h
src/gromacs/trajectoryanalysis/analysismodule-impl.h
src/gromacs/trajectoryanalysis/analysismodule.cpp
src/gromacs/trajectoryanalysis/analysismodule.h
src/gromacs/trajectoryanalysis/analysissettings.h
src/gromacs/trajectoryanalysis/cmdlinerunner.cpp
src/gromacs/trajectoryanalysis/cmdlinerunner.h
src/gromacs/trajectoryanalysis/modules.h

index 386c83048d01f25b5a1423e943b0c6dc32b1f693..26d3d083c6938358ae418154a294eba48ac413e5 100644 (file)
@@ -58,6 +58,9 @@
  * Internally, the module also defines a set of trajectory analysis modules that
  * can be instantiated using createTrajectoryAnalysisModule().
  *
+ * For an example of how to implement an analysis tool using the framework, see
+ * \ref template.cpp.
+ *
  * \author Teemu Murtola <teemu.murtola@cbr.su.se>
  */
 /*! \file
index 0d0f4973f6f89aee155d066de51d02863e2007e1..47539b5e8c86565c30b726efa64d5be4d5027a30 100644 (file)
@@ -52,31 +52,52 @@ class AbstractAnalysisData;
 class AnalysisData;
 class AnalysisDataHandle;
 
+/*! \internal \brief
+ * Private implementation class for TrajectoryAnalysisModuleData.
+ *
+ * \ingroup module_trajectoryanalysis
+ */
 class TrajectoryAnalysisModuleData::Impl
 {
     public:
+        //! Container that associates a data handle to its AnalysisData object.
         typedef std::map<const AnalysisData *, AnalysisDataHandle>
                 HandleContainer;
 
+        //! \copydoc TrajectoryAnalysisModuleData::TrajectoryAnalysisModuleData()
         Impl(TrajectoryAnalysisModule *module,
              const AnalysisDataParallelOptions &opt,
              const SelectionCollection &selections);
         ~Impl();
 
-        void finishHandles();
-
+        //! Keeps a data handle for each AnalysisData object.
         HandleContainer         _handles;
+        //! Stores thread-local selections.
         const SelectionCollection &_selections;
 };
 
+/*! \internal \brief
+ * Private implementation class for TrajectoryAnalysisModule.
+ *
+ * \ingroup module_trajectoryanalysis
+ */
 class TrajectoryAnalysisModule::Impl
 {
     public:
+        //! Container that associates a data set with its name.
         typedef std::map<std::string, AbstractAnalysisData *> DatasetContainer;
+        //! Container that associates a AnalysisData object with its name.
         typedef std::map<std::string, AnalysisData *> AnalysisDatasetContainer;
 
+        //! List of registered data set names.
         std::vector<std::string>        _datasetNames;
+        /*! \brief
+         * Keeps all registered data sets.
+         *
+         * This container also includes datasets from \a _analysisDatasets.
+         */
         DatasetContainer                _datasets;
+        //! Keeps registered AnalysisData objects.
         AnalysisDatasetContainer        _analysisDatasets;
 };
 
@@ -85,6 +106,8 @@ class TrajectoryAnalysisModule::Impl
  *
  * Most simple tools should only require data handles and selections to be
  * thread-local, so this class implements just that.
+ *
+ * \ingroup module_trajectoryanalysis
  */
 class TrajectoryAnalysisModuleDataBasic : public TrajectoryAnalysisModuleData
 {
index 3b9fe8f26c42c9a966c02c1b0879a2d70314f994..d2099899b399b4babf66d59dcc17eaa350cddaec 100644 (file)
@@ -72,18 +72,6 @@ TrajectoryAnalysisModuleData::Impl::~Impl()
 }
 
 
-void TrajectoryAnalysisModuleData::Impl::finishHandles()
-{
-    // FIXME: Call finishData() for all handles even if one throws
-    HandleContainer::iterator i;
-    for (i = _handles.begin(); i != _handles.end(); ++i)
-    {
-        i->second.finishData();
-    }
-    _handles.clear();
-}
-
-
 /********************************************************************
  * TrajectoryAnalysisModuleData
  */
@@ -104,7 +92,13 @@ TrajectoryAnalysisModuleData::~TrajectoryAnalysisModuleData()
 
 void TrajectoryAnalysisModuleData::finishDataHandles()
 {
-    _impl->finishHandles();
+    // FIXME: Call finishData() for all handles even if one throws
+    Impl::HandleContainer::iterator i;
+    for (i = _impl->_handles.begin(); i != _impl->_handles.end(); ++i)
+    {
+        i->second.finishData();
+    }
+    _impl->_handles.clear();
 }
 
 
@@ -238,6 +232,7 @@ AbstractAnalysisData &TrajectoryAnalysisModule::datasetFromName(const char *name
 void TrajectoryAnalysisModule::registerBasicDataset(AbstractAnalysisData *data,
                                                     const char *name)
 {
+    // TODO: Strong exception safety should be possible to implement.
     GMX_RELEASE_ASSERT(_impl->_datasets.find(name) == _impl->_datasets.end(),
                        "Duplicate data set name registered");
     _impl->_datasets[name] = data;
@@ -248,6 +243,7 @@ void TrajectoryAnalysisModule::registerBasicDataset(AbstractAnalysisData *data,
 void TrajectoryAnalysisModule::registerAnalysisDataset(AnalysisData *data,
                                                        const char *name)
 {
+    // TODO: Strong exception safety should be possible to implement.
     registerBasicDataset(data, name);
     _impl->_analysisDatasets[name] = data;
 }
index 32bc2c6452962be0caf610d9194daa9611437d7c..8226be1301c4a978822b07c14674d5afe9706c68 100644 (file)
@@ -67,7 +67,8 @@ class TrajectoryAnalysisSettings;
  *
  * Thread-local storage of data handles and selections is implemented in this
  * class; TrajectoryAnalysisModule instances can access the thread-local values
- * using dataHandle() and parallelSelection().
+ * in their TrajectoryAnalysisModule::analyzeFrame() method using dataHandle()
+ * and parallelSelection().
  *
  * \see TrajectoryAnalysisModule::startFrames()
  * \see TrajectoryAnalysisModule::analyzeFrame()
@@ -84,7 +85,13 @@ class TrajectoryAnalysisModuleData
         /*! \brief
          * Performs any finishing actions after all frames have been processed.
          *
-         * This function is called immediately before the destructor.
+         * \throws  unspecified Implementation may throw exceptions to indicate
+         *      errors.
+         *
+         * This function is called immediately before the destructor, after
+         * TrajectoryAnalysisModule::finishFrames().
+         * Derived classes should implement any final operations that need to
+         * be done after successful analysis.
          * All implementations should call finishDataHandles().
          */
         virtual void finish() = 0;
@@ -92,16 +99,37 @@ class TrajectoryAnalysisModuleData
         /*! \brief
          * Returns a data handle for a given dataset.
          *
-         * Allowed data sets are those that have been registered with
+         * \param[in] data  Analysis data object.
+         * \returns   Data handle for \p data stored in this thread-local data.
+         *
+         * \p data should have previously been registered with
          * TrajectoryAnalysisModule::registerAnalysisDataset().
+         *
+         * Does not throw.
          */
         AnalysisDataHandle dataHandle(const AnalysisData &data);
         /*! \brief
          * Returns a selection that corresponds to the given selection.
+         *
+         * \param[in] selection Global selection object.
+         * \returns   Selection object corresponding to this thread-local data.
+         *
+         * \p selection is the selection object that was obtained from
+         * SelectionOption.  The return value is the corresponding selection
+         * in the selection collection with which this data object was
+         * constructed with.
+         *
+         * Does not throw.
          */
         Selection parallelSelection(const Selection &selection);
         /*! \brief
          * Returns a set of selection that corresponds to the given selections.
+         *
+         * \throws std::bad_alloc if out of memory.
+         *
+         * Works as parallelSelection(), but for a list of selections at once.
+         *
+         * \see parallelSelection()
          */
         SelectionList parallelSelections(const SelectionList &selections);
 
@@ -112,6 +140,9 @@ class TrajectoryAnalysisModuleData
          * \param[in] module     Analysis module to use for data objects.
          * \param[in] opt        Data parallelization options.
          * \param[in] selections Thread-local selection collection.
+         * \throws  std::bad_alloc if out of memory.
+         * \throws  unspecified Can throw any exception thrown by
+         *      AnalysisData::startData().
          *
          * Calls AnalysisData::startData() on all data objects registered with
          * TrajectoryAnalysisModule::registerAnalysisDataset() in \p module.
@@ -124,6 +155,9 @@ class TrajectoryAnalysisModuleData
         /*! \brief
          * Calls finishData() on all data handles.
          *
+         * \throws  unspecified Can throw any exception thrown by
+         *      AnalysisDataHandle::finishData().
+         *
          * This function should be called from the implementation of finish()
          * in all subclasses.
          */
@@ -140,10 +174,10 @@ typedef gmx_unique_ptr<TrajectoryAnalysisModuleData>::type
         TrajectoryAnalysisModuleDataPointer;
 
 /*! \brief
- * Base class for trajectory analysis methods.
+ * Base class for trajectory analysis modules.
  *
  * Trajectory analysis methods should derive from this class and override the
- * necessary virtual functions to implement initialization (initOptions(),
+ * necessary virtual methods to implement initialization (initOptions(),
  * initOptionsDone(), initAnalysis(), initAfterFirstFrame()), per-frame analysis
  * (analyzeFrame()), and final processing (finishFrames(), finishAnalysis(),
  * writeOutput()).
@@ -155,6 +189,32 @@ typedef gmx_unique_ptr<TrajectoryAnalysisModuleData>::type
  * to the other methods.  The default implementation of startFrames() can be
  * used if only data handles and selections need to be thread-local.
  *
+ * To get the full benefit from this class,
+ * \ref module_analysisdata "analysis data objects" and
+ * \ref module_selection "selections" should be used in the implementation.
+ * See the corresponding modules' documentation for details of how they work.
+ *
+ * Typical way of using AnalysisData in derived classes is to have the
+ * AnalysisData object as a member variable and register it using
+ * registerAnalysisDataset().  Analysis modules are initialized in
+ * initAnalysis() and the processing chain is initialized.  If any of the
+ * modules is required, e.g., for post-processing in finishAnalysis(), it can
+ * be stored in a member variable.  To add data to the data object in
+ * analyzeFrame(), a data handle is obtained using
+ * TrajectoryAnalysisModuleData::dataHandle().
+ *
+ * Typical way of using selections in derived classes is to have the required
+ * \ref Selection objects (or ::SelectionList objects) as member variables, and
+ * add the required selection options in initOptions().  These member variables
+ * can be accessed in initAnalysis() to get general information about the
+ * selections.  In analyzeFrame(), these selection objects should not be used
+ * directly, but instead TrajectoryAnalysisModuleData::parallelSelection()
+ * should be used to obtain a selection object that works correctly also for
+ * parallel analysis.
+ *
+ * Derived classes should use exceptions to indicate errors in the virtual
+ * methods.
+ *
  * \inpublicapi
  * \ingroup module_trajectoryanalysis
  */
@@ -166,9 +226,17 @@ class TrajectoryAnalysisModule
         /*! \brief
          * Initializes options understood by the module.
          *
-         * In addition to initializing the options, this function can also
-         * provide information about its requirements using the \p settings
-         * object; see TrajectoryAnalysisSettings for more details.
+         * \param[in,out] settings Settings to pass to and from the module.
+         * \returns  Reference to options that are accepted by the module.
+         *
+         * Typical implementation returns a reference to a member variable
+         * after first adding necessary options to that object.  Output values
+         * from options (including selections) should be stored in member
+         * variables.
+         *
+         * In addition to initializing the options, this method can also
+         * provide information about the module's requirements using the
+         * \p settings object; see TrajectoryAnalysisSettings for more details.
          *
          * If settings depend on the option values provided by the user, see
          * initOptionsDone().
@@ -177,9 +245,15 @@ class TrajectoryAnalysisModule
         /*! \brief
          * Called after all option values have been set.
          *
+         * \param[in,out] settings Settings to pass to and from the module.
+         *
+         * This method is called after option values have been assigned (but
+         * interactive selection input has not yet been performed).
+         *
          * If the module needs to change settings that affect topology loading
-         * or selection initialization based on option values, this function
-         * has to be overridden.
+         * (can be done using the \p settings object) or selection
+         * initialization (can be done using SelectionOptionInfo) based on
+         * option values, this method has to be overridden.
          *
          * The default implementation does nothing.
          */
@@ -187,6 +261,9 @@ class TrajectoryAnalysisModule
         /*! \brief
          * Initializes the analysis.
          *
+         * \param[in]    settings Settings to pass to and from the module.
+         * \param[in]    top      Topology information.
+         *
          * When this function is called, selections have been initialized based
          * on user input, and a topology has been loaded if provided by the
          * user.  For dynamic selections, the selections have been evaluated to
@@ -242,14 +319,15 @@ class TrajectoryAnalysisModule
          * \param[in]     pbc    Periodic boundary conditions for \p fr.
          * \param[in,out] pdata  Data structure for frame-local data.
          *
-         * This function is called once for each frame to be analyzed,
-         * and should analyze the positions provided in \p sel.
+         * This method is called once for each frame to be analyzed, and should
+         * analyze the positions provided in the selections.  Data handles and
+         * selections should be obtained from the \p pdata structure.
          *
-         * For threaded analysis, this function is called asynchronously in
-         * different threads to analyze different frames. The \p pdata
+         * For threaded analysis, this method is called asynchronously in
+         * different threads to analyze different frames.  The \p pdata
          * structure is one of the structures created with startFrames(),
          * but no assumptions should be made about which of these data
-         * structures is used. It is guaranteed that two instances of
+         * structures is used.  It is guaranteed that two instances of
          * analyzeFrame() are not running concurrently with the same \p pdata
          * data structure.
          * Any access to data structures not stored in \p pdata should be
@@ -262,14 +340,14 @@ class TrajectoryAnalysisModule
          *
          * \param[in]  pdata    Data structure for thread-local data.
          *
-         * This function is called once for each call of startFrames(),
-         * with the data structure returned by the corresponding startFrames().
+         * This method is called once for each call of startFrames(), with the
+         * data structure returned by the corresponding startFrames().
          * The \p pdata object should be destroyed by the caller after this
          * function has been called.
          *
          * You only need to override this method if you need custom
          * operations to combine data from the frame-local data structures
-         * to get the final result. In such cases, the data should be
+         * to get the final result.  In such cases, the data should be
          * aggregated in this function and stored in a member attribute.
          *
          * The default implementation does nothing.
@@ -281,8 +359,12 @@ class TrajectoryAnalysisModule
         /*! \brief
          * Postprocesses data after frames have been read.
          *
+         * \param[in]  nframes  Total number of frames processed.
+         *
          * This function is called after all finishFrames() calls have been
          * called.
+         * \p nframes will equal the number of calls to analyzeFrame() that
+         * have occurred.
          */
         virtual void finishAnalysis(int nframes) = 0;
         /*! \brief
@@ -297,18 +379,24 @@ class TrajectoryAnalysisModule
 
         /*! \brief
          * Returns the number of datasets provided by the module.
+         *
+         * Does not throw.
          */
         int datasetCount() const;
         /*! \brief
-         * Returns a vector with the names of the datasets.
+         * Returns a vector with the names of datasets provided by the module.
+         *
+         * Does not throw.
          */
         const std::vector<std::string> &datasetNames() const;
         /*! \brief
          * Returns a pointer to the data set \p index.
          *
          * \param[in] index  Data set to query for.
-         * \returns   A pointer to the data set, or NULL if \p index is not
-         *      valid.
+         * \returns   Reference to the requested data set.
+         * \throws    APIError if \p index is not valid.
+         *
+         * \p index should be >= 0 and < datasetCount().
          *
          * The return value is not const to allow callers to add modules to the
          * data sets. However, the AbstractAnalysisData interface does not
@@ -320,8 +408,10 @@ class TrajectoryAnalysisModule
          * Returns a pointer to the data set with name \p name
          *
          * \param[in] name  Data set to query for.
-         * \returns   A pointer to the data set, or NULL if \p name is not
-         *      recognized.
+         * \returns   Reference to the requested data set.
+         * \throws    APIError if \p name is not valid.
+         *
+         * \p name should be one of the names returned by datasetNames().
          *
          * The return value is not const to allow callers to add modules to the
          * data sets. However, the AbstractAnalysisData interface does not
@@ -331,15 +421,42 @@ class TrajectoryAnalysisModule
         AbstractAnalysisData &datasetFromName(const char *name) const;
 
     protected:
-        //! Initializes the dataset registration mechanism.
+        /*! \brief
+         * Initializes the dataset registration mechanism.
+         *
+         * \throws std::bad_alloc if out of memory.
+         */
         TrajectoryAnalysisModule();
 
         /*! \brief
          * Registers a dataset that exports data.
+         *
+         * \param     data  Data object to register.
+         * \param[in] name  Name to register the dataset with.
+         * \throws    std::bad_alloc if out of memory.
+         *
+         * Registers \p data as a dataset that provides output from the
+         * analysis module.  Callers for the module can access the dataset
+         * with datasetFromName() using \p name as an AbstractAnalysisData
+         * object.  This allows them to add their own data modules to do extra
+         * processing.
+         *
+         * \p name must be unique across all calls within the same
+         * TrajectoryAnalysisModule instance.
          */
         void registerBasicDataset(AbstractAnalysisData *data, const char *name);
         /*! \brief
          * Registers a parallelized dataset that exports data.
+         *
+         * \param     data  AnalysisData object to register.
+         * \param[in] name  Name to register the dataset with.
+         * \throws    std::bad_alloc if out of memory.
+         *
+         * This method works as registerBasicDataset(), but additionally allows
+         * data handles for \p data to be accessed using
+         * TrajectoryAnalysisData.
+         *
+         * \see registerBasicDataset()
          */
         void registerAnalysisDataset(AnalysisData *data, const char *name);
 
index 7c116e8e4a0426cf2116374b2dfef80c8a4d4275..b005a9666f1f78c99452085694f78e4ce62e90a2 100644 (file)
@@ -56,13 +56,19 @@ class TrajectoryAnalysisRunnerCommon;
  *
  * This class is used by trajectory analysis modules to inform the caller
  * about the requirements they have on the input (e.g., whether a topology is
- * required, or whether PBC removal makes sense). It is also used to pass
+ * required, or whether PBC removal makes sense).  It is also used to pass
  * similar information back to the analysis module after parsing user input.
  *
  * Having this functionality as a separate class makes the
  * TrajectoryAnalysisModule interface much cleaner, and also reduces the need to
  * change existing code when new options are added.
  *
+ * Methods in this class do not throw, except for the constructor, which may
+ * throw an std::bad_alloc.
+ *
+ * \todo
+ * Remove plain flags from the public interface.
+ *
  * \inpublicapi
  * \ingroup module_trajectoryanalysis
  */
@@ -104,13 +110,6 @@ class TrajectoryAnalysisSettings
              * \see setRmPBC()
              */
             efNoUserRmPBC    = 1<<5,
-            /*! \brief
-             * Requests dumps of parsed and compiled selection trees.
-             *
-             * This flag is used by internal debugging tools to request
-             * the selection trees dumping to stderr.
-             */
-            efDebugSelection = 1<<16,
         };
 
         //! Initializes default settings.
@@ -225,6 +224,14 @@ class TrajectoryAnalysisSettings
 /*! \brief
  * Topology information passed to a trajectory analysis module.
  *
+ * This class is used to pass topology information to trajectory analysis
+ * modules and to manage memory for them.  Having a single wrapper object
+ * instead of passing each item separately makes TrajectoryAnalysisModule
+ * interface simpler, and also reduces the need to change existing code if
+ * additional information is added.
+ *
+ * Methods in this class do not throw if not explicitly stated.
+ *
  * \inpublicapi
  * \ingroup module_trajectoryanalysis
  */
@@ -246,6 +253,8 @@ class TopologyInformation
          *      (can be NULL, in which case it is not used).
          * \param[out] box   Box size from the topology file
          *      (can be NULL, in which case it is not used).
+         * \throws  APIError if topology coordinates are not available and
+         *      \p x is not NULL.
          *
          * If TrajectoryAnalysisSettings::efUseTopX has not been specified,
          * \p x should be NULL.
@@ -271,6 +280,9 @@ class TopologyInformation
 
         GMX_DISALLOW_COPY_AND_ASSIGN(TopologyInformation);
 
+        /*! \brief
+         * Needed to initialize the data.
+         */
         friend class TrajectoryAnalysisRunnerCommon;
 };
 
index 20bab5bde30150f45e23fc0adbd58e19166a79d5..d4fb89e40dcce10da04a95bb9e7d1c4cadf2dfd9 100644 (file)
@@ -150,17 +150,6 @@ TrajectoryAnalysisCommandLineRunner::Impl::parseOptions(
         return false;
     }
     _module->initOptionsDone(settings);
-    /*
-    if (rc != 0)
-    {
-        if (rc == eeInconsistentInput)
-        {
-            GMX_ERROR(rc, "Invalid options provided, "
-                          "see above for detailed error messages");
-        }
-        return rc;
-    }
-    */
 
     common->initIndexGroups(selections);
 
index a4cce2ecec2c3591cadc38955f9613db8450ce1e..a93c91f1caa23d1a3233f4182361a0b32c8a2579 100644 (file)
@@ -64,6 +64,9 @@ class TrajectoryAnalysisCommandLineRunner
         /*! \brief
          * Create a new runner with the provided module.
          *
+         * \param  module  Analysis module to run using the runner.
+         * \throws std::bad_alloc if out of memory.
+         *
          * The caller should ensure that the provided module is not destroyed
          * while the runner exists.
          */
@@ -75,11 +78,16 @@ class TrajectoryAnalysisCommandLineRunner
          *
          * This is intended only for use by internal debugging tools.
          *
+         * Does not throw.
+         *
          * \see SelectionCollection::setDebugLevel()
          */
         void setSelectionDebugLevel(int debuglevel);
         /*! \brief
          * Parses options from the given command line and runs the analysis.
+         *
+         * \throws  multiple  Exceptions are used to indicate errors.
+         * \returns Zero on success.
          */
         int run(int argc, char *argv[]);
 
index 170a73db92a03af6ca6260e0ed36afd6045d89c8..5659e89e034b853309e51d6ecbdec704dd326e1a 100644 (file)
@@ -52,7 +52,7 @@ namespace gmx
  * \returns  An allocated TrajectoryAnalysisModule object.
  * \throws   InvalidInputError if \p name is not recognized.
  *
- * This function should be used to instantiate selection methods defined in the
+ * This function should be used to instantiate analysis methods defined in the
  * library.
  *
  * In addition to recognizing exact matches on \p name, the function also