revp

revp.cpp

#include <string>
#include <iostream>
#include <fstream>
#include <vector>
#include <algorithm>
#include <cassert>
using namespace std;

//////////////////////////////////////////////////////
// Finding restriction sites (reverse palindromes)  //
//////////////////////////////////////////////////////

vector<string> parse_fasta(string file_path) 
{
    vector<string> fasta_out; 
    string line, first_char; 
    ifstream infile(file_path);  
    if (infile.is_open()) { 
        while (getline(infile, line)) {
            first_char = line[0]; 
            if (first_char != ">") { 
                fasta_out.push_back(line); 
            } 
        } 
    } else cerr << "OH NO! can't open file"; 
    infile.close(); 
    return fasta_out;
}


string complement(const string& s)
{
    string n, out_string;
    for (int i=s.length()-1; i>=0; i--) {
        n = toupper(s[i]);
        if (n == "A") {
            n = "T";
        } else if (n == "T") {
            n = "A";
        } else if (n == "C") {
            n = "G";
        } else if (n == "G") {
            n = "C";
        } else {
            cerr << "Unknown char: " << static_cast<int>(s[i]) << endl;
        }
        out_string += n;
    }
    assert(out_string.length() == s.length());
    return out_string;
}


void find_restrction_site(const vector<string> v, const int floor=4, const int ceiling=12)
{   
    // find reverse palindromes between length floor and ceiling
    const string s = string(v[0]);
    string win, first_half, second_half, c_second_half;

    // sliding window? ... sliding window!!
    for (int win_size = floor; win_size <= ceiling; win_size+=2) {
        for (int i=0; i<=(s.length()-win_size); ++i) {
            win = s.substr(i, win_size); // window sub-sequence

            // split sub-sequence into two halves
            first_half = win.substr(0, win.length()/2);
            second_half = win.substr(win.length()/2, win.length());
            // get reverse complement of second half
            c_second_half = complement(second_half);

            if (first_half == c_second_half) {
                // return index and window size
                cout << i+1 << "\t" << win_size << endl;
            }
        }
    }
    cout << endl;
}


int main()
{
    vector<string> fasta = parse_fasta("/home/scott/Downloads/rosalind_revp.txt");
    find_restrction_site(fasta);
    return 0;
}

compile with g++ -std=c++11 revp.cpp

output:

1   4
9   4
21  4
64  4
79  4
100 4
124 4
135 4
152 4
164 4
166 4
173 4
181 4
187 4
203 4
213 4
233 4
245 4
283 4
301 4
307 4
311 4
355 4
356 4
365 4
377 4
400 4
403 4
406 4
417 4
438 4
443 4
458 4
476 4
486 4
498 4
506 4
580 4
594 4
597 4
605 4
618 4
629 4
633 4
670 4
678 4
683 4
688 4
741 4
750 4
757 4
774 4
805 4
244 6
310 6
376 6
399 6
402 6
405 6
416 6
475 6
497 6
593 6
596 6
604 6
617 6
773 6
401 8
616 8
400 10
399 12