forked from OpenMS/OpenMS
-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathTutorial_TOPP.cpp
More file actions
175 lines (132 loc) · 5.73 KB
/
Copy pathTutorial_TOPP.cpp
File metadata and controls
175 lines (132 loc) · 5.73 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
//! [doxygen_snippet_TOPPexample]
// Copyright (c) 2002-present, OpenMS Inc. -- EKU Tuebingen, ETH Zurich, and FU Berlin
// SPDX-License-Identifier: BSD-3-Clause
// --------------------------------------------------------------------------
// $Maintainer: Oliver Alka $
// $Authors: Oliver Alka $
// This file is ONLY used for code snippets in the developer tutorial
// --------------------------------------------------------------------------
//! [doxygen_snippet_Includes]
#include <OpenMS/APPLICATIONS/TOPPBase.h>
#include <OpenMS/FORMAT/FASTAFile.h>
#include <OpenMS/FORMAT/FileHandler.h>
#include <OpenMS/FORMAT/FileTypes.h>
#include <OpenMS/METADATA/PeptideIdentification.h>
//! [doxygen_snippet_Includes]
using namespace OpenMS;
using namespace std;
//-------------------------------------------------------------
// Doxygen docu
//-------------------------------------------------------------
/**
@page TOPP_DatabaseFilter DatabaseFilter
@brief The DatabaseFilter tool filters a protein database in fasta format according to one or multiple filtering criteria.
The resulting database is written as output. Depending on the reporting method (method="whitelist" or "blacklist") only entries are retained that
passed all filters ("whitelist) or failed at least one filter ("blacklist").
Implemented filter criteria:
ID: Filter database according to the set of proteinIDs contained in an identification file (idXML, mzIdentML)
<B>The command line parameters of this tool are:</B>
@verbinclude TOPP_DatabaseFilter.cli
<B>INI file documentation of this tool:</B>
@htmlinclude TOPP_DatabaseFilter.html
*/
// We do not want this class to show up in the docu:
/// @cond TOPPCLASSES
class TOPPDatabaseFilter : public TOPPBase
{
public:
TOPPDatabaseFilter():
TOPPBase("DatabaseFilter", "Filters a protein database (FASTA format) based on identified proteins", false) // false: mark as unofficial tool
{
}
protected:
//! [doxygen_snippet_Register]
void registerOptionsAndFlags_() override
{
registerInputFile_("in", "<file>", "", "Input FASTA file, containing a protein database.");
setValidFormats_("in", {"fasta"});
registerInputFile_("id", "<file>", "", "Input file containing identified peptides and proteins.");
setValidFormats_("id", {"idXML", "mzid"});
registerStringOption_("method", "<choice>", "whitelist", "Switch between white-/blacklisting of protein IDs", false);
setValidStrings_("method", {"whitelist", "blacklist"});
registerOutputFile_("out", "<file>", "", "Output FASTA file where the reduced database will be written to.");
setValidFormats_("out", {"fasta"});
}
//! [doxygen_snippet_Register]
//! [doxygen_snippet_Functionality_1]
void filterByProteinAccessions_(const vector<FASTAFile::FASTAEntry>& db,
const vector<PeptideIdentification>& peptide_identifications,
bool whitelist,
vector<FASTAFile::FASTAEntry>& db_new)
{
set<String> id_accessions;
for (const auto& pep_id : peptide_identifications)
{
for (const auto& hit : pep_id.getHits())
{
for (const auto& ev : hit.getPeptideEvidences())
{
const String& id_accession = ev.getProteinAccession();
id_accessions.insert(id_accession);
}
}
}
//! [doxygen_snippet_Functionality_1]
OPENMS_LOG_INFO << "Number of Protein IDs: " << id_accessions.size() << endl;
//! [doxygen_snippet_Functionality_2]
for (const auto entry : db)
{
const String& fasta_accession = entry.identifier;
const bool found = id_accessions.find(fasta_accession) != id_accessions.end();
if ((found && whitelist) || (! found && ! whitelist)) // either found in the whitelist or not found in the blacklist
{
db_new.push_back(entry);
}
}
//! [doxygen_snippet_Functionality_2]
}
ExitCodes main_(int, const char**) override
{
//! [doxygen_snippet_InputParam]
//-------------------------------------------------------------
// parsing parameters
//-------------------------------------------------------------
String in(getStringOption_("in"));
String ids(getStringOption_("id"));
String method(getStringOption_("method"));
bool whitelist = (method == "whitelist");
String out(getStringOption_("out"));
//! [doxygen_snippet_InputParam]
//-------------------------------------------------------------
// reading input
//-------------------------------------------------------------
//! [doxygen_snippet_InputRead]
vector<FASTAFile::FASTAEntry> db;
FASTAFile().load(in, db);
//! [doxygen_snippet_InputRead]
vector<ProteinIdentification> protein_identifications;
vector<PeptideIdentification> peptide_identifications;
FileHandler().loadIdentifications(ids, protein_identifications, peptide_identifications);
OPENMS_LOG_INFO << "Identifications: " << ids.size() << endl;
// run filter
vector<FASTAFile::FASTAEntry> db_new;
filterByProteinAccessions_(db, peptide_identifications, whitelist, db_new);
//-------------------------------------------------------------
// writing output
//-------------------------------------------------------------
OPENMS_LOG_INFO << "Database entries (before / after): " << db.size() << " / " << db_new.size() << endl;
//! [doxygen_snippet_output]
FASTAFile().store(out, db_new);
//! [doxygen_snippet_output]
return EXECUTION_OK;
}
};
int main(int argc, const char** argv)
{
TOPPDatabaseFilter tool;
OPENMS_LOG_FATAL_ERROR << "THIS IS TEST CODE AND SHOULD NEVER BE RUN OUTSIDE OF TESTING" << endl;
tool.main(argc, argv);
return 0;
}
/// @endcond
//! [doxygen_snippet_TOPPexample]