2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2011,2012,2013,2014,2015, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
37 * Implements gmx::analysismodules::Angle.
39 * \author Teemu Murtola <teemu.murtola@gmail.com>
40 * \ingroup module_trajectoryanalysis
50 #include "gromacs/analysisdata/analysisdata.h"
51 #include "gromacs/analysisdata/modules/average.h"
52 #include "gromacs/analysisdata/modules/histogram.h"
53 #include "gromacs/analysisdata/modules/plot.h"
54 #include "gromacs/fileio/trx.h"
55 #include "gromacs/math/vec.h"
56 #include "gromacs/options/basicoptions.h"
57 #include "gromacs/options/filenameoption.h"
58 #include "gromacs/options/ioptionscontainer.h"
59 #include "gromacs/options/options.h"
60 #include "gromacs/pbcutil/pbc.h"
61 #include "gromacs/selection/selection.h"
62 #include "gromacs/selection/selectionoption.h"
63 #include "gromacs/trajectoryanalysis/analysissettings.h"
64 #include "gromacs/utility/arrayref.h"
65 #include "gromacs/utility/exceptions.h"
66 #include "gromacs/utility/gmxassert.h"
67 #include "gromacs/utility/stringutil.h"
72 namespace analysismodules
78 /********************************************************************
83 * Helper to encapsulate logic for looping over input selections.
85 * This class provides two-dimensional iteration:
86 * - Over _angle groups_, corresponding to an input selection. If the input
87 * selection list contains a single selection, that selection gets used
88 * for all angle groups.
89 * - Within a group, over _values_, each consisting of a fixed number of
90 * selection positions. If there is only a single value within a selection,
91 * that value is returned over and over again.
92 * This transparently provides the semantics of using a single selection/vector
93 * to compute angles against multiple selections/vectors as described in the
96 * This class isn't perferctly self-contained and requires the caller to know
97 * some of the internals to use it properly, but it serves its purpose for this
98 * single analysis tool by simplifying the loops.
99 * Some methods have also been tailored to allow the caller to use it a bit
102 * \ingroup module_trajectoryanalysis
104 class AnglePositionIterator
108 * Creates an iterator to loop over input selection positions.
110 * \param[in] selections List of selections.
111 * \param[in] posCountPerValue Number of selection positions that
112 * constitute a single value for the iteration.
114 * If \p selections is empty, and/or \p posCountPerValue is zero, the
115 * iterator can still be advanced and hasValue()/hasSingleValue()
116 * called, but values cannot be accessed.
118 AnglePositionIterator(const SelectionList &selections,
119 int posCountPerValue)
120 : selections_(selections), posCountPerValue_(posCountPerValue),
121 currentSelection_(0), nextPosition_(0)
125 //! Advances the iterator to the next group of angles.
128 if (selections_.size() > 1)
134 //! Advances the iterator to the next angle in the current group.
137 if (!hasSingleValue())
139 nextPosition_ += posCountPerValue_;
144 * Returns whether this iterator represents any values.
146 * If the return value is `false`, only nextGroup(), nextValue() and
147 * hasSingleValue() are allowed to be called.
149 bool hasValue() const
151 return !selections_.empty();
154 * Returns whether the current selection only contains a single value.
156 * Returns `false` if hasValue() returns false, which allows cutting
157 * some corners in consistency checks.
159 bool hasSingleValue() const
161 return hasValue() && currentSelection().posCount() == posCountPerValue_;
163 //! Returns whether the current selection is dynamic.
164 bool isDynamic() const
166 return currentSelection().isDynamic();
169 * Returns whether positions in the current value are either all
170 * selected or all unselected.
172 bool allValuesConsistentlySelected() const
174 if (posCountPerValue_ <= 1)
178 const bool bSelected = currentPosition(0).selected();
179 for (int i = 1; i < posCountPerValue_; ++i)
181 if (currentPosition(i).selected() != bSelected)
189 * Returns whether positions in the current value are selected.
191 * Only works reliably if allValuesConsistentlySelected() returns
194 bool currentValuesSelected() const
196 return selections_.empty() || currentPosition(0).selected();
199 //! Returns the currently active selection.
200 const Selection ¤tSelection() const
202 GMX_ASSERT(currentSelection_ < static_cast<int>(selections_.size()),
203 "Accessing an invalid selection");
204 return selections_[currentSelection_];
206 //! Returns the `i`th position for the current value.
207 SelectionPosition currentPosition(int i) const
209 return currentSelection().position(nextPosition_ + i);
212 * Extracts all coordinates corresponding to the current value.
214 * \param[out] x Array to which the positions are extracted.
216 * \p x should contain at minimum the number of positions per value
217 * passed to the constructor.
219 void getCurrentPositions(rvec x[]) const
221 GMX_ASSERT(posCountPerValue_ > 0,
222 "Accessing positions for an invalid angle type");
223 GMX_ASSERT(nextPosition_ + posCountPerValue_ <= currentSelection().posCount(),
224 "Accessing an invalid position");
225 for (int i = 0; i < posCountPerValue_; ++i)
227 copy_rvec(currentPosition(i).x(), x[i]);
232 const SelectionList &selections_;
233 const int posCountPerValue_;
234 int currentSelection_;
237 GMX_DISALLOW_COPY_AND_ASSIGN(AnglePositionIterator);
240 /********************************************************************
241 * Actual analysis module
244 class Angle : public TrajectoryAnalysisModule
249 virtual void initOptions(IOptionsContainer *options,
250 TrajectoryAnalysisSettings *settings);
251 virtual void optionsFinished(Options *options,
252 TrajectoryAnalysisSettings *settings);
253 virtual void initAnalysis(const TrajectoryAnalysisSettings &settings,
254 const TopologyInformation &top);
256 virtual void analyzeFrame(int frnr, const t_trxframe &fr, t_pbc *pbc,
257 TrajectoryAnalysisModuleData *pdata);
259 virtual void finishAnalysis(int nframes);
260 virtual void writeOutput();
263 void initFromSelections(const SelectionList &sel1,
264 const SelectionList &sel2);
265 void checkSelections(const SelectionList &sel1,
266 const SelectionList &sel2) const;
270 SelectionOptionInfo *sel1info_;
271 SelectionOptionInfo *sel2info_;
272 std::string fnAverage_;
274 std::string fnHistogram_;
280 AnalysisData angles_;
281 AnalysisDataFrameAverageModulePointer averageModule_;
282 AnalysisDataSimpleHistogramModulePointer histogramModule_;
284 std::vector<int> angleCount_;
287 std::vector<std::vector<RVec> > vt0_;
289 // Copy and assign disallowed by base.
293 : TrajectoryAnalysisModule(AngleInfo::name, AngleInfo::shortDescription),
294 sel1info_(NULL), sel2info_(NULL), binWidth_(1.0), natoms1_(0), natoms2_(0)
296 averageModule_.reset(new AnalysisDataFrameAverageModule());
297 angles_.addModule(averageModule_);
298 histogramModule_.reset(new AnalysisDataSimpleHistogramModule());
299 angles_.addModule(histogramModule_);
301 registerAnalysisDataset(&angles_, "angle");
302 registerBasicDataset(averageModule_.get(), "average");
303 registerBasicDataset(&histogramModule_->averager(), "histogram");
308 Angle::initOptions(IOptionsContainer *options, TrajectoryAnalysisSettings *settings)
310 static const char *const desc[] = {
311 "[THISMODULE] computes different types of angles between vectors.",
312 "It supports both vectors defined by two positions and normals of",
313 "planes defined by three positions.",
314 "The z axis or the local normal of a sphere can also be used as",
315 "one of the vectors.",
316 "There are also convenience options 'angle' and 'dihedral' for",
317 "calculating bond angles and dihedrals defined by three/four",
319 "The type of the angle is specified with [TT]-g1[tt] and [TT]-g2[tt].",
320 "If [TT]-g1[tt] is [TT]angle[tt] or [TT]dihedral[tt], [TT]-g2[tt]",
321 "should not be specified.",
322 "In this case, [TT]-group1[tt] should specify one or more selections,",
323 "and each should contain triplets or quartets of positions that define",
324 "the angles to be calculated.[PAR]",
325 "If [TT]-g1[tt] is [TT]vector[tt] or [TT]plane[tt], [TT]-group1[tt]",
326 "should specify selections that contain either pairs ([TT]vector[tt])",
327 "or triplets ([TT]plane[tt]) of positions. For vectors, the positions",
328 "set the endpoints of the vector, and for planes, the three positions",
329 "are used to calculate the normal of the plane. In both cases,",
330 "[TT]-g2[tt] specifies the other vector to use (see below).[PAR]",
331 "With [TT]-g2 vector[tt] or [TT]-g2 plane[tt], [TT]-group2[tt] should",
332 "specify another set of vectors. [TT]-group1[tt] and [TT]-group2[tt]",
333 "should specify the same number of selections. It is also allowed to",
334 "only have a single selection for one of the options, in which case",
335 "the same selection is used with each selection in the other group.",
336 "Similarly, for each selection in [TT]-group1[tt], the corresponding",
337 "selection in [TT]-group2[tt] should specify the same number of",
338 "vectors or a single vector. In the latter case, the angle is",
339 "calculated between that single vector and each vector from the other",
341 "With [TT]-g2 sphnorm[tt], each selection in [TT]-group2[tt] should",
342 "specify a single position that is the center of the sphere.",
343 "The second vector is calculated as the vector from the center to the",
344 "midpoint of the positions specified by [TT]-group1[tt].[PAR]",
345 "With [TT]-g2 z[tt], [TT]-group2[tt] is not necessary, and angles",
346 "between the first vectors and the positive Z axis are calculated.[PAR]",
347 "With [TT]-g2 t0[tt], [TT]-group2[tt] is not necessary, and angles",
348 "are calculated from the vectors as they are in the first frame.[PAR]",
349 "There are three options for output:",
350 "[TT]-oav[tt] writes an xvg file with the time and the average angle",
352 "[TT]-oall[tt] writes all the individual angles.",
353 "[TT]-oh[tt] writes a histogram of the angles. The bin width can be",
354 "set with [TT]-binw[tt].",
355 "For [TT]-oav[tt] and [TT]-oh[tt], separate average/histogram is",
356 "computed for each selection in [TT]-group1[tt]."
358 static const char *const cGroup1TypeEnum[] =
359 { "angle", "dihedral", "vector", "plane" };
360 static const char *const cGroup2TypeEnum[] =
361 { "none", "vector", "plane", "t0", "z", "sphnorm" };
363 settings->setHelpText(desc);
365 options->addOption(FileNameOption("oav").filetype(eftPlot).outputFile()
366 .store(&fnAverage_).defaultBasename("angaver")
367 .description("Average angles as a function of time"));
368 options->addOption(FileNameOption("oall").filetype(eftPlot).outputFile()
369 .store(&fnAll_).defaultBasename("angles")
370 .description("All angles as a function of time"));
371 options->addOption(FileNameOption("oh").filetype(eftPlot).outputFile()
372 .store(&fnHistogram_).defaultBasename("anghist")
373 .description("Histogram of the angles"));
375 options->addOption(StringOption("g1").enumValue(cGroup1TypeEnum)
376 .defaultEnumIndex(0).store(&g1type_)
377 .description("Type of analysis/first vector group"));
378 options->addOption(StringOption("g2").enumValue(cGroup2TypeEnum)
379 .defaultEnumIndex(0).store(&g2type_)
380 .description("Type of second vector group"));
381 options->addOption(DoubleOption("binw").store(&binWidth_)
382 .description("Binwidth for -oh in degrees"));
384 sel1info_ = options->addOption(SelectionOption("group1")
385 .required().dynamicMask().storeVector(&sel1_)
387 .description("First analysis/vector selection"));
388 sel2info_ = options->addOption(SelectionOption("group2")
389 .dynamicMask().storeVector(&sel2_)
391 .description("Second analysis/vector selection"));
396 Angle::optionsFinished(Options *options, TrajectoryAnalysisSettings * /* settings */)
398 const bool bSingle = (g1type_[0] == 'a' || g1type_[0] == 'd');
400 if (bSingle && g2type_[0] != 'n')
402 GMX_THROW(InconsistentInputError("Cannot use a second group (-g2) with "
403 "-g1 angle or dihedral"));
405 if (bSingle && options->isSet("group2"))
407 GMX_THROW(InconsistentInputError("Cannot provide a second selection "
408 "(-group2) with -g1 angle or dihedral"));
410 if (!bSingle && g2type_[0] == 'n')
412 GMX_THROW(InconsistentInputError("Should specify a second group (-g2) "
413 "if the first group is not an angle or a dihedral"));
416 // Set up the number of positions per angle.
419 case 'a': natoms1_ = 3; break;
420 case 'd': natoms1_ = 4; break;
421 case 'v': natoms1_ = 2; break;
422 case 'p': natoms1_ = 3; break;
424 GMX_THROW(InternalError("invalid -g1 value"));
428 case 'n': natoms2_ = 0; break;
429 case 'v': natoms2_ = 2; break;
430 case 'p': natoms2_ = 3; break;
431 case 't': natoms2_ = 0; break;
432 case 'z': natoms2_ = 0; break;
433 case 's': natoms2_ = 1; break;
435 GMX_THROW(InternalError("invalid -g2 value"));
437 if (natoms2_ == 0 && options->isSet("group2"))
439 GMX_THROW(InconsistentInputError("Cannot provide a second selection (-group2) with -g2 t0 or z"));
441 // TODO: If bSingle is not set, the second selection option should be
447 Angle::initFromSelections(const SelectionList &sel1,
448 const SelectionList &sel2)
450 const int angleGroups = std::max(sel1.size(), sel2.size());
451 const bool bHasSecondSelection = natoms2_ > 0;
453 if (bHasSecondSelection && sel1.size() != sel2.size()
454 && std::min(sel1.size(), sel2.size()) != 1)
456 GMX_THROW(InconsistentInputError(
457 "-group1 and -group2 should specify the same number of selections"));
460 AnglePositionIterator iter1(sel1, natoms1_);
461 AnglePositionIterator iter2(sel2, natoms2_);
462 for (int g = 0; g < angleGroups; ++g, iter1.nextGroup(), iter2.nextGroup())
464 const int posCount1 = iter1.currentSelection().posCount();
465 if (natoms1_ > 1 && posCount1 % natoms1_ != 0)
467 GMX_THROW(InconsistentInputError(formatString(
468 "Number of positions in selection %d in the first group not divisible by %d",
469 static_cast<int>(g + 1), natoms1_)));
471 const int angleCount1 = posCount1 / natoms1_;
472 int angleCount = angleCount1;
474 if (bHasSecondSelection)
476 const int posCount2 = iter2.currentSelection().posCount();
477 if (natoms2_ > 1 && posCount2 % natoms2_ != 0)
479 GMX_THROW(InconsistentInputError(formatString(
480 "Number of positions in selection %d in the second group not divisible by %d",
481 static_cast<int>(g + 1), natoms2_)));
483 if (g2type_[0] == 's' && posCount2 != 1)
485 GMX_THROW(InconsistentInputError(
486 "The second group should contain a single position with -g2 sphnorm"));
489 const int angleCount2 = posCount2 / natoms2_;
490 angleCount = std::max(angleCount1, angleCount2);
491 if (angleCount1 != angleCount2
492 && std::min(angleCount1, angleCount2) != 1)
494 GMX_THROW(InconsistentInputError(
495 "Number of vectors defined by the two groups are not the same"));
498 angleCount_.push_back(angleCount);
504 Angle::checkSelections(const SelectionList &sel1,
505 const SelectionList &sel2) const
507 AnglePositionIterator iter1(sel1, natoms1_);
508 AnglePositionIterator iter2(sel2, natoms2_);
509 for (size_t g = 0; g < angleCount_.size(); ++g, iter1.nextGroup(), iter2.nextGroup())
511 if (iter1.isDynamic() || (iter2.hasValue() && iter2.isDynamic()))
513 for (int n = 0; n < angleCount_[g]; ++n, iter1.nextValue(), iter2.nextValue())
516 if (!iter1.allValuesConsistentlySelected())
520 if (!iter2.allValuesConsistentlySelected())
524 if (angleCount_[g] > 1)
526 if (iter1.hasSingleValue() && !iter1.currentValuesSelected())
530 if (iter2.hasSingleValue() && !iter2.currentValuesSelected())
536 && (angleCount_[g] == 1
537 || (!iter1.hasSingleValue() && !iter2.hasSingleValue()))
538 && iter1.currentValuesSelected() != iter2.currentValuesSelected())
544 std::string message =
545 formatString("Dynamic selection %d does not select "
546 "a consistent set of angles over the frames",
547 static_cast<int>(g + 1));
548 GMX_THROW(InconsistentInputError(message));
557 Angle::initAnalysis(const TrajectoryAnalysisSettings &settings,
558 const TopologyInformation & /* top */)
560 initFromSelections(sel1_, sel2_);
562 // checkSelections() ensures that both selection lists have the same size.
563 angles_.setDataSetCount(angleCount_.size());
564 for (size_t i = 0; i < angleCount_.size(); ++i)
566 angles_.setColumnCount(i, angleCount_[i]);
568 double histogramMin = (g1type_ == "dihedral" ? -180.0 : 0);
569 histogramModule_->init(histogramFromRange(histogramMin, 180.0)
570 .binWidth(binWidth_).includeAll());
574 vt0_.resize(sel1_.size());
575 for (size_t g = 0; g < sel1_.size(); ++g)
577 vt0_[g].resize(sel1_[g].posCount() / natoms1_);
581 if (!fnAverage_.empty())
583 AnalysisDataPlotModulePointer plotm(
584 new AnalysisDataPlotModule(settings.plotSettings()));
585 plotm->setFileName(fnAverage_);
586 plotm->setTitle("Average angle");
587 plotm->setXAxisIsTime();
588 plotm->setYLabel("Angle (degrees)");
589 // TODO: Consider adding information about the second selection,
590 // and/or a subtitle describing what kind of angle this is.
591 for (size_t g = 0; g < sel1_.size(); ++g)
593 plotm->appendLegend(sel1_[g].name());
595 averageModule_->addModule(plotm);
600 AnalysisDataPlotModulePointer plotm(
601 new AnalysisDataPlotModule(settings.plotSettings()));
602 plotm->setFileName(fnAll_);
603 plotm->setTitle("Angle");
604 plotm->setXAxisIsTime();
605 plotm->setYLabel("Angle (degrees)");
606 // TODO: Add legends? (there can be a massive amount of columns)
607 angles_.addModule(plotm);
610 if (!fnHistogram_.empty())
612 AnalysisDataPlotModulePointer plotm(
613 new AnalysisDataPlotModule(settings.plotSettings()));
614 plotm->setFileName(fnHistogram_);
615 plotm->setTitle("Angle histogram");
616 plotm->setXLabel("Angle (degrees)");
617 plotm->setYLabel("Probability");
618 // TODO: Consider adding information about the second selection,
619 // and/or a subtitle describing what kind of angle this is.
620 for (size_t g = 0; g < sel1_.size(); ++g)
622 plotm->appendLegend(sel1_[g].name());
624 histogramModule_->averager().addModule(plotm);
629 //! Helper method to calculate a vector from two or three positions.
631 calc_vec(int natoms, rvec x[], t_pbc *pbc, rvec xout, rvec cout)
638 pbc_dx(pbc, x[1], x[0], xout);
642 rvec_sub(x[1], x[0], xout);
644 svmul(0.5, xout, cout);
645 rvec_add(x[0], cout, cout);
652 pbc_dx(pbc, x[1], x[0], v1);
653 pbc_dx(pbc, x[2], x[0], v2);
657 rvec_sub(x[1], x[0], v1);
658 rvec_sub(x[2], x[0], v2);
661 rvec_add(x[0], x[1], cout);
662 rvec_add(cout, x[2], cout);
663 svmul(1.0/3.0, cout, cout);
667 GMX_RELEASE_ASSERT(false, "Incorrectly initialized number of atoms");
673 Angle::analyzeFrame(int frnr, const t_trxframe &fr, t_pbc *pbc,
674 TrajectoryAnalysisModuleData *pdata)
676 AnalysisDataHandle dh = pdata->dataHandle(angles_);
677 const SelectionList &sel1 = pdata->parallelSelections(sel1_);
678 const SelectionList &sel2 = pdata->parallelSelections(sel2_);
680 checkSelections(sel1, sel2);
682 dh.startFrame(frnr, fr.time);
684 AnglePositionIterator iter1(sel1, natoms1_);
685 AnglePositionIterator iter2(sel2, natoms2_);
686 for (size_t g = 0; g < angleCount_.size(); ++g, iter1.nextGroup(), iter2.nextGroup())
691 // v2 & c2 are conditionally set in the switch statement below, and conditionally
692 // used in a different switch statement later. Apparently the clang static analyzer
693 // thinks there are cases where they can be used uninitialzed (which I cannot find),
694 // but to avoid trouble if we ever change just one of the switch statements it
695 // makes sense to clear them outside the first switch.
706 copy_rvec(sel2_[g].position(0).x(), c2);
714 for (int n = 0; n < angleCount_[g]; ++n, iter1.nextValue(), iter2.nextValue())
717 // x[] will be assigned below based on the number of atoms used to initialize iter1,
718 // which in turn should correspond perfectly to g1type_[0] (which determines how many we read),
719 // but unsurprisingly the static analyzer chokes a bit on that.
723 // checkSelections() ensures that this reflects all the involved
725 const bool bPresent =
726 iter1.currentValuesSelected() && iter2.currentValuesSelected();
727 iter1.getCurrentPositions(x);
733 pbc_dx(pbc, x[0], x[1], v1);
734 pbc_dx(pbc, x[2], x[1], v2);
738 rvec_sub(x[0], x[1], v1);
739 rvec_sub(x[2], x[1], v2);
741 angle = gmx_angle(v1, v2);
748 pbc_dx(pbc, x[0], x[1], dx[0]);
749 pbc_dx(pbc, x[2], x[1], dx[1]);
750 pbc_dx(pbc, x[2], x[3], dx[2]);
754 rvec_sub(x[0], x[1], dx[0]);
755 rvec_sub(x[2], x[1], dx[1]);
756 rvec_sub(x[2], x[3], dx[2]);
758 cprod(dx[0], dx[1], v1);
759 cprod(dx[1], dx[2], v2);
760 angle = gmx_angle(v1, v2);
761 real ipr = iprod(dx[0], v2);
770 calc_vec(natoms1_, x, pbc, v1, c1);
775 iter2.getCurrentPositions(x);
776 calc_vec(natoms2_, x, pbc, v2, c2);
779 // FIXME: This is not parallelizable.
782 copy_rvec(v1, vt0_[g][n]);
784 copy_rvec(vt0_[g][n], v2);
787 c1[XX] = c1[YY] = 0.0;
792 pbc_dx(pbc, c1, c2, v2);
796 rvec_sub(c1, c2, v2);
800 GMX_THROW(InternalError("invalid -g2 value"));
802 angle = gmx_angle(v1, v2);
805 GMX_THROW(InternalError("invalid -g1 value"));
807 dh.setPoint(n, angle * RAD2DEG, bPresent);
815 Angle::finishAnalysis(int /*nframes*/)
817 AbstractAverageHistogram &averageHistogram = histogramModule_->averager();
818 averageHistogram.normalizeProbability();
819 averageHistogram.done();
830 const char AngleInfo::name[] = "gangle";
831 const char AngleInfo::shortDescription[] =
834 TrajectoryAnalysisModulePointer AngleInfo::create()
836 return TrajectoryAnalysisModulePointer(new Angle);
839 } // namespace analysismodules