/*
Implementation of Nussinov Maximum Matching
Ronny Lorenz
*/
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <ctype.h>
#include <string.h>
#include "ViennaRNA/utils.h"
#include "ViennaRNA/energy_par.h"
#include "ViennaRNA/fold_vars.h"
#include "ViennaRNA/pair_mat.h"
#include "ViennaRNA/params.h"
/* the encoded string MUST have the length of the sequence at position 0!!! */
PUBLIC unsigned int maximumMatching(const char *string){
unsigned int i, j, l, length, max = 0;
unsigned int *mm; /* holds maximum matching on subsequence [i,j] */
short *encodedString = encode_sequence(string, 0);
int *iindx = vrna_idx_row_wise((unsigned) encodedString[0]);
make_pair_matrix();
length = (unsigned int)encodedString[0];
mm = (unsigned int *) vrna_alloc(sizeof(unsigned int)*((length*(length+1))/2+2));
for(j = 1; j<=length; j++)
for(i=(j>TURN?(j-TURN):1); i<j; i++)
mm[iindx[i]-j] = 0;
for(i=length-TURN-1;i>0; i--)
for(j=i+TURN+1; j<= length; j++){
max = mm[iindx[i]-j+1];
for(l=j-TURN-1; l>=i; l--)
if(pair[encodedString[l]][encodedString[j]])
max = MAX2(max, ((l>i) ? mm[iindx[i]-l+1] : 0) + 1 + mm[iindx[l+1]-j+1]);
mm[iindx[i]-j] = max;
}
max = mm[iindx[1]-length];
free(mm);
free(iindx);
free(encodedString);
return max;
}
/* the encoded string MUST have the length of the sequence at position 0!!! */
PUBLIC unsigned int *maximumMatchingConstraint(const char *string, short *ptable){
unsigned int i, j, l, length, max = 0;
unsigned int *mm; /* holds maximum matching on subsequence [i,j] */
short *encodedString = encode_sequence(string, 0);
int *iindx = vrna_idx_row_wise((unsigned) encodedString[0]);
make_pair_matrix();
length = (unsigned int)encodedString[0];
mm = (unsigned int *) vrna_alloc(sizeof(unsigned int)*((length*(length+1))/2+2));
for(j = 1; j<=length; j++)
for(i=(j>TURN?(j-TURN):1); i<j; i++)
mm[iindx[i]-j] = 0;
for(i=length-TURN-1;i>0; i--)
for(j=i+TURN+1; j<= length; j++){
max = mm[iindx[i]-j+1];
for(l=j-TURN-1; l>=i; l--){
if(pair[encodedString[l]][encodedString[j]]){
if(ptable[l] != j)
max = MAX2(max, ((l>i) ? mm[iindx[i]-l+1] : 0) + 1 + mm[iindx[l+1]-j+1]);
}
}
mm[iindx[i]-j] = max;
}
free(iindx);
free(encodedString);
return mm;
}
/* the encoded string MUST have the length of the sequence at position 0!!! */
PUBLIC unsigned int *maximumMatching2Constraint(const char *string, short *ptable, short *ptable2){
unsigned int i, j, l, length, max = 0;
unsigned int *mm; /* holds maximum matching on subsequence [i,j] */
short *encodedString = encode_sequence(string, 0);
int *iindx = vrna_idx_row_wise((unsigned) encodedString[0]);
make_pair_matrix();
length = (unsigned int)encodedString[0];
mm = (unsigned int *) vrna_alloc(sizeof(unsigned int)*((length*(length+1))/2+2));
for(j = 1; j<=length; j++)
for(i=(j>TURN?(j-TURN):1); i<j; i++)
mm[iindx[i]-j] = 0;
for(i=length-TURN-1;i>0; i--)
for(j=i+TURN+1; j<= length; j++){
max = mm[iindx[i]-j+1];
for(l=j-TURN-1; l>=i; l--){
if(pair[encodedString[l]][encodedString[j]]){
if(ptable[l] != j && ptable2[l] != j)
max = MAX2(max, ((l>i) ? mm[iindx[i]-l+1] : 0) + 1 + mm[iindx[l+1]-j+1]);
}
}
mm[iindx[i]-j] = max;
}
free(iindx);
free(encodedString);
return mm;
}