bpp-seq-omics  2.4.1
MafParser.cpp
Go to the documentation of this file.
1 //
2 // File: MafParser.cpp
3 // Authors: Julien Dutheil
4 // Created: Tue Apr 27 2010
5 //
6 
7 /*
8 Copyright or © or Copr. Bio++ Development Team, (2010)
9 
10 This software is a computer program whose purpose is to provide classes
11 for sequences analysis.
12 
13 This software is governed by the CeCILL license under French law and
14 abiding by the rules of distribution of free software. You can use,
15 modify and/ or redistribute the software under the terms of the CeCILL
16 license as circulated by CEA, CNRS and INRIA at the following URL
17 "http://www.cecill.info".
18 
19 As a counterpart to the access to the source code and rights to copy,
20 modify and redistribute granted by the license, users are provided only
21 with a limited warranty and the software's author, the holder of the
22 economic rights, and the successive licensors have only limited
23 liability.
24 
25 In this respect, the user's attention is drawn to the risks associated
26 with loading, using, modifying and/or developing or reproducing the
27 software by the user in light of its specific status of free software,
28 that may mean that it is complicated to manipulate, and that also
29 therefore means that it is reserved for developers and experienced
30 professionals having in-depth computer knowledge. Users are therefore
31 encouraged to load and test the software's suitability as regards their
32 requirements in conditions enabling the security of their systems and/or
33 data to be ensured and, more generally, to use and operate it in the
34 same conditions as regards security.
35 
36 The fact that you are presently reading this means that you have had
37 knowledge of the CeCILL license and that you accept its terms.
38 */
39 
40 #include "MafParser.h"
43 #include <Bpp/Text/TextTools.h>
44 #include <Bpp/Text/KeyvalTools.h>
45 
46 #include <algorithm>
47 
48 using namespace std;
49 using namespace bpp;
50 
51 MafBlock* MafParser::analyseCurrentBlock_()
52 {
53  MafBlock* block = 0;
54 
55  string line;
56  bool test = true;
57  unique_ptr<MafSequence> currentSequence;
58 
59  while (test)
60  {
61  if (stream_->eof()) {
62  break;
63  }
64  getline(*stream_, line, '\n');
65  if (TextTools::isEmpty(line))
66  {
67  if (firstBlock_)
68  continue;
69  if (currentSequence) {
70  //Add previous sequence:
71  block->addSequence(*currentSequence); //The sequence is copied in the container.
72  currentSequence.reset();
73  }
74 
75  //end of paragraph
76  test = false;
77  }
78  else if (line[0] == 'a')
79  {
80  if (currentSequence) {
81  //Add previous sequence:
82  block->addSequence(*currentSequence); //The sequence is copied in the container.
83  currentSequence.reset();
84  }
85 
86  //New block.
87  block = new MafBlock();
88  firstBlock_ = false;
89 
90  map<string, string> args;
91  if (line.size() > 2)
92  {
93  KeyvalTools::multipleKeyvals(line.substr(2), args, " ");
94 
95  if (args.find("score") != args.end())
96  if (args["score"] != "NA")
97  block->setScore(TextTools::toDouble(args["score"]));
98 
99  if (args.find("pass") != args.end())
100  block->setPass(TextTools::to<unsigned int>(args["pass"]));
101  }
102  }
103  else if (line[0] == 's')
104  {
105  StringTokenizer st(line);
106  st.nextToken(); //The 's' tag
107  if (! st.hasMoreToken())
108  throw IOException("Sequence description should include a source field.");
109  string src = st.nextToken();
110  if (! st.hasMoreToken())
111  throw IOException("Sequence description should include a start field.");
112  unsigned int start = TextTools::to<unsigned int>(st.nextToken());
113  if (! st.hasMoreToken())
114  throw IOException("Sequence description should include a size field.");
115  unsigned int size = TextTools::to<unsigned int>(st.nextToken());
116  if (! st.hasMoreToken())
117  throw IOException("Sequence description should include a strand field.");
118  string tmp = st.nextToken();
119  if (tmp.size() != 1)
120  throw Exception("MafAlignmentParser::nextBlock. Strand specification is incorrect, should be only one character long, found " + TextTools::toString(tmp.size()) + ".");
121  char strand = tmp[0];
122 
123  if (! st.hasMoreToken())
124  throw IOException("Sequence description should include a source size field.");
125  unsigned int srcSize = TextTools::to<unsigned int>(st.nextToken());
126  if (currentSequence) {
127  //Add previous sequence:
128  block->addSequence(*currentSequence); //The sequence is copied in the container.
129  currentSequence.reset();
130  }
131  if (! st.hasMoreToken())
132  throw IOException("Sequence description without a sequence.");
133  string seq = st.nextToken();
134  if (dotOption_ == DOT_ASGAP) {
135  std::replace(seq.begin(), seq.end(), '.', '-');
136  }
137  if (dotOption_ == DOT_ASUNRES) {
138  std::replace(seq.begin(), seq.end(), '.', 'N');
139  }
140  currentSequence.reset(new MafSequence(src, seq, start, strand, srcSize));
141  if (currentSequence->getGenomicSize() != size) {
142  if (checkSequenceSize_)
143  throw Exception("MafAlignmentParser::nextBlock. Sequence found (" + src + ") does not match specified size: " + TextTools::toString(currentSequence->getGenomicSize()) + ", should be " + TextTools::toString(size) + ".");
144  else {
145  if (verbose_) {
146  ApplicationTools::displayWarning("MafAlignmentParser::nextBlock. Sequence found (" + src + ") does not match specified size: " + TextTools::toString(currentSequence->getGenomicSize()) + ", should be " + TextTools::toString(size) + ".");
147  }
148  }
149  }
150  //Add mask:
151  if (mask_) {
152  vector<bool> mask(currentSequence->size());
153  for (size_t i = 0; i < mask.size(); ++i) {
154  mask[i] = cmAlphabet_.isMasked(seq[i]);
155  }
156  currentSequence->addAnnotation(new SequenceMask(mask));
157  }
158  }
159  else if (line[0] == 'q')
160  {
161  if (!currentSequence)
162  throw Exception("MafAlignmentParser::nextBlock(). Quality scores found, but there is currently no sequence!");
163  StringTokenizer st(line);
164  st.nextToken(); //The 'q' tag
165  string name = st.nextToken();
166  if (name != currentSequence->getName())
167  throw Exception("MafAlignmentParser::nextBlock(). Quality scores found, but with a different name from the previous sequence: " + name + ", should be " + currentSequence->getName() + ".");
168  string qstr = st.nextToken();
169  //Now parse the score string:
170  SequenceQuality* seqQual = new SequenceQuality(qstr.size());
171  for (size_t i = 0; i < qstr.size(); ++i) {
172  char c = qstr[i];
173  if (c == '-') {
174  seqQual->setScore(i, -1);
175  } else if (c == '0' || c == '1' || c == '2' || c== '3' || c == '4' || c == '5' || c == '6' || c == '7' || c == '8' || c == '9') {
176  seqQual->setScore(i, c - '0');
177  } else if (c == 'F' || c == 'f') { //Finished
178  seqQual->setScore(i, 10);
179  } else if (c == '?' || c == '.') {
180  seqQual->setScore(i, -2);
181  } else {
182  throw Exception("MafAlignmentParser::nextBlock(). Unvalid quality score: " + TextTools::toString(c) + ". Should be 0-9, F or '-'.");
183  }
184  }
185  currentSequence->addAnnotation(seqQual);
186  }
187  }
189 
190  //In case last line in not empty:
191  if (currentSequence) {
192  //Add previous sequence:
193  block->addSequence(*currentSequence); //The sequence is copied in the container.
194  currentSequence.reset();
195  }
196 
197  //Returning block:
198  return block;
199 }
200 
A synteny block data structure, the basic unit of a MAF alignement file.
Definition: MafBlock.h:57
void setScore(double score)
Definition: MafBlock.h:102
void setPass(unsigned int pass)
Definition: MafBlock.h:103
void addSequence(const MafSequence &sequence)
Definition: MafBlock.h:115
A sequence class which is used to store data from MAF files.
Definition: MafSequence.h:64
void setScore(size_t pos, int score)
const std::string & nextToken()
bool hasMoreToken() const