ProteoWizard
PrecursorMaskCodec.hpp
Go to the documentation of this file.
1 //
2 // $Id: PrecursorMaskCodec.hpp 10042 2016-09-27 22:27:10Z chambm $
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 _PRECURSORMASKCODEC_HPP
21 #define _PRECURSORMASKCODEC_HPP
22 
23 #include "IPrecursorMaskCodec.hpp"
24 
25 namespace pwiz{
26 namespace analysis{
27 
28  /// Implementation of the IPrecursorMaskCodec interface that is able to handle both overlapping MSX experiments. Some features that are only
29  /// applicable when overlapping is used without MSX or MSX used without overlapping are not implemented in this class. Such features and optimizations
30  /// are currently left to other more targeted implementations of IPrecursorMaskCodec. One missing feature is interpolation of overlap-only data to
31  /// optimize weighting of nearby spectra before demultiplexing.
33  {
34  public:
35 
36  /// Construct a PrecursorMaskCodec for interpreting overlapping and MSX experiments for demultiplexing.
37  /// @param[in] slPtr SpectrumList to demux
38  /// @param[in] variableFill Set to true if fill times are allowed to vary for each scan window
39  explicit PrecursorMaskCodec(msdata::SpectrumList_const_ptr slPtr, bool variableFill = false);
40  virtual ~PrecursorMaskCodec(){}
41 
42  /// \name IPrecursorMaskCodec interface
43  ///@{
44 
45  Eigen::VectorXd GetMask(msdata::Spectrum_const_ptr sPtr, double weight) const override;
46  void GetMask(msdata::Spectrum_const_ptr sPtr, DemuxTypes::MatrixType& m, size_t rowNum, double weight) const override;
47  void SpectrumToIndices(msdata::Spectrum_const_ptr spectrumPtr, std::vector<size_t>& indices) const override;
48  IsolationWindow GetIsolationWindow(size_t i) const override;
49  size_t GetNumDemuxWindows() const override;
50  int GetSpectraPerCycle() const override;
51  int GetPrecursorsPerSpectrum() const override;
52  int GetOverlapsPerCycle() const override;
53  size_t GetDemuxBlockSize() const override;
55  ///@}
56 
57  protected:
58 
59 
60  /// Simple container that is useful for breaking up DemuxWindows into their edges and resolving overlap
62  {
63  /// Constructs a DemuxBoundary from an m/z floating point value
64  explicit DemuxBoundary(double mz) : mz(mz), mzHash(IsoWindowHasher::Hash(mz)) {}
65 
66  double mz; ///< Full precision m/z value
67 
68  MZHash mzHash; ///< Hashed m/z value for fast and simple comparison operations
69 
70  /// DemuxBoundaries are sorted to the precision of their hash
71  bool operator<(const DemuxBoundary& rhs) const { return this->mzHash < rhs.mzHash; }
72 
73  /// DemuxBoundaries are equated only by their hashes
74  bool operator==(const DemuxBoundary& rhs) const { return this->mzHash == rhs.mzHash; }
75  };
76 
77  /// Interpret the experimental design of the multiplexed experiment and cache values for building the design matrix when later given spectra.
79 
80  /// Identifies the repeating scan pattern in the experiment and extracts features of the experimental design in order to interpret the intended demux
81  /// scheme. Note that an alternative to this would be to have the experimental design specified by the user or written in some form in metadata.
82  /// @param[in] spectrumList The SpectrumListPtr containing the multiplexed experiment
83  /// @param[out] demuxWindows The largest set of windows that the experiment can be demultiplexed into when not accounting for overlap.
84  /// \post demuxWindows is sorted and contains no duplicate elements
85  void IdentifyCycle(msdata::SpectrumList_const_ptr spectrumList, std::vector<IsolationWindow>& demuxWindows);
86 
87  /// Identifies any overlap in a DemuxWindow set and splits any overlapping regions such that a non-overlapping DemuxWindow set is produced.
88  /// @param[in] demuxWindows Set of possibly overlapping DemuxWindows.
89  /// @param[out] demuxWindows Set of non-overlapping DemuxWindows produced from the input demuxWindows.
90  /// \pre demuxWindows is sorted and has no duplicate elements.
91  /// \post demuxWindows is output as a vector of size equal to or greater than the input vector.
92  void IdentifyOverlap(std::vector<IsolationWindow>& demuxWindows);
93 
94  private:
95 
96  /// Template for retrieving masks from array-like objects with [] accessors. Rvalue references are used for compatibility with Eigen types.
97  /// See IPrecursorMaskCodec::GetMask
98  /// @param[out] arrayType Array-like object with [] accessor
99  /// @param[in] sPtr Multiplexed spectrum from which to extract precursor windows.
100  /// @param[in] weight Scalar value by which to weight the resulting design matrix vector.
101  /// This weighting is a simple scalar multiplication of the vector.
102  template <class T>
103  void GetMask(T&& arrayType, msdata::Spectrum_const_ptr sPtr, double weight) const;
104 
105  /// This is effectively the index to isolation window map for translating isolation windows to
106  std::vector<IsolationWindow> isolationWindows_;
107 
108  /// Number of spectra required to cover all precursor windows. This is the number of spectra required to make a fully determined
109  /// system of equations.
111 
112  /// Number of precursors (or isolation windows) per multiplexed spectrum. This is calculated in ReadDemuxScheme() and is assumed to be constant
113  /// for all spectra. An error is thrown if this is ever observed to change.
115 
116  /// Number of overlap windows per multiplexed spectrum. E.g. having one additional round of acquisition with an overlap offset would result in
117  /// two overlap regions per spectrum. Alternatively, in the case of no overlap the number of "overlap windows" is one. This is calculated in
118  /// ReadDemuxScheme() and assumed to be constant for all spectra. An error is thrown if this is ever observed to change.
120 
121  /// Whether this data acquired with variable fill times or not. This is set by the user.
123 
124  /// Cached processing method to return from GetProcessingMethod()
126  };
127 } // namespace analysis
128 } // namespace pwiz
129 #endif // _PRECURSORMASKCODEC_HPP
A method of hashing an isolation window to a unique long value mz is and m/z of a unique point in the...
size_t overlapsPerSpectrum_
Number of overlap windows per multiplexed spectrum.
msdata::ProcessingMethod processingMethod_
Cached processing method to return from GetProcessingMethod()
int GetSpectraPerCycle() const override
Returns the number of spectra required to cover all precursor isolation windows.
void SpectrumToIndices(msdata::Spectrum_const_ptr spectrumPtr, std::vector< size_t > &indices) const override
Identifies the precursor windows within a spectrum and returns the indices to the design matrix colum...
size_t spectraPerCycle_
Number of spectra required to cover all precursor windows.
int GetOverlapsPerCycle() const override
Returns the number of overlap repeats per cycle. So for no overlap, this returns 1. For an overlap that splits each precursor in two, this returns 2. Etc.
A container that wraps DemuxWindow to preserve the full precision window boundaries.
MZHash mzHash
Hashed m/z value for fast and simple comparison operations.
Eigen::VectorXd GetMask(msdata::Spectrum_const_ptr sPtr, double weight) const override
Generates a design matrix row describing which precursor isolation windows are present in the given s...
bool operator==(const DemuxBoundary &rhs) const
DemuxBoundaries are equated only by their hashes.
Implementation of the IPrecursorMaskCodec interface that is able to handle both overlapping MSX exper...
msdata::ProcessingMethod GetProcessingMethod() const override
Returns a descriptor of the processing done by this PrecursorMaskCodec.
Interface for generating and accessing precursor masks for a demultiplexing scheme.
IsolationWindow GetIsolationWindow(size_t i) const override
Returns the precursor window for a given index.
Description of the default peak processing method. This element describes the base method used in the...
Definition: MSData.hpp:253
bool variableFill_
Whether this data acquired with variable fill times or not. This is set by the user.
size_t GetNumDemuxWindows() const override
Returns the total number of demux&#39;d precursor windows. This is the number of possible indices returne...
bool operator<(const DemuxBoundary &rhs) const
DemuxBoundaries are sorted to the precision of their hash.
int GetPrecursorsPerSpectrum() const override
Returns the number of precursor isolations per spectrum. This is verified to be constant for all spec...
void IdentifyOverlap(std::vector< IsolationWindow > &demuxWindows)
Identifies any overlap in a DemuxWindow set and splits any overlapping regions such that a non-overla...
size_t precursorsPerSpectrum_
Number of precursors (or isolation windows) per multiplexed spectrum.
Matrix< DemuxScalar, Dynamic, Dynamic > MatrixType
Definition: DemuxTypes.hpp:38
boost::shared_ptr< const msdata::Spectrum > Spectrum_const_ptr
Definition: DemuxTypes.hpp:29
size_t GetDemuxBlockSize() const override
Returns the number of windows required to demultiplex.
std::vector< IsolationWindow > isolationWindows_
This is effectively the index to isolation window map for translating isolation windows to...
PrecursorMaskCodec(msdata::SpectrumList_const_ptr slPtr, bool variableFill=false)
Construct a PrecursorMaskCodec for interpreting overlapping and MSX experiments for demultiplexing...
void IdentifyCycle(msdata::SpectrumList_const_ptr spectrumList, std::vector< IsolationWindow > &demuxWindows)
Identifies the repeating scan pattern in the experiment and extracts features of the experimental des...
Simple container that is useful for breaking up DemuxWindows into their edges and resolving overlap...
DemuxBoundary(double mz)
Constructs a DemuxBoundary from an m/z floating point value.
boost::shared_ptr< const msdata::SpectrumList > SpectrumList_const_ptr
Definition: DemuxTypes.hpp:28
void ReadDemuxScheme(msdata::SpectrumList_const_ptr spectrumList)
Interpret the experimental design of the multiplexed experiment and cache values for building the des...