libpappsomspp
Library for mass spectrometry
tracedetectionzivy.cpp
Go to the documentation of this file.
1
2/*******************************************************************************
3 * Copyright (c) 2015 Olivier Langella <Olivier.Langella@moulon.inra.fr>.
4 *
5 * This file is part of the PAPPSOms++ library.
6 *
7 * PAPPSOms++ is free software: you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation, either version 3 of the License, or
10 * (at your option) any later version.
11 *
12 * PAPPSOms++ is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with PAPPSOms++. If not, see <http://www.gnu.org/licenses/>.
19 *
20 * Contributors:
21 * Olivier Langella <Olivier.Langella@moulon.inra.fr> - initial API and
22 *implementation
23 ******************************************************************************/
24
25#include "tracedetectionzivy.h"
26#include "tracepeak.h"
27#include "../../exception/exceptionoutofrange.h"
28#include "../../exception/exceptionnotpossible.h"
29
30#include <QObject>
31
32
33namespace pappso
34{
36 unsigned int smoothing_half_window_length,
37 unsigned int minmax_half_window_length,
38 unsigned int maxmin_half_window_length,
39 pappso_double detection_threshold_on_minmax,
40 pappso_double detection_threshold_on_maxmin)
41 : m_smooth(smoothing_half_window_length),
42 m_minMax(minmax_half_window_length),
43 m_maxMin(maxmin_half_window_length)
44{
45 m_detectionThresholdOnMaxMin = detection_threshold_on_maxmin;
46 m_detectionThresholdOnMinMax = detection_threshold_on_minmax;
47}
49{
50}
51
52void
54{
55 m_smooth = smooth;
56}
57void
59{
60 m_minMax = minMax;
61}
62void
64{
65 m_maxMin = maxMin;
66}
67
68void
70 double detectionThresholdOnMinMax)
71{
72 m_detectionThresholdOnMinMax = detectionThresholdOnMinMax;
73}
74void
76 double detectionThresholdOnMaxMin)
77{
78 m_detectionThresholdOnMaxMin = detectionThresholdOnMaxMin;
79}
80unsigned int
82{
84};
85unsigned int
87{
89};
90
91unsigned int
93{
95};
98{
100};
103{
105};
106
107
108void
111 bool remove_peak_base) const
112{
113
114 // detect peak positions on close curve : a peak is an intensity value
115 // strictly greater than the two surrounding values. In case of
116 // equality (very rare, can happen with some old old spectrometers) we
117 // take the last equal point to be the peak
118 qDebug();
119 std::size_t size = xic.size();
120 if(size < 4)
121 {
122 qDebug() << QObject::tr(
123 "The original XIC is too small to detect peaks (%1)")
124 .arg(xic.size());
125 return;
126 }
128 return;
130 return;
131
132 Trace xic_minmax(xic); //"close" courbe du haut
134 {
135 qDebug() << "f";
136 m_smooth.filter(xic_minmax);
137 }
138
139 qDebug();
140 Trace xic_maxmin(xic_minmax); //"open" courbe du bas
141
142 qDebug();
143
144 try
145 {
146 qDebug() << "f1";
147 m_minMax.filter(xic_minmax);
148 qDebug() << "f2";
149 m_maxMin.filter(xic_maxmin);
150 }
151 catch(const ExceptionOutOfRange &e)
152 {
153 qDebug() << QObject::tr(
154 "The original XIC is too small to detect peaks (%1) :\n%2")
155 .arg(xic.size())
156 .arg(e.qwhat());
157 return;
158 }
159 qDebug() << "a";
160
161 std::vector<DataPoint>::const_iterator previous_rt = xic_minmax.begin();
162 std::vector<DataPoint>::const_iterator current_rt = (previous_rt + 1);
163 std::vector<DataPoint>::const_iterator last_rt = (xic_minmax.end() - 1);
164
165 std::vector<DataPoint>::const_iterator current_rt_on_maxmin =
166 (xic_maxmin.begin() + 1);
167
168 std::vector<DataPoint>::const_iterator xic_position = xic.begin();
169 qDebug() << "b";
170 while(current_rt != last_rt)
171 // for (unsigned int i = 1, count = 0; i < xic_minmax.size() - 1; )
172 {
173 // conditions to have a peak
174 if((previous_rt->y < current_rt->y) &&
175 (current_rt->y > m_detectionThresholdOnMinMax) &&
176 (current_rt_on_maxmin->y > m_detectionThresholdOnMaxMin))
177 {
178 // here we test the last condition to have a peak
179
180 // no peak case
181 if(current_rt->y < (current_rt + 1)->y)
182 {
183 //++i;
184 previous_rt = current_rt;
185 current_rt++;
186 current_rt_on_maxmin++;
187 }
188 // there is a peak here ! case
189 else if(current_rt->y > (current_rt + 1)->y)
190 {
191 // peak.setMaxXicElement(*current_rt);
192
193 // find left boundary
194 std::vector<DataPoint>::const_iterator it_left =
195 moveLowerYLeftDataPoint(xic_minmax, current_rt);
196 // walk back
197 it_left = moveLowerYRigthDataPoint(xic_minmax, it_left);
198 // peak.setLeftBoundary(*it_left);
199
200 // find right boundary
201 std::vector<DataPoint>::const_iterator it_right =
202 moveLowerYRigthDataPoint(xic_minmax, current_rt);
203 // walk back
204 it_right = moveLowerYLeftDataPoint(xic_minmax, it_right);
205 // peak.setRightBoundary(*it_right);
206
207 // integrate peak surface :
208 auto it =
209 findFirstEqualOrGreaterX(xic_position, xic.end(), it_left->x);
210 xic_position =
211 findFirstEqualOrGreaterX(it, xic.end(), it_right->x) + 1;
212 // peak.setArea(areaTrace(it, xic_position));
213
214 // find the maximum :
215 // peak.setMaxXicElement(*maxYDataPoint(it, xic_position));
216
217 // areaTrace()
218 TracePeak peak(it, xic_position, remove_peak_base);
219 sink.setTracePeak(peak);
220 // }
221 //++i;
222 previous_rt = current_rt;
223 current_rt++;
224 current_rt_on_maxmin++;
225 }
226 // equality case, skipping equal points
227 else
228 {
229 // while (v_minmax[i] == v_minmax[i + 1]) {
230 //++i;
231 current_rt++;
232 current_rt_on_maxmin++;
233
234 //++count;
235 }
236 }
237 // no chance to have a peak at all, continue looping
238 else
239 {
240 //++i;
241 previous_rt = current_rt;
242 current_rt++;
243 current_rt_on_maxmin++;
244 }
245 } // end loop for peaks
246 qDebug();
247}
248} // namespace pappso
transform the trace with the maximum of the minimum equivalent of the erode filter for pictures
Definition: filtermorpho.h:138
Trace & filter(Trace &data_points) const override
std::size_t getMaxMinHalfEdgeWindows() const
mean filter apply mean of y values inside the window : this results in a kind of smoothing
Definition: filtermorpho.h:210
std::size_t getMeanHalfEdgeWindows() const
transform the trace with the minimum of the maximum equivalent of the dilate filter for pictures
Definition: filtermorpho.h:118
std::size_t getMinMaxHalfEdgeWindows() const
Trace & filter(Trace &data_points) const override
virtual Trace & filter(Trace &data_points) const override
virtual const QString & qwhat() const
virtual void setTracePeak(TracePeak &xic_peak)=0
void setFilterMorphoMinMax(const FilterMorphoMinMax &m_minMax)
void setDetectionThresholdOnMaxmin(double detectionThresholdOnMaxMin)
unsigned int getMinMaxHalfEdgeWindows() const
void detect(const Trace &xic, TraceDetectionSinkInterface &sink, bool remove_peak_base) const override
detect peaks on a trace
unsigned int getSmoothingHalfEdgeWindows() const
pappso_double getDetectionThresholdOnMinmax() const
pappso_double getDetectionThresholdOnMaxmin() const
void setFilterMorphoMaxMin(const FilterMorphoMaxMin &maxMin)
TraceDetectionZivy(unsigned int smoothing_half_window_length, unsigned int minmax_half_window_length, unsigned int maxmin_half_window_length, pappso_double detection_threshold_on_minmax, pappso_double detection_threshold_on_maxmin)
pappso_double m_detectionThresholdOnMaxMin
void setFilterMorphoMean(const FilterMorphoMean &smooth)
pappso_double m_detectionThresholdOnMinMax
unsigned int getMaxMinHalfEdgeWindows() const
void setDetectionThresholdOnMinmax(double detectionThresholdOnMinMax)
A simple container of DataPoint instances.
Definition: trace.h:148
tries to keep as much as possible monoisotopes, removing any possible C13 peaks and changes multichar...
Definition: aa.cpp:39
std::vector< DataPoint >::iterator findFirstEqualOrGreaterX(std::vector< DataPoint >::iterator begin, std::vector< DataPoint >::iterator end, const double &value)
find the first element in which X is equal or greater than the value searched important : it implies ...
Definition: trace.cpp:71
std::vector< DataPoint >::const_iterator moveLowerYLeftDataPoint(const Trace &trace, std::vector< DataPoint >::const_iterator begin)
Move left to the lower value.
Definition: trace.cpp:223
double pappso_double
A type definition for doubles.
Definition: types.h:49
std::vector< DataPoint >::const_iterator moveLowerYRigthDataPoint(const Trace &trace, std::vector< DataPoint >::const_iterator begin)
Move right to the lower value.
Definition: trace.cpp:205