binaries/src/ViennaRNA/lib/RNAstruct.c

   1 /*
   2                 parse and convert secondary structures
   3            Walter Fontana, Ivo L Hofacker, Peter F Stadler
   4                         Vienna RNA Package
   5 */
   6 /* Last changed Time-stamp: <2005-07-23 10:12:19 ivo> */
   7
   8
   9 #include <stdio.h>
  10 #include <stdlib.h>
  11 #include <string.h>
  12 #include <ctype.h>
  13 #include "utils.h"
  14 #include "RNAstruct.h"
  15
  16 #define PRIVATE  static
  17 #define PUBLIC
  18
  19 #define MAXLEN    10000
  20
  21 static char rcsid[] = "$Id: RNAstruct.c,v 1.8 2006/05/09 20:44:33 ivo Exp $";
  22
  23
  24 PRIVATE char *aux_struct(const char *structure);
  25
  26 /* on return from parse_structure(), b2C() or b2Shapiro() ... */
  27 PUBLIC int    loop_size[STRUC];       /* contains loop sizes of a structure */
  28 PUBLIC int    helix_size[STRUC];      /* contains helix sizes of a structure */
  29 PUBLIC int    loop_degree[STRUC];     /* contains loop degrees of a structure */
  30 PUBLIC int    loops;                  /* n of loops and stacks in a structure */
  31 PUBLIC int    unpaired, pairs;        /* n of unpaired digits and pairs */
  32
  33 /*---------------------------------------------------------------------------*/
  34
  35 PRIVATE char *aux_struct(const char* structure )
  36 {
  37    short        *match_paren;
  38    int          i, o, p;
  39    char        *string;
  40
  41    string = (char *) space(sizeof(char)*(strlen(structure)+1));
  42    match_paren = (short *) space(sizeof(short)*(strlen(structure)/2+1));
  43    strcpy(string, structure);
  44
  45    i = o = 0;
  46    while (string[i]) {
  47       switch (string[i]) {
  48        case '.': break;
  49        case '(':
  50          match_paren[++o]=i;
  51          break;
  52        case ')':
  53          p=i;
  54          while ((string[p+1]==')')&&(match_paren[o-1]==match_paren[o]-1)) {
  55             p++; o--;
  56          }
  57          string[p]=']';
  58          i=p;
  59          string[match_paren[o]]='[';
  60          o--;
  61          break;
  62        default:
  63          nrerror("Junk in structure at aux_structure\n");
  64       }
  65       i++;
  66    }
  67    free(match_paren);
  68    return(string);
  69 }
  70
  71 /*---------------------------------------------------------------------------*/
  72
  73 PUBLIC char *b2HIT(const char *structure)
  74 {
  75
  76    int            i, u, p, l;
  77    char          *string, *temp, *HIT, tt[10];
  78
  79    temp = (char *) space(strlen(structure)*4+4);
  80    string = aux_struct( structure );
  81
  82    strcpy(temp,"(");
  83    i=p=u=0; l=1;
  84    while (string[i]) {
  85       switch(string[i]) {
  86        case '.':
  87          u++; break;
  88        case '[':
  89          if (u>0) {
  90             sprintf(tt, "(U%d)" , u);
  91             strcat(temp+l, tt);
  92             l+=strlen(tt);
  93             u=0;
  94          }
  95          strcat(temp+l, "("); l++;
  96          break;
  97        case ')':
  98          if (u>0) {
  99             sprintf(tt, "(U%d)" , u);
 100             strcat(temp+l, tt);
 101             l+=strlen(tt);
 102             u=0;
 103          }
 104          p++;
 105          break;
 106        case ']':
 107          if (u>0) {
 108             sprintf(tt, "(U%d)" , u);
 109             strcat(temp+l, tt);
 110             l+=strlen(tt);
 111             u=0;
 112          }
 113          sprintf(tt,"P%d)", p+1);
 114          strcat(temp+l, tt);
 115          l+=strlen(tt);
 116          p=0;
 117          break;
 118       }
 119       i++;
 120    }
 121    if (u>0) {
 122       sprintf(tt, "(U%d)" , u);
 123       strcat(temp+l, tt);
 124       l+=strlen(tt);
 125    }
 126    strcat(temp+l, "R)");
 127
 128    free( string );
 129
 130    HIT = (char *) space(sizeof(char)*(strlen(temp)+2));
 131    strcpy(HIT, temp);
 132    free(temp);
 133    return(HIT);
 134 }
 135
 136 /*---------------------------------------------------------------------------*/
 137
 138 PUBLIC char *b2C(const char *structure )
 139 {
 140    short *bulge, *loop;
 141
 142    int    i, lp, p, l;
 143    char  *string, *Coarse, *temp;
 144
 145    bulge = (short *) space(sizeof(short)*(strlen(structure)/3+1));
 146    loop = (short *) space(sizeof(short)*(strlen(structure)/3+1));
 147    temp = (char *) space(4*strlen(structure)+2);
 148
 149    for (i = 0; i < STRUC; i++) {
 150       loop_size[i] = helix_size[i] = 0;
 151    }
 152    loop_degree[0]=0;         /* open structure has degree 0 */
 153    pairs = unpaired = loops = lp = 0;
 154    loop[0]=0;
 155
 156    string = aux_struct( structure );
 157
 158    i=p=l=0;
 159    temp[l++] = '(';
 160    while (string[i]) {
 161       switch(string[i]) {
 162        case '.':
 163          loop_size[loop[lp]]++;
 164          break;
 165        case '[':
 166          temp[l++]='(';
 167          if ((i>0)&&(string[i-1]=='(')) bulge[lp]=1;
 168          lp++;
 169          loop_degree[++loops]=1;
 170          loop[lp]=loops;
 171          bulge[lp]=0;
 172          break;
 173        case ')':
 174          if (string[i-1]==']') bulge[lp]=1;
 175          p++;
 176          break;
 177        case ']':
 178          if (string[i-1]==']') bulge[lp]=1;
 179          switch (loop_degree[loop[lp]]) {
 180           case 1:  temp[l++]='H'; break;           /* hairpin */
 181           case 2:
 182             if (bulge[lp]==1)
 183                temp[l++] = 'B';                    /* bulge */
 184             else
 185                temp[l++] = 'I';                    /* internal loop */
 186             break;
 187           default: temp[l++] = 'M';                /* multiloop */
 188          }
 189          temp[l++] = ')';
 190          pairs+=p+1;
 191          p=0;
 192          loop_degree[loop[--lp]]++;
 193          break;
 194       }
 195       i++;
 196    }
 197    temp[l++] = 'R';
 198    temp[l++] = ')';
 199    temp[l]='\0';
 200    free(string);
 201    Coarse = (char *) space(sizeof(char)*(strlen(temp)+2));
 202    strcpy(Coarse, temp);
 203    free(temp);
 204    free(bulge); free(loop);
 205    return(Coarse);
 206 }
 207
 208 /*---------------------------------------------------------------------------*/
 209
 210 PUBLIC char *b2Shapiro(const char *structure )
 211 {
 212
 213    short *bulge, *loop;
 214
 215    int            i, lp, p, l, k;
 216    char          *string, *Shapiro, *temp, tt[10];
 217
 218    bulge = (short *) space(sizeof(short)*(strlen(structure)/3+1));
 219    loop = (short *) space(sizeof(short)*(strlen(structure)/3+1));
 220    temp = (char *) space(4*strlen(structure)+3);
 221
 222    for (i = 0; i < STRUC; i++) {
 223       loop_size[i] = helix_size[i] = 0;
 224    }
 225    loop_degree[0]=0;         /* open structure has degree 0 */
 226    pairs = unpaired = loops = lp = 0;
 227    loop[0]=0;
 228
 229    string = aux_struct( structure );
 230
 231    i=p=l=0;
 232    temp[l++] = '(';    /* root */
 233    while (string[i]) {
 234       switch(string[i]) {
 235        case '.':
 236          unpaired++;
 237          loop_size[loop[lp]]++;
 238          break;
 239        case '[':
 240          temp[l++]='(';
 241          temp[l++]='(';
 242          if ((i>0)&&(string[i-1]=='(' || string[i-1]=='['))
 243            bulge[lp]=1;
 244          lp++;
 245          loop_degree[++loops]=1;
 246          loop[lp]=loops;
 247          bulge[lp]=0;
 248          break;
 249        case ')':
 250          if (string[i-1]==']') bulge[lp]=1;
 251          p++;
 252          break;
 253        case ']':
 254          if (string[i-1]==']') bulge[lp]=1;
 255          switch (loop_degree[loop[lp]]) {
 256           case 1:  temp[l++]='H'; break;           /* hairpin */
 257           case 2:
 258             if (bulge[lp]==1)
 259                temp[l++] = 'B';                    /* bulge */
 260             else
 261                temp[l++] = 'I';                    /* internal loop */
 262             break;
 263           default: temp[l++] = 'M';                /* multiloop */
 264          }
 265          helix_size[loop[lp]]=p+1;
 266
 267          sprintf(tt, "%d)" , loop_size[loop[lp]]);
 268          for(k=0; k<strlen(tt); k++) temp[l++] = tt[k];
 269          sprintf(tt, "S%d)" , helix_size[loop[lp]]);
 270          for(k=0; k<strlen(tt); k++) temp[l++] = tt[k];
 271
 272          pairs+=p+1;
 273          p=0;
 274          loop_degree[loop[--lp]]++;
 275          break;
 276       }
 277       i++;
 278    }
 279
 280    *tt = '\0';
 281    if (loop_size[0]) sprintf(tt, "E%d)" , loop_size[0]);
 282    strcat(tt,"R)");
 283    temp[l]='\0';
 284    strcat(temp, tt);
 285    Shapiro = (char *) space(sizeof(char)*(strlen(temp)+2));
 286    if (loop_size[0]) {
 287       Shapiro[0]='(';
 288       strcpy(Shapiro+1, temp);
 289    } else strcpy(Shapiro, temp);
 290    free(string);
 291    free(temp);
 292    free(loop); free(bulge);
 293    return Shapiro;
 294 }
 295
 296
 297
 298
 299 /*---------------------------------------------------------------------------*/
 300
 301 PUBLIC void parse_structure(const char *structure)
 302
 303 /*-----------------------------------------------------------------------------
 304
 305     upon return from parse_structure():
 306
 307     loops    ....................... number of loops or stacks in structure.
 308     loop_size[1 <= i <= loops] ..... size of i-th loop.
 309     loop_size[0] ................... number of external digits.
 310     loop_degree[1 <= i <= loops] ... degree (branches) of i-th loop.
 311     loop_degree[0] ................. number of components.
 312     helix_size[1 <= i <= loops] .... size of i-th stack.
 313     unpaired ....................... n of unpaired digits.
 314     pairs .......................... n of base pairs.
 315
 316 -----------------------------------------------------------------------------*/
 317
 318 {
 319    short  *bulge, *loop;
 320
 321    int            i, lp, p;
 322    char          *string, *temp;
 323
 324    temp = (char *)  space(strlen(structure)*4+2);
 325    bulge = (short *) space(sizeof(short)*(strlen(structure)/3+1));
 326    loop = (short *) space(sizeof(short)*(strlen(structure)/3+1));
 327
 328    for (i = 0; i < STRUC; i++) {
 329       loop_size[i] = helix_size[i] = 0;
 330    }
 331    loop[0] = loop_degree[0]=0;         /* open structure has degree 0 */
 332    pairs = unpaired = loops = lp = 0;
 333    *temp='\0';
 334
 335    string = aux_struct(structure);
 336
 337    i=p=0;
 338    while (string[i]) {
 339       switch(string[i]) {
 340        case '.':
 341          unpaired++;
 342          loop_size[loop[lp]]++;
 343          break;
 344        case '[':
 345          if ((i>0)&&(string[i-1]=='(')) bulge[lp]=1;
 346          lp++;
 347          loop_degree[++loops]=1;
 348          loop[lp]=loops;
 349          bulge[lp]=0;
 350          break;
 351        case ')':
 352          if (string[i-1]==']') bulge[lp]=1;
 353          p++;
 354          break;
 355        case ']':
 356          if (string[i-1]==']') bulge[lp]=1;
 357          helix_size[loop[lp]]=p+1;
 358          pairs+=p+1;
 359          p=0;
 360          loop_degree[loop[--lp]]++;
 361          break;
 362       }
 363       i++;
 364    }
 365    free(string);
 366    free(bulge); free(loop);
 367    free(temp);
 368 }
 369
 370 /*---------------------------------------------------------------------------*/
 371
 372 PUBLIC char *add_root(const char *structure)
 373 {
 374     char *xS;
 375     xS = (char *) space(sizeof(char)*(strlen(structure)+4));
 376     xS[0] = '(';
 377     strcat(xS,structure);
 378     strcat(xS,"R)");
 379     return xS;
 380 }
 381
 382
 383 /*---------------------------------------------------------------------------*/
 384
 385 PUBLIC char *expand_Shapiro(const char *structure)
 386 {
 387    char  *xS, *temp;
 388    int  i, l;
 389
 390    temp = (char *) space(4*strlen(structure)+2);
 391
 392    i = 1;
 393    l = 1;
 394    temp[0] = '(';
 395    while (i<strlen(structure)-1) {
 396       temp[l++] = structure[i];
 397       if      (structure[i] == '(') temp[l++] = '(';
 398       else if (structure[i] == ')') {
 399          temp[l++] = 'S';
 400          temp[l++] = ')';
 401       }
 402       i++;
 403    }
 404    temp[l++] = ')';
 405    temp[l] = '\0';
 406
 407    xS = (char *) space(sizeof(char)*(strlen(temp)+1));
 408    strcpy(xS, temp);
 409    free(temp);
 410    return (xS);
 411 }
 412
 413 /*---------------------------------------------------------------------------*/
 414
 415 PUBLIC char *expand_Full(const char *structure)
 416 {
 417     char *xF, *temp;
 418     int  i, l;
 419
 420     temp = (char *) space(4*strlen(structure)+2);
 421
 422     i = 0;
 423     l = 0;
 424     while (structure[i]) {
 425         if      (structure[i] == '(') temp[l++] = '(';
 426         else if (structure[i] == ')') {
 427             temp[l++] = 'P';
 428             temp[l++] = ')';
 429         }
 430         else {
 431             temp[l++] = '(';
 432             temp[l++] = 'U';
 433             temp[l++] = ')';
 434         }
 435         i++;
 436      }
 437      temp[l] = '\0';
 438
 439      xF = (char *) space(sizeof(char)*(l+5));
 440      strcpy(xF, "(");
 441      strcat(xF, temp);
 442      strcat(xF, "R)");
 443      free(temp);
 444      return (xF);
 445 }
 446
 447 /*---------------------------------------------------------------------------*/
 448
 449 PUBLIC char *unexpand_Full(const char *structure)
 450 {
 451    short        *match_paren;
 452    char id[10], *full, *temp;
 453    int    i, j, k, l, o, w;
 454
 455    temp = (char *) space(4*strlen(structure)+2);
 456    match_paren = (short *) space(sizeof(short)*(strlen(structure)/2+1));
 457
 458    i = strlen(structure)-1;
 459    l = o = 0; k=9;
 460    id[9]='\0';
 461    while (i>=0) {
 462      switch (structure[i]) {
 463      case '(':
 464        for (j=0; j<match_paren[o]; j++) temp[l++]='(';
 465        match_paren[o--] = 0;
 466        break;
 467      case 'U':
 468        w=1;
 469        sscanf(id+k, "%d", &w);
 470        for (j=0; j<w; j++) temp[l++]='.';
 471        k=9;
 472        break;
 473      case 'P':
 474            w=1;
 475        sscanf(id+k, "%d", &w);
 476        for (j=0; j<w; j++) temp[l++]=')';
 477        match_paren[o]=w;
 478        k=9;
 479        break;
 480      case 'R':
 481        break;
 482      case ')':
 483        o++;
 484        break;
 485      default:
 486        id[--k]=structure[i];
 487      }
 488      i--;
 489    }
 490
 491    temp[l] = '\0';
 492    full = (char *) space(sizeof(char)*(l+1));
 493    for (i=0; i<l; i++) full[i]=temp[l-i-1];
 494    full[l]='\0';
 495    free(temp);
 496    free(match_paren);
 497    return full;
 498 }
 499
 500
 501 /*---------------------------------------------------------------------------*/
 502
 503 PUBLIC char *unweight(const char *structure)
 504 {
 505    int i,l;
 506    char *full, *temp;
 507
 508    temp = (char *) space(4*strlen(structure)+1);
 509
 510    i=l=0;
 511    while (structure[i]) {
 512       if (!isdigit((int)structure[i])) temp[l++]=structure[i];
 513       i++;
 514    }
 515    temp[l]='\0';
 516    full = (char *) space(sizeof(char)*(l+1));
 517    strcpy(full, temp);
 518    free(temp);
 519    return full;
 520 }
 521
 522 /*---------------------------------------------------------------------------*/
 523
 524 PUBLIC void unexpand_aligned_F(char *align[2])
 525 {
 526    char *t0, *t1;
 527    int i,l;
 528
 529    t0 = (char *) space(strlen(align[0])+1);
 530    t1 = (char *) space(strlen(align[0])+1);
 531
 532    for (i=0, l=0; i<strlen(align[0]); i++) {
 533       switch (align[0][i]) {
 534        case '(':
 535        case ')':
 536          t0[l] = align[0][i];
 537          t1[l++]=align[1][i];
 538          break;
 539        case 'U':
 540          switch (align[1][i]) {
 541           case 'U':
 542             t0[l-1]=t1[l-1]='.';
 543             break;
 544           case '_':
 545             t0[l-1]='.';
 546             t1[l-1]='_';
 547             break;
 548           case 'P':
 549             t0[l-1]='_'; t0[l]='.';
 550             t1[l-1]='('; t1[l]=')'; l++;
 551          }
 552          while (align[0][i]!=')') i++;
 553          break;
 554        case '_':
 555          switch (align[1][i]) {
 556           case '(':
 557           case ')':
 558             t0[l] = align[0][i];
 559             t1[l++]=align[1][i];
 560             break;
 561           case 'U':
 562             while (align[1][i]!=')') i++;
 563             t1[l-1]='.';
 564             t0[l-1]='_';
 565             break;
 566          }
 567        case 'P':
 568          if (align[1][i]=='U') {
 569             t1[l-1]='_'; t1[l]='.'; t0[l++]=')';
 570             while (align[0][i]!=')') i++;
 571          }
 572          break;
 573       }
 574    }
 575    t0[l-1]=t1[l-1]='\0';
 576    strcpy(align[0], t0+1);
 577    strcpy(align[1], t1+1);
 578    free(t0); free(t1);
 579 }