ProteoWizard
OverlapDemultiplexer.hpp
Go to the documentation of this file.
1 //
2 // $Id$
3 //
4 //
5 // Original author: Jarrett Egertson <jegertso .@. uw.edu>
6 //
7 // Licensed under the Apache License, Version 2.0 (the "License");
8 // you may not use this file except in compliance with the License.
9 // You may obtain a copy of the License at
10 //
11 // http://www.apache.org/licenses/LICENSE-2.0
12 //
13 // Unless required by applicable law or agreed to in writing, software
14 // distributed under the License is distributed on an "AS IS" BASIS,
15 // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
16 // See the License for the specific language governing permissions and
17 // limitations under the License.
18 //
19 
20 #ifndef _OVERLAPDEMULTIPLEXER_HPP
21 #define _OVERLAPDEMULTIPLEXER_HPP
22 
23 #include "IDemultiplexer.hpp"
24 #include "DemuxHelpers.hpp"
25 
26 namespace pwiz {
27 namespace analysis {
28 
29  /// Implementation of the IDemultiplexer interface that is able to handle overlap experiments.
31  {
32  public:
33 
34  /// User-defined options for demultiplexing
35  struct Params
36  {
37  Params() :
38 
39  massError(10.0, pwiz::chemistry::MZTolerance::PPM),
40  applyWeighting(true),
41  variableFill(false)
42  {}
43 
44  /// Mass error for extracting MS/MS peaks
46 
47  /// Weight the spectra nearby to the input spectrum more heavily in the solve
48  /// than the outer ones
50 
51  /// Set to true if fill times are allowed to vary for each scan window
53  };
54  /// Constructs an OverlapDemultiplexer with optional user-specified parameters
55  /// @param p Options to use in demultiplexing (see Params for available options)
56  explicit OverlapDemultiplexer(Params p = Params());
57 
58  virtual ~OverlapDemultiplexer();
59 
60  /// \name IDemultiplexer interface
61  ///@{
62 
64  void BuildDeconvBlock(size_t index,
65  const std::vector<size_t>& muxIndices,
66  DemuxTypes::MatrixPtr& masks,
67  DemuxTypes::MatrixPtr& signal) override;
68  void GetMatrixBlockIndices(size_t indexToDemux, std::vector<size_t>& muxIndices, double demuxBlockExtra) const override;
69  const std::vector<size_t>& SpectrumIndices() const override;
70  ///@}
71 
72  protected:
73 
74  /// Performs interpolation on a matrix of intensities using a vector of scanTimes and outputs them to a row vector of interpolated intensities
75  /// @param[out] interpolatedIntensities The row vector of interpolated intensities mapped from timeToInterpolate.
76  /// The Ref template convinces Eigen that this can be column or row.
77  /// @param[in] timeToInterpolate The time corresponding to the scan to be demuxed
78  /// @param[in] intensities A matrix of intensities of a number of nearby spectra. Each spectrum should have a row of transitions.
79  /// @param[in] scanTimes A vector of scanTimes. Because the same scan times are used for every transition in a given spectrum it is only
80  /// necessary to pass a vector of scanTimes rather than a matrix. This vector is reused for interpolating every transition.
81  static void InterpolateMuxRegion(
82  Eigen::Ref<Eigen::MatrixXd, 0, Eigen::Stride<Eigen::Dynamic, Eigen::Dynamic> > interpolatedIntensities,
83  double timeToInterpolate,
84  Eigen::Ref<const Eigen::MatrixXd> intensities,
85  Eigen::Ref<const Eigen::VectorXd> scanTimes);
86 
87  /// Takes two vectors of equal length and solves an interpolation for the given point.
88  /// @param pointToInterpolate Independent variable to interpolate
89  /// @param points The independent values as a monotonically increasing series
90  /// @param values The dependent values
91  /// @return Returns the solved interpolation value
92  /// \pre points must be sorted in order of increasing value with no duplicates
93  /// \pre points and values must be of the same size
94  static double InterpolateMatrix(double pointToInterpolate, Eigen::Ref<const Eigen::VectorXd> points, Eigen::Ref<const Eigen::VectorXd> values);
95 
96  private:
97 
98  /// The number of mux spectra nearby the spectrum to demux (in both retention time and m/z space) to use for demuxing
100 
101  /// The number of spectra with identical isolation parameters to use for interpolation
103 
104  /// A SpectrumList that provides access to the spectra specified in the muxIndices list provided to BuildDeconvBlock()
106 
107  /// An IPrecursorMaskCodec that provides information about the experiment's scheme and can generate the masks for given mux spectra
109 
110  /// A set of user-defined options
112 
113  /// A cache of the indices provided by SpectrumIndices()
114  std::vector<size_t> spectrumIndices_;
115  };
116 } // namespace analysis
117 } // namespace pwiz
118 #endif // _OVERLAPDEMULTIPLEXER_HPP
const std::vector< size_t > & SpectrumIndices() const override
Returns the indices to the demultiplexed windows in the solution matrix corresponding to the windows ...
Params params_
A set of user-defined options.
User-defined options for demultiplexing.
static void InterpolateMuxRegion(Eigen::Ref< Eigen::MatrixXd, 0, Eigen::Stride< Eigen::Dynamic, Eigen::Dynamic > > interpolatedIntensities, double timeToInterpolate, Eigen::Ref< const Eigen::MatrixXd > intensities, Eigen::Ref< const Eigen::VectorXd > scanTimes)
Performs interpolation on a matrix of intensities using a vector of scanTimes and outputs them to a r...
Implementation of the IDemultiplexer interface that is able to handle overlap experiments.
OverlapDemultiplexer(Params p=Params())
Constructs an OverlapDemultiplexer with optional user-specified parameters.
size_t overlapRegionsInApprox_
The number of mux spectra nearby the spectrum to demux (in both retention time and m/z space) to use ...
boost::shared_ptr< MatrixType > MatrixPtr
Definition: DemuxTypes.hpp:39
Interface for calculating demultiplexing scheme.
bool applyWeighting
Weight the spectra nearby to the input spectrum more heavily in the solve than the outer ones...
IPrecursorMaskCodec::const_ptr pmc_
An IPrecursorMaskCodec that provides information about the experiment&#39;s scheme and can generate the m...
Helper functions for demultiplexing Helper functions include nice methods of accessing CV parameters ...
bool variableFill
Set to true if fill times are allowed to vary for each scan window.
pwiz::chemistry::MZTolerance massError
Mass error for extracting MS/MS peaks.
size_t cyclesInBlock_
The number of spectra with identical isolation parameters to use for interpolation.
boost::shared_ptr< const IPrecursorMaskCodec > const_ptr
Constant shared pointer definition.
void BuildDeconvBlock(size_t index, const std::vector< size_t > &muxIndices, DemuxTypes::MatrixPtr &masks, DemuxTypes::MatrixPtr &signal) override
Translates a spectrum into a set of matrices to be solved by NNLS.
std::vector< size_t > spectrumIndices_
A cache of the indices provided by SpectrumIndices()
void Initialize(msdata::SpectrumList_const_ptr sl, IPrecursorMaskCodec::const_ptr pmc) override
Initializes the demultiplexer using the demux scheme provided by an IPrecursorMaskCodec.
struct for expressing m/z tolerance in either amu or ppm
Definition: MZTolerance.hpp:38
msdata::SpectrumList_const_ptr sl_
A SpectrumList that provides access to the spectra specified in the muxIndices list provided to Build...
static double InterpolateMatrix(double pointToInterpolate, Eigen::Ref< const Eigen::VectorXd > points, Eigen::Ref< const Eigen::VectorXd > values)
Takes two vectors of equal length and solves an interpolation for the given point.
boost::shared_ptr< const msdata::SpectrumList > SpectrumList_const_ptr
Definition: DemuxTypes.hpp:28
void GetMatrixBlockIndices(size_t indexToDemux, std::vector< size_t > &muxIndices, double demuxBlockExtra) const override
Figures out which spectra to include in the system of equations to demux.