add: sortdiamond and splitfasta
This commit is contained in:
parent
fdc588efc7
commit
e090370586
3 changed files with 236 additions and 0 deletions
12
CMakeLists.txt
Normal file
12
CMakeLists.txt
Normal file
|
@ -0,0 +1,12 @@
|
|||
cmake_minimum_required(VERSION 3.0)
|
||||
|
||||
project(RGB EPP)
|
||||
SET( CMAKE_EXPORT_COMPILE_COMMANDS ON )
|
||||
|
||||
set(CMAKE_CXX_STANDARD 17)
|
||||
|
||||
add_executable(splitfasta splitfasta.cpp)
|
||||
add_executable(sortdiamond sortdiamond.cpp)
|
||||
|
||||
target_compile_options(splitfasta PRIVATE -Wall -Wextra -pedantic -O3)
|
||||
target_compile_options(sortdiamond PRIVATE -Wall -Wextra -pedantic -O3)
|
143
sortdiamond.cpp
Normal file
143
sortdiamond.cpp
Normal file
|
@ -0,0 +1,143 @@
|
|||
#include <algorithm>
|
||||
#include <fstream>
|
||||
#include <iostream>
|
||||
#include <map>
|
||||
#include <sstream>
|
||||
#include <vector>
|
||||
|
||||
using namespace std;
|
||||
|
||||
// Function to generate reverse complement of a DNA sequence
|
||||
string revcomp(const string &seq) {
|
||||
string revseq;
|
||||
for (auto it = seq.rbegin(); it != seq.rend(); ++it) {
|
||||
switch (*it) {
|
||||
case 'A':
|
||||
revseq += 'T';
|
||||
break;
|
||||
case 'T':
|
||||
revseq += 'A';
|
||||
break;
|
||||
case 'C':
|
||||
revseq += 'G';
|
||||
break;
|
||||
case 'G':
|
||||
revseq += 'C';
|
||||
break;
|
||||
case 'R':
|
||||
revseq += 'Y';
|
||||
break;
|
||||
case 'Y':
|
||||
revseq += 'R';
|
||||
break;
|
||||
case 'S':
|
||||
revseq += 'S';
|
||||
break;
|
||||
case 'W':
|
||||
revseq += 'W';
|
||||
break;
|
||||
case 'K':
|
||||
revseq += 'M';
|
||||
break;
|
||||
case 'M':
|
||||
revseq += 'K';
|
||||
break;
|
||||
case 'B':
|
||||
revseq += 'V';
|
||||
break;
|
||||
case 'D':
|
||||
revseq += 'H';
|
||||
break;
|
||||
case 'H':
|
||||
revseq += 'D';
|
||||
break;
|
||||
case 'V':
|
||||
revseq += 'B';
|
||||
break;
|
||||
case 'N':
|
||||
revseq += 'N';
|
||||
break;
|
||||
default:
|
||||
break;
|
||||
}
|
||||
}
|
||||
return revseq;
|
||||
}
|
||||
|
||||
int main(int argc, char *argv[]) {
|
||||
if (argc != 3) {
|
||||
cerr << "Usage: " << argv[0] << " <input_file> <output_file>"
|
||||
<< endl;
|
||||
return 1;
|
||||
}
|
||||
|
||||
string in_name = argv[1];
|
||||
string ot_name = argv[2];
|
||||
|
||||
ifstream infile(in_name);
|
||||
if (!infile) {
|
||||
cerr << "Error opening input file: " << in_name << endl;
|
||||
return 1;
|
||||
}
|
||||
|
||||
ofstream outfile(ot_name);
|
||||
if (!outfile) {
|
||||
cerr << "Error opening output file: " << ot_name << endl;
|
||||
return 1;
|
||||
}
|
||||
|
||||
map<string, pair<double, string>>
|
||||
max_map; // Key: sseqid, Value: pair<score, line>
|
||||
string line;
|
||||
|
||||
while (getline(infile, line)) {
|
||||
istringstream iss(line);
|
||||
string fields[20]; // Adjust size if needed
|
||||
int i = 0;
|
||||
|
||||
while (iss >> fields[i] && i < 20) {
|
||||
i++;
|
||||
}
|
||||
|
||||
string key = fields[1];
|
||||
double value = stod(fields[11]);
|
||||
|
||||
if (max_map.find(key) != max_map.end()) {
|
||||
if (value > max_map[key].first) {
|
||||
max_map[key] = make_pair(value, line);
|
||||
}
|
||||
} else {
|
||||
max_map[key] = make_pair(value, line);
|
||||
}
|
||||
}
|
||||
|
||||
vector<pair<string, string>> result;
|
||||
|
||||
for (const auto &entry : max_map) {
|
||||
istringstream iss(entry.second.second);
|
||||
string fields[20]; // Adjust size if needed
|
||||
int i = 0;
|
||||
|
||||
while (iss >> fields[i] && i < 20) {
|
||||
i++;
|
||||
}
|
||||
|
||||
if (stoi(fields[6]) > stoi(fields[7])) {
|
||||
fields[17] = revcomp(fields[17]);
|
||||
}
|
||||
|
||||
result.push_back(make_pair(">" + fields[1], fields[17]));
|
||||
}
|
||||
|
||||
sort(result.begin(), result.end());
|
||||
|
||||
for (const auto &entry : result) {
|
||||
outfile << entry.first << "\n" << entry.second << endl;
|
||||
}
|
||||
|
||||
infile.close();
|
||||
outfile.close();
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
81
splitfasta.cpp
Normal file
81
splitfasta.cpp
Normal file
|
@ -0,0 +1,81 @@
|
|||
// #include <sys/stat.h>
|
||||
|
||||
#include <filesystem>
|
||||
#include <fstream>
|
||||
#include <iostream>
|
||||
// #include <sstream>
|
||||
#include <string>
|
||||
// #include <vector>
|
||||
|
||||
namespace fs = std::filesystem;
|
||||
|
||||
std::string removeExtension(const std::string& filename) {
|
||||
size_t lastdot = filename.find_last_of(".");
|
||||
if (lastdot != std::string::npos) {
|
||||
return filename.substr(0, lastdot);
|
||||
}
|
||||
return filename;
|
||||
}
|
||||
|
||||
void splitFasta(const std::string& input_fasta) {
|
||||
std::ifstream infile(input_fasta);
|
||||
if (!infile) {
|
||||
std::cerr << "Error: Unable to open input file " << input_fasta
|
||||
<< std::endl;
|
||||
return;
|
||||
}
|
||||
|
||||
std::string line;
|
||||
std::string dir_name;
|
||||
std::ofstream outfile;
|
||||
bool in_sequence = false;
|
||||
|
||||
while (std::getline(infile, line)) {
|
||||
if (line.empty()) continue;
|
||||
|
||||
if (line[0] == '>') {
|
||||
// New sequence header
|
||||
if (in_sequence) {
|
||||
outfile.close();
|
||||
}
|
||||
dir_name = line.substr(1); // Remove '>'
|
||||
// mkdir(dir_name.c_str(), 0777); // Create
|
||||
// directory
|
||||
fs::create_directories(dir_name);
|
||||
fs::path output_file = fs::path(dir_name) / input_fasta;
|
||||
outfile.open(output_file);
|
||||
outfile << ">" << removeExtension(input_fasta)
|
||||
<< std::endl;
|
||||
in_sequence = true;
|
||||
} else if (in_sequence) {
|
||||
// Inside sequence content
|
||||
outfile << line << std::endl;
|
||||
}
|
||||
}
|
||||
|
||||
if (in_sequence) {
|
||||
outfile.close();
|
||||
}
|
||||
|
||||
if (infile.eof()) {
|
||||
std::cout << "Sequences have been split into individual files."
|
||||
<< std::endl;
|
||||
} else {
|
||||
std::cerr << "Error occurred while reading file." << std::endl;
|
||||
}
|
||||
|
||||
infile.close();
|
||||
}
|
||||
|
||||
int main(int argc, char* argv[]) {
|
||||
if (argc != 2) {
|
||||
std::cerr << "Usage: " << argv[0] << " <input_fasta>"
|
||||
<< std::endl;
|
||||
return 1;
|
||||
}
|
||||
|
||||
splitFasta(argv[1]);
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
Loading…
Reference in a new issue