BGI/SOAPdenovo2/standardPregraph/loadPath.c

   1 /*
   2  * loadPath.c
   3  *
   4  * Copyright (c) 2008-2012 BGI-Shenzhen <soap at genomics dot org dot cn>.
   5  *
   6  * This file is part of SOAPdenovo.
   7  *
   8  * SOAPdenovo is free software: you can redistribute it and/or modify
   9  * it under the terms of the GNU General Public License as published by
  10  * the Free Software Foundation, either version 3 of the License, or
  11  * (at your option) any later version.
  12  *
  13  * SOAPdenovo is distributed in the hope that it will be useful,
  14  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  16  * GNU General Public License for more details.
  17  *
  18  * You should have received a copy of the GNU General Public License
  19  * along with SOAPdenovo.  If not, see <http://www.gnu.org/licenses/>.
  20  *
  21  */
  22
  23 #include "stdinc.h"
  24 #include "newhash.h"
  25 #include "kmerhash.h"
  26 #include "extfunc.h"
  27 #include "extvab.h"
  28
  29 /*************************************************
  30 Function:
  31     add1marker2edge
  32 Description:
  33     Records the id of read which crosses the edge.
  34 Input:
  35     1. edgeno:      the edge index
  36     2. readid:      the read id
  37 Output:
  38     None.
  39 Return:
  40     None.
  41 *************************************************/
  42 static void add1marker2edge ( unsigned int edgeno, long long readid )
  43 {
  44         if ( edge_array[edgeno].multi == 255 )
  45         {
  46                 return;
  47         }
  48
  49         unsigned int bal_ed = getTwinEdge ( edgeno );
  50         unsigned char counter = edge_array[edgeno].multi++;
  51         edge_array[edgeno].markers[counter] = readid;
  52         counter = edge_array[bal_ed].multi++;
  53         edge_array[bal_ed].markers[counter] = -readid;
  54 }
  55
  56 /*************************************************
  57 Function:
  58     loadPath
  59 Description:
  60     1. Loads the path info.
  61     2. Records the ids of reads crossing edges.
  62 Input:
  63     1. graphfile:       the input prefix
  64 Output:
  65     None.
  66 Return:
  67     None.
  68 *************************************************/
  69 boolean loadPath ( char * graphfile )
  70 {
  71         FILE * fp;
  72         char name[256], line[1024];
  73         unsigned int i, bal_ed, num1, edgeno, num2;
  74         long long markCounter = 0, readid = 0;
  75         char * seg;
  76         sprintf ( name, "%s.markOnEdge", graphfile );
  77         fp = fopen ( name, "r" );
  78
  79         if ( !fp )
  80         {
  81                 return 0;
  82         }
  83
  84         for ( i = 1; i <= num_ed; i++ )
  85         {
  86                 edge_array[i].multi = 0;
  87         }
  88
  89         for ( i = 1; i <= num_ed; i++ )
  90         {
  91                 fscanf ( fp, "%d", &num1 );
  92
  93                 if ( EdSmallerThanTwin ( i ) )
  94                 {
  95                         fscanf ( fp, "%d", &num2 );
  96                         bal_ed = getTwinEdge ( i );
  97
  98                         if ( num1 + num2 >= 255 )
  99                         {
 100                                 edge_array[i].multi = 255;
 101                                 edge_array[bal_ed].multi = 255;
 102                         }
 103                         else
 104                         {
 105                                 edge_array[i].multi = num1 + num2;
 106                                 edge_array[bal_ed].multi = num1 + num2;
 107                                 markCounter += 2 * ( num1 + num2 );
 108                         }
 109
 110                         i++;
 111                 }
 112                 else
 113                 {
 114                         if ( 2 * num1 >= 255 )
 115                         {
 116                                 edge_array[i].multi = 255;
 117                         }
 118                         else
 119                         {
 120                                 edge_array[i].multi = 2 * num1;
 121                                 markCounter += 2 * num1;
 122                         }
 123                 }
 124         }
 125
 126         fclose ( fp );
 127         fprintf ( stderr, "%lld markers overall.\n", markCounter );
 128         markersArray = ( long long * ) ckalloc ( markCounter * sizeof ( long long ) );
 129         markCounter = 0;
 130
 131         for ( i = 1; i <= num_ed; i++ )
 132         {
 133                 if ( edge_array[i].multi == 255 )
 134                 {
 135                         continue;
 136                 }
 137
 138                 edge_array[i].markers = markersArray + markCounter;
 139                 markCounter += edge_array[i].multi;
 140                 edge_array[i].multi = 0;
 141         }
 142
 143         sprintf ( name, "%s.path", graphfile );
 144         fp = fopen ( name, "r" );
 145
 146         if ( !fp )
 147         {
 148                 return 0;
 149         }
 150
 151         while ( fgets ( line, sizeof ( line ), fp ) != NULL )
 152         {
 153                 //printf("%s",line);
 154                 readid++;
 155                 seg = strtok ( line, " " );
 156
 157                 while ( seg )
 158                 {
 159                         edgeno = atoi ( seg );
 160                         //printf("%s, %d\n",seg,edgeno);
 161                         add1marker2edge ( edgeno, readid );
 162                         seg = strtok ( NULL, " " );
 163                 }
 164         }
 165
 166         fclose ( fp );
 167         markCounter = 0;
 168
 169         for ( i = 1; i <= num_ed; i++ )
 170         {
 171                 if ( edge_array[i].multi == 255 )
 172                 {
 173                         continue;
 174                 }
 175
 176                 markCounter += edge_array[i].multi;
 177         }
 178
 179         fprintf ( stderr, "%lld marks loaded.\n", markCounter );
 180         return 1;
 181 }
 182
 183 static int  comp ( const void * a , const void * b )
 184 {
 185         long long m , n ;
 186         m = * ( long long * ) a > 0 ? * ( long long * ) a : -* ( long long * ) a;
 187         n = * ( long long * ) b > 0 ? * ( long long * ) b : -* ( long long * ) b;
 188
 189         //      return (int)(m-n);
 190         if ( m > n )
 191                 { return 1; }
 192         else if ( m < n )
 193                 { return -1; }
 194         else
 195                 { return 0; }
 196 }
 197
 198 /*************************************************
 199 Function:
 200     loadPathBin
 201 Description:
 202     1. Loads the path info.
 203     2. Records the ids of reads crossing edges.
 204 Input:
 205     1. graphfile:       the input prefix
 206 Output:
 207     None.
 208 Return:
 209     0 if it's fail to load the path.
 210 *************************************************/
 211
 212 boolean loadPathBin ( char * graphfile )
 213 {
 214         FILE * fp;
 215         char name[256];
 216         unsigned int i, bal_ed, num1, num2;
 217         long long markCounter = 0, readid = 0;
 218         unsigned char seg, ch;
 219         unsigned int * freadBuf;
 220         sprintf ( name, "%s.markOnEdge", graphfile );
 221         fp = fopen ( name, "r" );
 222
 223         if ( !fp )
 224         {
 225                 return 0;
 226         }
 227
 228         for ( i = 1; i <= num_ed; i++ )
 229         {
 230                 edge_array[i].multi = 0;
 231                 edge_array[i].markers = NULL;
 232         }
 233
 234         for ( i = 1; i <= num_ed; i++ )
 235         {
 236                 fscanf ( fp, "%d", &num1 );
 237
 238                 if ( EdSmallerThanTwin ( i ) )
 239                 {
 240                         fscanf ( fp, "%d", &num2 );
 241                         bal_ed = getTwinEdge ( i );
 242
 243                         if ( num1 + num2 >= 255 )
 244                         {
 245                                 edge_array[i].multi = 255;
 246                                 edge_array[bal_ed].multi = 255;
 247                         }
 248                         else
 249                         {
 250                                 edge_array[i].multi = num1 + num2;
 251                                 edge_array[bal_ed].multi = num1 + num2;
 252                                 markCounter += 2 * ( num1 + num2 );
 253                         }
 254
 255                         i++;
 256                 }
 257                 else
 258                 {
 259                         if ( 2 * num1 >= 255 )
 260                         {
 261                                 edge_array[i].multi = 255;
 262                         }
 263                         else
 264                         {
 265                                 edge_array[i].multi = 2 * num1;
 266                                 markCounter += 2 * num1;
 267                         }
 268                 }
 269         }
 270
 271         fclose ( fp );
 272         fprintf ( stderr, "%lld markers overall.\n", markCounter );
 273         markersArray = ( long long * ) ckalloc ( markCounter * sizeof ( long long ) );
 274         markCounter = 0;
 275
 276         for ( i = 1; i <= num_ed; i++ )
 277         {
 278                 if ( edge_array[i].multi == 255 )
 279                 {
 280                         continue;
 281                 }
 282
 283                 edge_array[i].markers = markersArray + markCounter;
 284                 markCounter += edge_array[i].multi;
 285                 edge_array[i].multi = 0;
 286         }
 287
 288         sprintf ( name, "%s.path", graphfile );
 289         fp = fopen ( name, "rb" );
 290
 291         if ( !fp )
 292         {
 293                 return 0;
 294         }
 295
 296         freadBuf = ( unsigned int * ) ckalloc ( ( maxReadLen - overlaplen + 1 ) * sizeof ( unsigned int ) );
 297
 298         while ( fread ( &ch, sizeof ( char ), 1, fp ) == 1 )
 299         {
 300                 //printf("%s",line);
 301                 if ( fread ( freadBuf, sizeof ( unsigned int ), ch, fp ) != ch )
 302                 {
 303                         break;
 304                 }
 305
 306                 readid++;
 307
 308                 for ( seg = 0; seg < ch; seg++ )
 309                 {
 310                         add1marker2edge ( freadBuf[seg], readid );
 311                 }
 312         }
 313
 314         fclose ( fp );
 315         markCounter = 0;
 316
 317         for ( i = 1; i <= num_ed; i++ )
 318         {
 319                 if ( edge_array[i].multi == 255 )
 320                 {
 321                         continue;
 322                 }
 323
 324                 markCounter += edge_array[i].multi;
 325         }
 326
 327         for ( i = 0; i <= num_ed; i++ )
 328         {
 329                 if ( edge_array[i].multi >= 2 && edge_array[i].multi != 255 )
 330                         { qsort ( edge_array[i].markers, ( int ) edge_array[i].multi, sizeof ( long long ), comp ); }
 331         }
 332
 333         fprintf ( stderr, "%lld markers loaded.\n", markCounter );
 334         free ( ( void * ) freadBuf );
 335         return 1;
 336 }