MNE-CPP  beta 1.0
fixdictmp.cpp
Go to the documentation of this file.
1 //=============================================================================================================
40 //*************************************************************************************************************
41 //=============================================================================================================
42 // INCLUDES
43 //=============================================================================================================
44 
45 #include <vector>
46 #include "fixdictmp.h"
47 
48 //*************************************************************************************************************
49 //=============================================================================================================
50 // USED NAMESPACES
51 //=============================================================================================================
52 
53 using namespace UTILSLIB;
54 
55 //*************************************************************************************************************
56 //=============================================================================================================
57 // DEFINE MEMBER METHODS
58 //=============================================================================================================
59 
61 : it(0)
62 , signal_energy(0)
63 , current_energy(0){}
64 
66 
68 
70 
71 //*************************************************************************************************************
72 
73 void FixDictMp::matching_pursuit(MatrixXd signal, qint32 max_iterations, qreal epsilon, qint32 boost, QString path, qreal delta)
74 {
75  Eigen::initParallel();
76 
77  std::cout << "\nFixDict Matching Pursuit Algorithm started...\n";
78 
79  QList<Dictionary> parsed_dicts;
80 
81  qint32 sample_count = signal.rows();
82  qint32 channel_count = signal.cols();
83 
84  signal_energy = 0;
85  qreal residuum_energy = 0;
86  qreal energy_threshold = 0;
87  qreal last_energy = 0;
88  bool sample_count_mismatch = false;
89 
90  this->residuum = signal;
91  parsed_dicts = parse_xml_dict(path);
92 
93  //calculate signal_energy
94  for(qint32 channel = 0; channel < channel_count; channel++)
95  {
96  for(qint32 sample = 0; sample < sample_count; sample++)
97  signal_energy += (signal(sample, channel) * signal(sample, channel));
98 
99  energy_threshold = 0.01 * epsilon * signal_energy;
100  residuum_energy = signal_energy;
101  }
102 
103  std::cout << "absolute energy of signal: " << residuum_energy << "\n";
104 
105  while(it < max_iterations && energy_threshold < residuum_energy)
106  {
107  FixDictAtom global_best_matching;
108 
109  QList<find_best_matching> list_of_best;
110  for(qint32 i = 0; i < parsed_dicts.length(); i++)
111  {
112  find_best_matching current_best_matching;
113  current_best_matching.pdict = parsed_dicts.at(i);
114  current_best_matching.current_resid = this->residuum;
115  current_best_matching.boost = boost;
116  list_of_best.append(current_best_matching);
117  }
118 
119  QFuture<FixDictAtom> mapped_best_matchings = QtConcurrent::mapped(list_of_best, &find_best_matching::parallel_correlation);// parse_threads;
120  mapped_best_matchings.waitForFinished();
121 
122  QFuture<FixDictAtom>::const_iterator i;
123  for (i = mapped_best_matchings.constBegin(); i != mapped_best_matchings.constEnd(); i++)
124  {
125  //FixDictAtom current_best_matching = correlation(parsed_dicts.at(i));
126  if(i == mapped_best_matchings.constBegin())
127  global_best_matching = *i;
128 
129  else if(abs(i->max_scalar_product) > abs(global_best_matching.max_scalar_product))
130  global_best_matching = *i;
131  }
132 
133  global_best_matching.display_text = create_display_text(global_best_matching);
134 
135  VectorXd fitted_atom = VectorXd::Zero(this->residuum.rows());
136  VectorXd resized_atom = VectorXd::Zero(this->residuum.rows());
137 
138  if(global_best_matching.atom_samples.rows() > this->residuum.rows())
139  for(qint32 k = 0; k < this->residuum.rows(); k++)
140  resized_atom[k] = global_best_matching.atom_samples[k + floor(global_best_matching.atom_samples.rows() / 2) - floor(this->residuum.rows() / 2)];
141  else resized_atom = global_best_matching.atom_samples;
142 
143 
144  for(qint32 k = 0; k < resized_atom.rows(); k++)
145  if((k + global_best_matching.translation - floor(resized_atom.rows() / 2)) >= 0
146  && (k + global_best_matching.translation - floor(resized_atom.rows() / 2)) < fitted_atom.rows())
147  fitted_atom[(k + global_best_matching.translation - floor(resized_atom.rows() / 2))] += resized_atom[k];
148 
149  //normalization
150  qreal norm = 0;
151  norm = fitted_atom.norm();
152  if(norm != 0) fitted_atom /= norm;
153 
154  for(qint32 chn = 0; chn < this->residuum.cols(); chn++)
155  {
156  qreal scalarproduct = 0;
157  for(qint32 sample = 0; sample < this->residuum.rows(); sample++)
158  scalarproduct += (this->residuum(sample, chn) * fitted_atom[sample]);
159 
160  global_best_matching.max_scalar_list.append(scalarproduct);//residuum(global_best_matching.translation, chn) / fitted_atom[global_best_matching.translation]);
161 
162  for(qint32 k = 0; k < fitted_atom.rows(); k++)
163  {
164  this->residuum(k,chn) -= global_best_matching.max_scalar_list.at(chn) * fitted_atom[k];
165  global_best_matching.energy += pow(global_best_matching.max_scalar_list.at(chn) * fitted_atom[k], 2); // * global_best_matching.max_scalar_list.at(chn) * fitted_atom[k];
166  }
167  }
168 
169  global_best_matching.atom_samples = fitted_atom;
170 
171 
172  last_energy = residuum_energy;
173  /*residuum_energy = 0;
174  for(qint32 channel = 0; channel < residuum.cols(); channel++)
175  {
176  for(qint32 sample = 0; sample < residuum.rows(); sample++)
177  residuum_energy += pow(residuum(sample, channel), 2);
178  }*/
179 
180  residuum_energy -= global_best_matching.energy;
181  current_energy += global_best_matching.energy;
182 
183  fix_dict_list.append(global_best_matching);
184 
185  if(global_best_matching.sample_count != signal.rows() && !sample_count_mismatch)
186  {
187  std::cout << "\n=============================================\n"
188  << " INFO\n\n(part-)dictionary does not fit the signal!\n\nResults might be better by creating dictionaries\ncontaining atoms of the same length as the signal\n"
189  << "\nsignal samples: " << signal.rows() << "\n"
190  << "atom samples: " << global_best_matching.sample_count << "\n"
191  << "=============================================\n";
192  sample_count_mismatch = true;
193  }
194 
195  std::cout << "\n" << "===============" << it + 1 <<"th atom found" << "===============" << ":\n\n"<<
196  qPrintable(global_best_matching.display_text) << "\n\n" << "sample_count: " << global_best_matching.sample_count <<
197  " source_dict: " << qPrintable(global_best_matching.dict_source) << " Atom ID: " << global_best_matching.id <<
198  "\nsclr_prdct: " << global_best_matching.max_scalar_product << "\n" <<
199  "\nabsolute energy of residue: " << residuum_energy << "\n";
200 
201  it++;
202 
203  emit current_result(it, max_iterations, current_energy, signal_energy, residuum, adaptive_list, fix_dict_list);
204 
205 
206  if(((last_energy * 100 / signal_energy) - (residuum_energy * 100 / signal_energy)) < delta)
207  {
208  std::cout << "\n=============================================\n"
209  << " ALGORITHM ABORTED\n\ndictionary excludes atoms to reduce further residual energy\n"
210  << "=============================================\n";
211  emit send_warning(1);
212  break;
213  }
214  if( QThread::currentThread()->isInterruptionRequested())
215  {
216  send_warning(10);
217  break;
218  }
219 
220  }//end while iterations
221 
222  std::cout << "\nFixDict Matching Pursuit Algorithm finished.\n";
223  emit finished_calc();
224  //return atom_list;
225 }
226 
227 //*************************************************************************************************************
228 
229 // calc scalarproduct of Atom and Signal
230 FixDictAtom FixDictMp::correlation(Dictionary current_pdict, MatrixXd current_resid, qint32 boost)
231 {
232  qint32 channel_count = current_resid.cols() * (boost / 100.0); //reducing the number of observed channels in the algorithm to increase speed performance
233  if(boost == 0 || channel_count == 0)
234  channel_count = 1;
235 
236  Eigen::FFT<double> fft;
237  std::ptrdiff_t max_index;
238  VectorXcd fft_atom = VectorXcd::Zero(current_resid.rows());
239 
240  FixDictAtom best_matching;
241  qreal max_scalar_product = 0;
242 
243  for(qint32 i = 0; i < current_pdict.atoms.length(); i++)
244  {
245  VectorXd fitted_atom = VectorXd::Zero(current_resid.rows());
246  qint32 p = floor(current_resid.rows() / 2);//translation
247 
248  VectorXd resized_atom = VectorXd::Zero(current_resid.rows());
249 
250  if(current_pdict.atoms.at(i).atom_samples.rows() > current_resid.rows())
251  for(qint32 k = 0; k < current_resid.rows(); k++)
252  resized_atom[k] = current_pdict.atoms.at(i).atom_samples[k + floor(current_pdict.atoms.at(i).atom_samples.rows() / 2) - floor(current_resid.rows() / 2)];
253  else resized_atom = current_pdict.atoms.at(i).atom_samples;
254 
255  if(resized_atom.rows() < current_resid.rows())
256  for(qint32 k = 0; k < resized_atom.rows(); k++)
257  fitted_atom[(k + p - floor(resized_atom.rows() / 2))] = resized_atom[k];
258  else fitted_atom = resized_atom;
259 
260  //normalization
261  qreal norm = 0;
262  norm = fitted_atom.norm();
263  if(norm != 0) fitted_atom /= norm;
264 
265  fft.fwd(fft_atom, fitted_atom);
266 
267  for(qint32 chn = 0; chn < channel_count; chn++)
268  {
269  p = floor(current_resid.rows() / 2);//translation
270  VectorXd corr_coeffs = VectorXd::Zero(current_resid.rows());
271  VectorXcd fft_signal = VectorXcd::Zero(current_resid.rows());
272  VectorXcd fft_sig_atom = VectorXcd::Zero(current_resid.rows());
273 
274  fft.fwd(fft_signal, current_resid.col(chn));
275 
276  for( qint32 m = 0; m < current_resid.rows(); m++)
277  fft_sig_atom[m] = fft_signal[m] * conj(fft_atom[m]);
278 
279  fft.inv(corr_coeffs, fft_sig_atom);
280 
281  //find index of maximum correlation-coefficient to use in translation
282  max_scalar_product = corr_coeffs.maxCoeff(&max_index);
283 
284  if(i == 0)
285  {
286  best_matching = current_pdict.atoms.at(i);
287  best_matching.max_scalar_product = max_scalar_product;
288 
289  //adapting translation p to create atomtranslation correctly
290  if(max_index >= p && current_resid.rows() % (2) == 0) p = max_index - p;
291  else if(max_index >= p && current_resid.rows() % (2) != 0) p = max_index - p - 1;
292  else p = max_index + p;
293 
294  best_matching.translation = p;
295  }
296 
297  else if(abs(max_scalar_product) > abs(best_matching.max_scalar_product))
298  {
299  best_matching = current_pdict.atoms.at(i);
300  best_matching.max_scalar_product = max_scalar_product;
301 
302  //adapting translation p to create atomtranslation correctly
303  if(max_index >= p && current_resid.rows() % (2) == 0) p = max_index - p;
304  else if(max_index >= p && current_resid.rows() % (2) != 0) p = max_index - p - 1;
305  else p = max_index + p;
306 
307  best_matching.translation = p;
308  }
309  }
310  }
311  best_matching.atom_formula = current_pdict.atom_formula;
312  best_matching.dict_source = current_pdict.source;
313  best_matching.type = current_pdict.type;
314  best_matching.sample_count = current_pdict.sample_count;
315 
316  return best_matching;
317 }
318 //*************************************************************************************************************
319 
320 QList<Dictionary> FixDictMp::parse_xml_dict(QString path)
321 {
322  QFile current_dict(path);
323  std::cout << "\nparsing dictionary, please be patient...";
324 
325  QDomDocument dictionary;
326  dictionary.setContent(&current_dict);
327  QList<Dictionary> parsed_dict;
328  QDomElement current_element = dictionary.documentElement();
329  QDomNodeList node_list = current_element.childNodes();
330  QList<parse_node> pdict_nodes;
331  QList<QDomNode> nodes_listed;
332  bool is_emitted = false;
333 
334  for(qint32 i = 0; i < node_list.length(); i++)
335  {
336  parse_node temp_parse;
337  temp_parse.node = node_list.at(i);
338  pdict_nodes.append(temp_parse);
339  nodes_listed.append(node_list.at(i));
340  }
341 
342  //multithreading
343  QFuture<Dictionary> mapped_dicts = QtConcurrent::mapped(pdict_nodes, &parse_node::fill_dict_in_map);// parse_threads;
344  mapped_dicts.waitForFinished();
345 
346  QFuture<Dictionary>::const_iterator i;
347 
348  for (i = mapped_dicts.constBegin(); i != mapped_dicts.constEnd(); i++)
349  {
350  //find_best_matching current_best_matching;
351  //current_best_matching.pdict = *i;
352  //current_best_matching.current_resid = this->residuum;
353  if(i->sample_count != this->residuum.rows() && !is_emitted)
354  {
355  is_emitted = true;
356  emit send_warning(2);
357  }
358 
359  parsed_dict.append(*i);
360  }
361 
362  std::cout << " done.\n\n";
363 
364  return parsed_dict;
365 }
366 //*************************************************************************************************************
367 
368 Dictionary FixDictMp::fill_dict(const QDomNode &pdict)
369 {
370  Dictionary current_dict;
371  QDomElement atom = pdict.firstChildElement("ATOM");
372 
373  current_dict.source = pdict.toElement().attribute("source_dict");
374  current_dict.atom_count();// = pdict.toElement().attribute("atom_count").toInt();
375  current_dict.atom_formula = pdict.toElement().attribute("formula");
376  current_dict.sample_count = pdict.toElement().attribute("sample_count").toInt();
377 
378  if(current_dict.atom_formula == QString("Gaboratom"))
379  {
380  current_dict.type = GABORATOM;
381 
382  while(!atom.isNull())
383  {
384  FixDictAtom current_atom;
385 
386  current_atom.id = atom.attribute("ID").toInt();
387  current_atom.gabor_atom.scale = atom.attribute("scale").toDouble();
388  current_atom.gabor_atom.modulation = atom.attribute("modu").toDouble();
389  current_atom.gabor_atom.phase = atom.attribute("phase").toDouble();
390 
391  if(atom.hasChildNodes())
392  {
393  QString sample_string = atom.firstChild().toElement().attribute("samples");
394  QStringList sample_list = sample_string.split(":", QString::SkipEmptyParts);
395  current_atom.atom_samples = VectorXd::Zero(sample_list.length());
396  for(qint32 i = 0; i < sample_list.length(); i++)
397  current_atom.atom_samples[i] = sample_list.at(i).toDouble();
398  }
399  current_dict.atoms.append(current_atom);
400  atom = atom.nextSiblingElement("ATOM");
401  }
402  }
403  else if(current_dict.atom_formula == QString("Chirpatom"))
404  {
405  current_dict.type = CHIRPATOM;
406 
407  while(!atom.isNull())
408  {
409  FixDictAtom current_atom;
410 
411  current_atom.id = atom.attribute("id").toInt();
412  current_atom.chirp_atom.scale = atom.attribute("scale").toDouble();
413  current_atom.chirp_atom.modulation = atom.attribute("modu").toDouble();
414  current_atom.chirp_atom.phase = atom.attribute("phase").toDouble();
415  current_atom.chirp_atom.chirp = atom.attribute("chirp").toDouble();
416 
417  if(atom.hasChildNodes())
418  {
419  QString sample_string = atom.firstChild().toElement().attribute("samples");
420  QStringList sample_list = sample_string.split(":", QString::SkipEmptyParts);
421  current_atom.atom_samples = VectorXd::Zero(sample_list.length());
422  for(qint32 i = 0; i < sample_list.length(); i++)
423  current_atom.atom_samples[i] = sample_list.at(i).toDouble();
424  }
425  current_dict.atoms.append(current_atom);
426  atom = atom.nextSiblingElement("ATOM");
427  }
428  }
429  else
430  {
431  current_dict.type = FORMULAATOM;
432 
433  while(!atom.isNull())
434  {
435  FixDictAtom current_atom;
436 
437  current_atom.id = atom.attribute("id").toInt();
438  current_atom.formula_atom.a = atom.attribute("a").toDouble();
439  current_atom.formula_atom.b = atom.attribute("b").toDouble();
440  current_atom.formula_atom.c = atom.attribute("c").toDouble();
441  current_atom.formula_atom.d = atom.attribute("d").toDouble();
442  current_atom.formula_atom.e = atom.attribute("e").toDouble();
443  current_atom.formula_atom.f = atom.attribute("f").toDouble();
444  current_atom.formula_atom.g = atom.attribute("g").toDouble();
445  current_atom.formula_atom.h = atom.attribute("h").toDouble();
446 
447  if(atom.hasChildNodes())
448  {
449  QString sample_string = atom.firstChild().toElement().attribute("samples");
450  QStringList sample_list = sample_string.split(":", QString::SkipEmptyParts);
451  current_atom.atom_samples = VectorXd::Zero(sample_list.length());
452  for(qint32 i = 0; i < sample_list.length(); i++)
453  current_atom.atom_samples[i] = sample_list.at(i).toDouble();
454  }
455  current_dict.atoms.append(current_atom);
456  atom = atom.nextSiblingElement("ATOM");
457  }
458  }
459 
460  return current_dict;
461 }
462 
463 //*************************************************************************************************************
464 
465 QString FixDictMp::create_display_text(FixDictAtom global_best_matching)
466 {
467  QString display_text = "";
468  if(global_best_matching.type == GABORATOM)
469  {
470  display_text = QString("Gaboratom: scale: %0, translation: %1, modulation: %2, phase: %3")
471  .arg(QString::number(global_best_matching.gabor_atom.scale, 'f', 2))
472  .arg(QString::number(global_best_matching.translation, 'f', 2))
473  .arg(QString::number(global_best_matching.gabor_atom.modulation, 'f', 2))
474  .arg(QString::number(global_best_matching.gabor_atom.phase, 'f', 2));
475  }
476  else if(global_best_matching.type == CHIRPATOM)
477  {
478  display_text = QString("Chripatom: scale: %0, translation: %1, modulation: %2, phase: %3, chirp: %4")
479  .arg(QString::number(global_best_matching.chirp_atom.scale, 'f', 2))
480  .arg(QString::number(global_best_matching.translation, 'f', 2))
481  .arg(QString::number(global_best_matching.chirp_atom.modulation, 'f', 2))
482  .arg(QString::number(global_best_matching.chirp_atom.phase, 'f', 2))
483  .arg(QString::number(global_best_matching.chirp_atom.chirp, 'f', 2));
484  }
485  else if(global_best_matching.type == FORMULAATOM)
486  {
487  display_text = QString("%0: transl: %1 a: %2, b: %3 c: %4, d: %5, e: %6, f: %7, g: %8, h: %9")
488  .arg(global_best_matching.atom_formula)
489  .arg(QString::number(global_best_matching.translation, 'f', 2))
490  .arg(QString::number(global_best_matching.formula_atom.a, 'f', 2))
491  .arg(QString::number(global_best_matching.formula_atom.b, 'f', 2))
492  .arg(QString::number(global_best_matching.formula_atom.c, 'f', 2))
493  .arg(QString::number(global_best_matching.formula_atom.d, 'f', 2))
494  .arg(QString::number(global_best_matching.formula_atom.e, 'f', 2))
495  .arg(QString::number(global_best_matching.formula_atom.f, 'f', 2))
496  .arg(QString::number(global_best_matching.formula_atom.g, 'f', 2))
497  .arg(QString::number(global_best_matching.formula_atom.h, 'f', 2));
498  }
499 
500  return display_text;
501 }
502 
503 //*************************************************************************************************************
504 
505 void FixDictMp::recieve_input(MatrixXd signal, qint32 max_iterations, qreal epsilon, qint32 boost, QString path, qreal delta)
506 {
507  matching_pursuit(signal, max_iterations, epsilon, boost, path, delta);
508 }
509 
510 //*************************************************************************************************************
511 
512 qint32 Dictionary::atom_count()
513 {
514  return atoms.length();
515 }
516 
517 //*************************************************************************************************************
518 
519  void Dictionary::clear()
520  {
521  this->atoms.clear();
522  this->atom_formula = "";
523  this->sample_count = 0;
524  this->source = "";
525  }
526 
527  //*************************************************************************************************************
528 
529 /*
530 void FixDictMp::create_tree_dict(QString save_path)
531 {
532  QFile file(save_path);
533  QDomDocument atom_xml_file;
534  atom_xml_file.setContent(&file);
535  GaborAtom* gabor_Atom = new GaborAtom;
536 
537  qint32 sample_count = 0;
538  qreal scale = 0;
539  quint32 translation = 0;
540  qreal modulation = 0;
541  qreal phase = 0;
542 
543  QDomElement xml_element = atom_xml_file.documentElement();
544 
545  QString root_tag = xml_element.tagName();
546 
547  if(root_tag == "builtAtomsTreebasedMP")
548  sample_count = xml_element.attribute("sample_count").toInt();
549 
550  VectorXd compare_atom = VectorXd::Zero(sample_count);
551 
552  QString write_save_path = save_path;
553  write_save_path.replace(".tbd", ".dict");
554  QFile write_file(write_save_path);
555  write_file.open(QIODevice::WriteOnly);
556 
557  QXmlStreamWriter write_molecules_to_xml(&write_file);
558  write_molecules_to_xml.setAutoFormatting(true);
559 
560  write_molecules_to_xml.writeStartDocument();
561  write_molecules_to_xml.writeStartElement("TreebasedStructureFor_TBMP");
562  write_molecules_to_xml.writeAttribute("sample_count", QString::number(sample_count));
563 
564  QDomNodeList node_list = xml_element.elementsByTagName("Atom");
565  qint32 max_it = node_list.count();
566  qint32 number_of_molecs = 0;//only debugout
567 
568  for(qint32 i = 0; i < max_it; i++)
569  {
570  QDomElement current_element = node_list.at(0).toElement();
571 
572  if(!current_element.isNull())
573  {
574  scale = (current_element.attribute("scale", current_element.text())).toDouble();
575  translation = (current_element.attribute("translation", current_element.text())).toInt();
576  modulation = (current_element.attribute("modulation", current_element.text())).toDouble();
577  phase = (current_element.attribute("phase", current_element.text())).toDouble();
578  }
579 
580  compare_atom = gabor_Atom->create_real(sample_count,scale,translation,modulation,phase);
581  qreal threshold = 0.8 * compare_atom.dot(compare_atom); //vary this to find optimum of seperation
582 
583  QList<qint32> similar_atoms;
584  VectorXd temp_atom = VectorXd::Zero(sample_count);
585 
586  //finding atoms with low differences
587  for (qint32 next = 0; next < node_list.count(); next++)//;32)
588  {
589  //try the next atom
590  current_element = node_list.at(next).toElement();
591 
592  if(!current_element.isNull())
593  {
594  scale = (current_element.attribute("scale", current_element.text())).toDouble();
595  translation = (current_element.attribute("translation", current_element.text())).toInt();
596  modulation = (current_element.attribute("modulation", current_element.text())).toDouble();
597  phase = (current_element.attribute("phase", current_element.text())).toDouble();
598  }
599 
600  //set atoms out of xml to compare with compare_atom
601  temp_atom = gabor_Atom->create_real(sample_count,scale,translation,modulation,phase);
602 
603  //fill list of similar atoms to save as molecule
604  if( compare_atom.dot(temp_atom) > threshold)
605  similar_atoms.append(current_element.attribute("ID", current_element.text()).toInt());
606 
607  //if (similar_atoms.length() == 32)//vary this for testing
608  // break;
609  }
610 
611  qreal molec_scale = 0;
612  quint32 molec_translation = 0;
613  qreal molec_modulation = 0;
614  qreal molec_phase = 0;
615 
616  //calc molecule params
617 
618  for(qint32 j = 0; j < similar_atoms.length(); j++)
619  {
620  current_element = node_list.at(0).toElement();
621 
622  while(!current_element.isNull())
623  {
624  if(current_element.attribute("ID", current_element.text()).toInt() == similar_atoms.at(j))
625  {
626  molec_scale += (current_element.attribute("scale", current_element.text())).toDouble();
627  molec_translation += (current_element.attribute("translation", current_element.text())).toInt();
628  molec_modulation += (current_element.attribute("modulation", current_element.text())).toDouble();
629  molec_phase += (current_element.attribute("phase", current_element.text())).toDouble();
630  }
631  current_element = current_element.nextSibling().toElement();
632  }
633  }
634 
635  molec_scale /= similar_atoms.length();
636  molec_translation /= similar_atoms.length();
637  molec_modulation /= similar_atoms.length();
638  molec_phase /= similar_atoms.length();
639 
640  write_molecules_to_xml.writeStartElement("Molecule");
641  write_molecules_to_xml.writeAttribute("level", "0");
642  write_molecules_to_xml.writeAttribute("ID", QString::number(i));
643  write_molecules_to_xml.writeAttribute("scale", QString::number(molec_scale));
644  write_molecules_to_xml.writeAttribute("translation", QString::number(molec_translation));
645  write_molecules_to_xml.writeAttribute("modulation", QString::number(molec_modulation));
646  write_molecules_to_xml.writeAttribute("phase", QString::number(molec_phase));
647 
648  // write atoms to molecule
649  for(qint32 k = 0; k < similar_atoms.length(); k++)
650  {
651  current_element = node_list.at(0).toElement();
652 
653  while(!current_element.isNull())
654  {
655  if(current_element.attribute("ID", current_element.text()).toInt() == similar_atoms.at(k))
656  {
657  write_molecules_to_xml.writeStartElement("Atom");
658  write_molecules_to_xml.writeAttribute("ID", current_element.attribute("ID", current_element.text()));
659  write_molecules_to_xml.writeAttribute("scale", current_element.attribute("scale", current_element.text()));
660  write_molecules_to_xml.writeAttribute("translation", current_element.attribute("translation", current_element.text()));
661  write_molecules_to_xml.writeAttribute("modulation", current_element.attribute("modulation", current_element.text()));
662  write_molecules_to_xml.writeAttribute("phase", current_element.attribute("phase", current_element.text()));
663 
664  write_molecules_to_xml.writeEndElement();//atom
665  }
666  current_element = current_element.nextSibling().toElement();
667  }
668  }
669  write_molecules_to_xml.writeEndElement();//molecule
670 
671  //remove found atoms from node list
672  for(qint32 n = 0; n < similar_atoms.length(); n++)
673  {
674  current_element = node_list.at(0).toElement();
675  QDomNode to_remove = current_element.parentNode();
676 
677  while(!current_element.isNull())
678  {
679  if(current_element.attribute("ID", current_element.text()).toInt() == similar_atoms.at(n))
680  to_remove.removeChild(current_element);
681 
682  current_element = current_element.nextSibling().toElement();
683  }
684  }
685 
686  if(node_list.count() == 0)
687  break;
688 
689  number_of_molecs++;
690  }//for max_it of nodelist
691 
692  write_molecules_to_xml.writeEndElement();//header
693  write_molecules_to_xml.writeEndDocument();
694 
695  cout << "...start to build tree recursive... number of built molecules: " << number_of_molecs + 1 << "\n";
696 
697  //close reader and writer
698  file.close();
699  write_file.close();
700  delete gabor_Atom;
701 
702  //ToDo: recursion
703  build_molecule_xml_file(0);
704 
705  cout << "finished treebuilding\n";
706 
707 
708 }
709 
710 //******************************************************************************************************************
711 
712 void FixDictMp::build_molecule_xml_file(qint32 level_counter)
713 {
714  //qint32 i = 0;
715  if(level_counter < 10)
716  {
717  QString save_path = QString("Matching-Pursuit-Toolbox/my_first_16atoms.dict");//ToDo: make it flexible
718  QFile file(save_path);
719  QString temp_path("Matching-Pursuit-Toolbox/__temp.dict");
720  file.copy(temp_path);
721  QFile temp_file(temp_path);
722 
723  QDomDocument atom_xml_file;
724  atom_xml_file.setContent(&file);
725  QDomElement xml_element = atom_xml_file.documentElement();
726  xml_element = atom_xml_file.documentElement();
727 
728  GaborAtom* gabor_Atom = new GaborAtom;
729  qint32 sample_count = 0;
730  qreal scale = 0;
731  quint32 translation = 0;
732  qreal modulation = 0;
733  qreal phase = 0;
734  //qint32 level_counter = 0;
735 
736  QString root_tag = xml_element.tagName();
737  root_tag = xml_element.tagName();
738 
739  if(root_tag == "TreebasedStructureFor_TBMP")
740  sample_count = xml_element.attribute("sample_count").toInt();
741 
742  VectorXd compare_molec = VectorXd::Zero(sample_count);
743 
744  QDomNodeList node_list = xml_element.elementsByTagName("Molecule");
745 
746  qint32 max_it = node_list.count();
747  QDomElement current_element;
748 
749  temp_file.open(QIODevice::WriteOnly);
750 
751  QXmlStreamWriter write_molecules_to_xml(&temp_file);
752  write_molecules_to_xml.setAutoFormatting(true);
753 
754  write_molecules_to_xml.writeStartDocument();
755  write_molecules_to_xml.writeStartElement("TreebasedStructureFor_TBMP");
756  write_molecules_to_xml.writeAttribute("sample_count", QString::number(sample_count));
757 
758  for(qint32 i = 0; i < max_it; i++)
759  {
760  current_element = node_list.at(0).toElement();
761  QList<qint32> similar_molecs;
762 
763  if(current_element.attribute("level").toInt() == level_counter && current_element.nodeName() == "Molecule")
764  {
765  if(!current_element.isNull())
766  {
767  scale = (current_element.attribute("scale", current_element.text())).toDouble();
768  translation = (current_element.attribute("translation", current_element.text())).toInt();
769  modulation = (current_element.attribute("modulation", current_element.text())).toDouble();
770  phase = (current_element.attribute("phase", current_element.text())).toDouble();
771  }
772  compare_molec = gabor_Atom->create_real(sample_count,scale,translation,modulation,phase);
773  qreal threshold = 0.6 * compare_molec.dot(compare_molec); //vary this to find optimum of seperation
774 
775  VectorXd temp_molec = VectorXd::Zero(sample_count);
776 
777  //finding molecules with low differences
778  for (qint32 next = 0; next < node_list.count(); next++)//;32)
779  {
780  current_element = node_list.at(next).toElement();
781  //try the next molecule
782  if(current_element.attribute("level").toInt() == level_counter && current_element.nodeName() == "Molecule")
783  { //current_element = node_list.at(next).toElement();
784 
785  if(!current_element.isNull())
786  {
787  scale = (current_element.attribute("scale", current_element.text())).toDouble();
788  translation = (current_element.attribute("translation", current_element.text())).toInt();
789  modulation = (current_element.attribute("modulation", current_element.text())).toDouble();
790  phase = (current_element.attribute("phase", current_element.text())).toDouble();
791  }
792  //set molecules out of xml to compare with molecule
793  temp_molec = gabor_Atom->create_real(sample_count,scale,translation,modulation,phase);
794 
795  //fill list of similar molecules to save as new molecule
796  if( compare_molec.dot(temp_molec) > threshold)
797  similar_molecs.append(current_element.attribute("ID", current_element.text()).toInt());
798 
799  //if (similar_molecs.length() == 32)//vary this for testing
800  // break;
801  }
802  }
803 
804  qreal molec_scale = 0;
805  quint32 molec_translation = 0;
806  qreal molec_modulation = 0;
807  qreal molec_phase = 0;
808 
809  //calc molecule params
810  for(qint32 j = 0; j < similar_molecs.length(); j++)
811  {
812  current_element = node_list.at(0).toElement();
813 
814  while(!current_element.isNull())
815  {
816  if(current_element.attribute("ID", current_element.text()).toInt() == similar_molecs.at(j) && current_element.attribute("level").toInt() == level_counter)
817  {
818  molec_scale += (current_element.attribute("scale", current_element.text())).toDouble();
819  molec_translation += (current_element.attribute("translation", current_element.text())).toInt();
820  molec_modulation += (current_element.attribute("modulation", current_element.text())).toDouble();
821  molec_phase += (current_element.attribute("phase", current_element.text())).toDouble();
822  }
823  current_element = current_element.nextSiblingElement("Molecule").toElement();
824  }
825  }
826 
827  molec_scale /= similar_molecs.length();
828  molec_translation /= similar_molecs.length();
829  molec_modulation /= similar_molecs.length();
830  molec_phase /= similar_molecs.length();
831 
832  write_molecules_to_xml.writeStartElement("Molecule");
833  write_molecules_to_xml.writeAttribute("level", QString::number(level_counter + 1));
834  write_molecules_to_xml.writeAttribute("ID", QString::number(i));
835  write_molecules_to_xml.writeAttribute("scale", QString::number(molec_scale));
836  write_molecules_to_xml.writeAttribute("translation", QString::number(molec_translation));
837  write_molecules_to_xml.writeAttribute("modulation", QString::number(molec_modulation));
838  write_molecules_to_xml.writeAttribute("phase", QString::number(molec_phase));
839 
840  //treedepth building
841 
842  current_element = node_list.at(0).toElement();
843  qint32 root_counter = 0;
844 
845  for(qint32 k = 0; k < similar_molecs.length(); k++)
846  {
847  //current_element = node_list.at(current_element.attribute("ID", current_element.text()).toInt() == similar_molecs.at(k) && current_element.attribute("level").toInt() == level_counter).toElement();
848  current_element = node_list.at(0).toElement();
849 
850  while(!current_element.isNull())
851  {
852  root_counter = 0;
853 
854  if(current_element.attribute("ID").toInt() == similar_molecs.at(k) && current_element.attribute("level").toInt() == level_counter)
855  {
856  QDomElement temp_element = current_element;
857  qint32 end_element_counter = 0;
858  qint32 child_counter = 0;//ToDo: set right child for right position in tree to write all atoms and molecules
859  root_counter = 0;
860 
861  while(temp_element.hasChildNodes())
862  {
863  end_element_counter = 0;
864  child_counter = 0;
865 
866  while(!temp_element.isNull())
867  {
868  write_molecules_to_xml.writeStartElement(temp_element.nodeName());
869 
870  if(temp_element.nodeName() == "Molecule")
871  write_molecules_to_xml.writeAttribute("level", temp_element.attribute("level", temp_element.text()));
872 
873  write_molecules_to_xml.writeAttribute("ID", temp_element.attribute("ID", temp_element.text()));
874  write_molecules_to_xml.writeAttribute("scale", temp_element.attribute("scale", temp_element.text()));
875  write_molecules_to_xml.writeAttribute("translation", temp_element.attribute("translation", temp_element.text()));
876  write_molecules_to_xml.writeAttribute("modulation", temp_element.attribute("modulation", temp_element.text()));
877  write_molecules_to_xml.writeAttribute("phase", temp_element.attribute("phase", temp_element.text()));
878 
879  if(temp_element.nodeName() == "Atom")
880  write_molecules_to_xml.writeEndElement();//close atoms immediately
881  else
882  end_element_counter++;//keep in mind how many molecules are opened
883 
884  if(temp_element.hasChildNodes())//if there are children, that means the current element is not an atom: go to next child
885  {
886  temp_element = temp_element.firstChild().toElement();
887  child_counter++;//keep in mind how many children exist
888  }
889  else
890  temp_element = temp_element.nextSibling().toElement();//else there are no children, there must be an atom , so take the next atom
891 
892  }
893  root_counter++;
894  write_molecules_to_xml.writeEndElement();
895 
896  temp_element = current_element;
897 
898  //for(qint32 node_depth = 2; node_depth < child_counter; node_depth++)
899  // temp_element = temp_element.firstChildElement();
900 
901  for(qint32 node_depth = 2; node_depth < child_counter; node_depth++)
902  temp_element = temp_element.childNodes().at(root_counter).toElement();
903 
904  }
905  //for(qint32 close_root_molecs = end_element_counter; close_root_molecs > 1; close_root_molecs--)
906  for(qint32 close_root_molecs = level_counter; close_root_molecs > 0; close_root_molecs--)
907  write_molecules_to_xml.writeEndElement();//end elements for one submolecule in a "parentmolecule"
908 
909  }
910 
911  //for(qint32 close_parent_molecs = root_counter; close_parent_molecs > 1; close_parent_molecs--)
912  // write_molecules_to_xml.writeEndElement();//close all submolecules inside a new parentmolecule
913 
914  current_element = current_element.nextSiblingElement("Molecule").toElement();
915  }
916 
917  }//for all similar molecules
918 
919  //for(qint32 close_parent_molecs = root_counter; close_parent_molecs > 1; close_parent_molecs--)
920  write_molecules_to_xml.writeEndElement();
921 
922  //remove found molecules from node list
923  for(qint32 n = 0; n < similar_molecs.length(); n++)
924  {
925  current_element = node_list.at(0).toElement();
926  QDomNode to_remove = current_element.parentNode();
927 
928  while(!current_element.isNull())
929  {
930  if(current_element.attribute("ID").toInt() == similar_molecs.at(n) && current_element.attribute("level").toInt() == level_counter)
931  to_remove.removeChild(current_element);
932 
933  current_element = current_element.nextSiblingElement("Molecule").toElement();
934  }
935  }
936  }//if
937 
938  cout << node_list.count() <<"\n";
939 
940  if(node_list.count() == 0)
941  break;
942  }//for nodelist
943 
944  write_molecules_to_xml.writeEndElement();//header
945  write_molecules_to_xml.writeEndDocument();
946 
947  file.flush();
948  temp_file.flush();
949  file.close();
950  temp_file.close();
951 
952  file.remove();
953 
954  temp_file.copy(save_path);
955  temp_file.remove();
956 
957  level_counter++;
958 
959  build_molecule_xml_file(level_counter);
960  }
961 }
962 */
FixDictAtom used in fix dict MP Algorithm.
Definition: atom.h:123