MNE-CPP  beta 1.0
atom.cpp
Go to the documentation of this file.
1 //=============================================================================================================
42 //*************************************************************************************************************
43 //=============================================================================================================
44 // INCLUDES
45 //=============================================================================================================
46 
47 #include "atom.h"
48 
49 //*************************************************************************************************************
50 //=============================================================================================================
51 // STL INCLUDES
52 //=============================================================================================================
53 
54 #include <iostream>
55 #include <algorithm> // std::sort
56 #include <vector> // std::vector
57 
58 //*************************************************************************************************************
59 //=============================================================================================================
60 // QT INCLUDES
61 //=============================================================================================================
62 
63 #include <QFile>
64 #include <QDebug>
65 
66 
67 //*************************************************************************************************************
68 //=============================================================================================================
69 // USED NAMESPACES
70 //=============================================================================================================
71 
72 using namespace UTILSLIB;
73 
74 //*************************************************************************************************************
75 //=============================================================================================================
76 // DEFINE MEMBER METHODS
77 //=============================================================================================================
78 
79 //*************************************************************************************************************
80 
81 FixDictAtom::FixDictAtom(qint32 _id, qint32 _sample_count, QString _dict_source)
82 {
83  id = _id;
84  sample_count = _sample_count;
85  dict_source = _dict_source;
86  energy = 0;
87  max_scalar_product = 0;
88 }
89 
90 //*************************************************************************************************************
91 
93 {
94  energy = 0;
95  max_scalar_product = 0;
96 }
97 
98 //*************************************************************************************************************
99 
101 
102 //*************************************************************************************************************
103 
105 {
106  energy = 0;
107  max_scalar_product = 0;
108 }
109 
110 //*************************************************************************************************************
111 
113 
114 //*************************************************************************************************************
115 
116 ChirpAtom::ChirpAtom() {}
117 
118 //*************************************************************************************************************
119 
120 ChirpAtom::~ChirpAtom() {}
121 
122 //************************************************************************************************************
123 
124 QStringList GaborAtom::create_string_values(qint32 sample_count, qreal scale, qint32 translation, qreal modulation, qreal phase)
125 {
126  VectorXd atomValues = create_real(sample_count, scale, translation, modulation, phase);
127  QStringList atomStringValues;
128  for(qint32 i = 0; i < atomValues.rows(); i++)
129  atomStringValues.append(QString("%1").arg(atomValues[i]));
130  return atomStringValues;
131 }
132 
133 //*************************************************************************************************************
134 
135 VectorXd GaborAtom::gauss_function (qint32 sample_count, qreal scale, quint32 translation)
136 {
137  VectorXd gauss = VectorXd::Zero(sample_count);
138 
139  for(qint32 n = 0; n < sample_count; n++)
140  {
141  qreal t = (qreal(n) - translation) / scale;
142  gauss[n] = exp(-PI * pow(t, 2))*pow(sqrt(scale),(-1))*pow(qreal(2),(0.25));
143  }
144 
145  return gauss;
146 }
147 
148 //*************************************************************************************************************
149 
150 VectorXcd GaborAtom::create_complex(qint32 sample_count, qreal scale, quint32 translation, qreal modulation)
151 {
152  VectorXcd complex_atom(sample_count);
153  qreal norm_atom = 0;
154 
155  //if scale == signalLength and translation == middle of signal there is no window or envelope necessary
156  //thats a simplification to save calculation time and reduces the windowing effects
157  if(scale == sample_count && translation == floor(sample_count / 2))
158  for(qint32 i = 0; i < sample_count; i++)
159  complex_atom[i] = std::polar(1 / sqrt(qreal(sample_count)), 2 * PI * modulation / qreal(sample_count) * qreal(i));
160 
161  //else scale is smaler than signalLength and translation is not in the middle an envelopement is required
162  else
163  for(qint32 i = 0; i < sample_count; i++)
164  {
165  qreal t = (qreal(i) - qreal(translation)) / qreal(scale);
166  complex_atom[i] = std::polar(1 / sqrt(sample_count) * pow(2.0, 0.25) * exp( -PI * pow(t, 2.0))
167  , 2 * PI * modulation / qreal(sample_count) * qreal(i));
168  }
169 
170  //normalization
171  for(qint32 i = 0; i < sample_count; i++)
172  norm_atom += pow(abs(complex_atom[i]), 2.0);
173 
174  norm_atom = sqrt(norm_atom);
175 
176  if(norm_atom != 0)
177  for(qint32 i = 0; i < sample_count; i++)
178  complex_atom[i] = complex_atom[i] / norm_atom;
179 
180  return complex_atom;
181 }
182 
183 //*************************************************************************************************************
184 
185 VectorXd GaborAtom::create_real(qint32 sample_count, qreal scale, quint32 translation, qreal modulation, qreal phase)
186 {
187  VectorXd real_atom(sample_count);
188  qreal norm = 0;
189 
190  if(scale == sample_count)
191  for(qint32 i = 0; i < sample_count; i++)
192  real_atom[i] = 1 / sqrt(sample_count) * cos( 2 * PI * modulation / sample_count * qreal(i) + phase);
193  else
194  {
195  VectorXd envelope = GaborAtom::gauss_function(sample_count, scale, translation);
196  for(qint32 i = 0; i < sample_count; i++)
197  real_atom[i] = envelope[i] * cos(2 * PI * modulation / sample_count * qreal(i) + phase);
198  }
199 
200  //normalization
201  norm = real_atom.norm();
202  if(norm != 0) real_atom /= norm;
203 
204  return real_atom; //length of the vector realAtom is 1 after normalization
205 }
206 
207 //*************************************************************************************************************
208 
209 VectorXd ChirpAtom::gauss_function (qint32 sample_count, qreal scale, quint32 translation)
210 {
211  VectorXd gauss = VectorXd::Zero(sample_count);
212 
213  for(qint32 n = 0; n < sample_count; n++)
214  {
215  qreal t = (qreal(n) - translation) / scale;
216  gauss[n] = exp(-PI * pow(t, 2)) * pow(sqrt(scale), (-1)) * pow(qreal(2),(0.25));
217  }
218 
219  return gauss;
220 }
221 
222 //*************************************************************************************************************
223 
224 VectorXd ChirpAtom::create_real(qint32 sample_count, qreal scale, quint32 translation, qreal modulation, qreal phase, qreal chirp)
225 {
226  VectorXd real_atom(sample_count);
227  qreal norm = 0;
228 
229 
230  if(scale == sample_count)
231  for(qint32 i = 0; i < sample_count; i++)
232  real_atom[i] = 1 / sqrt(sample_count) * cos( 2 * PI * modulation / sample_count * qreal(i) +(chirp / (sample_count * qreal(2))) * pow(qreal(i) - qreal(translation), 2) + phase);
233  else
234  {
235  VectorXd envelope = GaborAtom::gauss_function(sample_count, scale, translation);
236  for(qint32 i = 0; i < sample_count; i++)
237  real_atom[i] = envelope[i] * cos(2 * PI * modulation / sample_count * qreal(i) + (chirp / (sample_count * qreal(2))) * pow(qreal(i) - qreal(translation), 2) + phase);
238  }
239 
240  //normalization
241  norm = real_atom.norm();
242  if(norm != 0) real_atom /= norm;
243 
244  return real_atom; //length of the vector realAtom is 1 after normalization
245 }
246 
247 //*************************************************************************************************************
248 
249 QStringList ChirpAtom::create_string_values(qint32 sample_count, qreal scale, quint32 translation, qreal modulation, qreal phase, qreal chirp)
250 {
251  VectorXd atomValues = create_real(sample_count, scale, translation, modulation, phase, chirp);
252  QStringList atomStringValues;
253  for(qint32 i = 0; i < atomValues.rows(); i++)
254  atomStringValues.append(QString("%1").arg(atomValues[i]));
255  return atomStringValues;
256 }
257 
258 //*************************************************************************************************************
259 
260 
static VectorXd gauss_function(qint32 sample_count, qreal scale, quint32 translation)
Definition: atom.cpp:135
VectorXcd create_complex(qint32 sample_count, qreal scale, quint32 translation, qreal modulation)
Definition: atom.cpp:150
VectorXd create_real(qint32 sample_count, qreal scale, quint32 translation, qreal modulation, qreal phase)
Definition: atom.cpp:185
QStringList create_string_values(qint32 sample_count, qreal scale, qint32 translation, qreal modulation, qreal phase)
Definition: atom.cpp:124
ATOM class declaration, providing core features and parameters of Atoms used in Matching Pursiut Algo...