code_c/Maasha/src/lib/ucsc.c

   1 /* Martin Asser Hansen (mail@maasha.dk) Copyright (C) 2008 - All right reserved */
   2
   3
   4 #include "common.h"
   5 #include "mem.h"
   6 #include "filesys.h"
   7 #include "list.h"
   8 #include "strings.h"
   9 #include "ucsc.h"
  10
  11
  12 bed_entry *bed_entry_new( const int cols )
  13 {
  14     /* Martin A. Hansen, September 2008 */
  15
  16     /* Returns a new BED entry with memory allocated for */
  17     /* a given number of columns. */
  18
  19     bed_entry *entry = mem_get( sizeof( bed_entry ) );
  20
  21     entry->cols    = cols;
  22     entry->chr     = mem_get( BED_CHR_MAX );
  23     entry->chr_beg = 0;
  24     entry->chr_end = 0;
  25
  26     if ( cols == 3 ) {
  27         return entry;
  28     }
  29
  30     entry->q_id = mem_get( BED_QID_MAX );
  31
  32     if ( cols == 4 ) {
  33         return entry;
  34     }
  35
  36     entry->score = 0;
  37
  38     if ( cols == 5 ) {
  39         return entry;
  40     }
  41
  42     entry->strand = 0;
  43
  44     if ( cols == 6 ) {
  45         return entry;
  46     }
  47
  48     entry->thick_beg  = 0;
  49     entry->thick_end  = 0;
  50     entry->itemrgb    = mem_get( BED_ITEMRGB_MAX );
  51     entry->blockcount = 0;
  52     entry->blocksizes = mem_get( BED_BLOCKSIZES_MAX );
  53     entry->q_begs     = mem_get( BED_QBEGS_MAX );
  54
  55     return entry;
  56 }
  57
  58
  59 void bed_entry_destroy( bed_entry *entry )
  60 {
  61     /* Martin A. Hansen, September 2008 */
  62
  63     /* Free memory for a BED entry. */
  64
  65     int cols = entry->cols;
  66
  67     if ( cols > 6 )
  68     {
  69         free( entry->itemrgb );
  70         free( entry->blocksizes );
  71         free( entry->q_begs );
  72         free( entry->q_id );
  73         free( entry->chr );
  74     }
  75     else if ( cols > 3 )
  76     {
  77         free( entry->q_id );
  78         free( entry->chr );
  79     }
  80     else
  81     {
  82         free( entry->chr );
  83     }
  84
  85     mem_free( &entry );
  86 }
  87
  88
  89 bed_entry *bed_entry_get( FILE *fp, int cols )
  90 {
  91     /* Martin A. Hansen, September 2008 */
  92
  93     /* Get next BED entry of a given number of columns from a file pointer. */
  94
  95     bed_entry *entry = bed_entry_new( cols );
  96     char       buffer[ BED_BUFFER ];
  97
  98     assert( cols == 0 || cols == 3 || cols == 4 || cols == 5 || cols == 6 || cols == 12 );
  99
 100     if ( fgets( buffer, sizeof( buffer ), fp ) != NULL )
 101     {
 102         if ( ! cols )
 103         {
 104             cols = 1 + strchr_total( buffer, '\t' );
 105             entry->cols = cols;
 106         }
 107
 108         if ( cols == 3 )
 109         {
 110             sscanf(
 111                     buffer,
 112                     "%s\t%u\t%u",
 113                     entry->chr,
 114                     &entry->chr_beg,
 115                     &entry->chr_end
 116             );
 117
 118             return entry;
 119         }
 120
 121         if ( cols == 4 )
 122         {
 123             sscanf(
 124                     buffer,
 125                     "%s\t%u\t%u\t%s",
 126                     entry->chr,
 127                     &entry->chr_beg,
 128                     &entry->chr_end,
 129                     entry->q_id
 130             );
 131
 132             return entry;
 133         }
 134
 135         if ( cols == 5 )
 136         {
 137             sscanf(
 138                     buffer,
 139                     "%s\t%u\t%u\t%s\t%i",
 140                     entry->chr,
 141                     &entry->chr_beg,
 142                     &entry->chr_end,
 143                     entry->q_id,
 144                     &entry->score
 145             );
 146
 147             return entry;
 148         }
 149
 150         if ( cols == 6 )
 151         {
 152             sscanf(
 153                     buffer,
 154                     "%s\t%u\t%u\t%s\t%i\t%c",
 155                     entry->chr,
 156                     &entry->chr_beg,
 157                     &entry->chr_end,
 158                     entry->q_id,
 159                     &entry->score,
 160                     &entry->strand
 161             );
 162
 163             return entry;
 164         }
 165
 166         if ( cols == 12 )
 167         {
 168             sscanf(
 169                     buffer,
 170                     "%s\t%u\t%u\t%s\t%i\t%c\t%u\t%u\t%s\t%u\t%s\t%s",
 171                     entry->chr,
 172                     &entry->chr_beg,
 173                     &entry->chr_end,
 174                     entry->q_id,
 175                     &entry->score,
 176                     &entry->strand,
 177                     &entry->thick_beg,
 178                     &entry->thick_end,
 179                     entry->itemrgb,
 180                     &entry->blockcount,
 181                     entry->blocksizes,
 182                     entry->q_begs
 183             );
 184
 185             return entry;
 186         }
 187     }
 188
 189     return NULL;
 190 }
 191
 192
 193 list_sl *bed_entries_get( char *path, const int cols )
 194 {
 195     /* Martin A. Hansen, September 2008 */
 196
 197     /* Get a singly linked list with all BED entries (of a given number of coluns */
 198     /* from a specified file. */
 199
 200     list_sl   *list     = list_sl_new();
 201     node_sl   *node     = node_sl_new();
 202     node_sl   *old_node = NULL;
 203     bed_entry *entry    = NULL;
 204     FILE      *fp       = NULL;
 205
 206     fp = read_open( path );
 207
 208     if ( ( entry = bed_entry_get( fp, cols ) ) != NULL )
 209     {
 210         node->val = entry;
 211
 212         list_sl_add_beg( &list, &node );
 213
 214         old_node = node;
 215     }
 216
 217     while ( ( entry = bed_entry_get( fp, cols ) ) != NULL )
 218     {
 219         node = node_sl_new();
 220
 221         node->val = entry;
 222
 223         list_sl_add_after( &old_node, &node );
 224
 225         old_node = node;
 226     }
 227
 228     close_stream( fp );
 229
 230     return list;
 231 }
 232
 233
 234 void bed_entry_put( bed_entry *entry, int cols )
 235 {
 236     /* Martin A. Hansen, September 2008 */
 237
 238     /* Output a given number of columns from a BED entry to stdout. */
 239
 240     if ( ! cols ) {
 241         cols = entry->cols;
 242     }
 243
 244     if ( cols == 3 )
 245     {
 246         printf(
 247             "%s\t%u\t%u\n",
 248             entry->chr,
 249             entry->chr_beg,
 250             entry->chr_end
 251         );
 252     }
 253     else if ( cols == 4 )
 254     {
 255         printf(
 256             "%s\t%u\t%u\t%s\n",
 257             entry->chr,
 258             entry->chr_beg,
 259             entry->chr_end,
 260             entry->q_id
 261         );
 262     }
 263     else if ( cols == 5 )
 264     {
 265         printf(
 266             "%s\t%u\t%u\t%s\t%i\n",
 267             entry->chr,
 268             entry->chr_beg,
 269             entry->chr_end,
 270             entry->q_id,
 271             entry->score
 272         );
 273     }
 274     else if ( cols == 6 )
 275     {
 276         printf(
 277             "%s\t%u\t%u\t%s\t%i\t%c\n",
 278             entry->chr,
 279             entry->chr_beg,
 280             entry->chr_end,
 281             entry->q_id,
 282             entry->score,
 283             entry->strand
 284         );
 285     }
 286     else if ( cols == 12 )
 287     {
 288         printf(
 289             "%s\t%u\t%u\t%s\t%i\t%c\t%u\t%u\t%s\t%u\t%s\t%s\n",
 290             entry->chr,
 291             entry->chr_beg,
 292             entry->chr_end,
 293             entry->q_id,
 294             entry->score,
 295             entry->strand,
 296             entry->thick_beg,
 297             entry->thick_end,
 298             entry->itemrgb,
 299             entry->blockcount,
 300             entry->blocksizes,
 301             entry->q_begs
 302         );
 303     }
 304     else
 305     {
 306         fprintf( stderr, "ERROR: Wrong number of columns in bed_entry_put: %d\n", cols );
 307
 308         abort();
 309     }
 310 }
 311
 312
 313 void bed_entries_put( list_sl *entries, int cols )
 314 {
 315     /* Martin A. Hansen, September 2008 */
 316
 317     /* Output a given number of columns from all BED entries */
 318     /* in a singly linked list. */
 319
 320     node_sl *node = NULL;
 321
 322     for ( node = entries->first; node != NULL; node = node->next ) {
 323         bed_entry_put( ( bed_entry * ) node->val, cols );
 324     }
 325 }
 326
 327
 328 void bed_entries_destroy( list_sl **entries_ppt )
 329 {
 330     /* Martin A. Hansen, September 2008 */
 331
 332     /* Free memory for all BED entries and list nodes. */
 333
 334     list_sl *entries = *entries_ppt;
 335     node_sl *node    = NULL;
 336     node_sl *next    = NULL;
 337
 338     next = entries->first;
 339
 340     while ( next != NULL )
 341     {
 342         node = next;
 343
 344         bed_entry_destroy( ( bed_entry * ) node->val );
 345
 346         next = node->next;
 347
 348         mem_free( &node );
 349     }
 350
 351     mem_free( &entries );
 352
 353     *entries_ppt = NULL;
 354 }
 355
 356
 357 void bed_file_sort_beg( char *path, int cols )
 358 {
 359     /* Martin A. Hansen, September 2008 */
 360
 361     /* Given a path to a BED file, read the given number of cols */
 362     /* according to the begin position. The result is written to stdout. */
 363
 364     list_sl *entries = NULL;
 365
 366     entries = bed_entries_get( path, cols );
 367
 368     list_sl_sort( &entries, cmp_bed_sort_beg );
 369
 370     bed_entries_put( entries, cols );
 371
 372     bed_entries_destroy( &entries );
 373 }
 374
 375
 376 void bed_file_sort_strand_beg( char *path, int cols )
 377 {
 378     /* Martin A. Hansen, September 2008 */
 379
 380     /* Given a path to a BED file, read the given number of cols */
 381     /* according to the strand AND begin position. The result is written to stdout. */
 382
 383     assert( cols >= 6 );
 384
 385     list_sl *entries = NULL;
 386
 387     entries = bed_entries_get( path, cols );
 388
 389     list_sl_sort( &entries, cmp_bed_sort_strand_beg );
 390
 391     bed_entries_put( entries, cols );
 392
 393     bed_entries_destroy( &entries );
 394 }
 395
 396
 397 void bed_file_sort_chr_beg( char *path, int cols )
 398 {
 399     /* Martin A. Hansen, September 2008 */
 400
 401     /* Given a path to a BED file, read the given number of cols */
 402     /* according to the chromosome AND begin position. The result is written to stdout. */
 403
 404     list_sl *entries = NULL;
 405
 406     entries = bed_entries_get( path, cols );
 407
 408     list_sl_sort( &entries, cmp_bed_sort_chr_beg );
 409
 410     bed_entries_put( entries, cols );
 411
 412     bed_entries_destroy( &entries );
 413 }
 414
 415
 416 void bed_file_sort_chr_strand_beg( char *path, int cols )
 417 {
 418     /* Martin A. Hansen, September 2008 */
 419
 420     /* Given a path to a BED file, read the given number of cols */
 421     /* according to the chromosome AND strand AND begin position. The result is written to stdout. */
 422
 423     assert( cols >= 6 );
 424
 425     list_sl *entries = NULL;
 426
 427     entries = bed_entries_get( path, cols );
 428
 429     list_sl_sort( &entries, cmp_bed_sort_chr_strand_beg );
 430
 431     bed_entries_put( entries, cols );
 432
 433     bed_entries_destroy( &entries );
 434 }
 435
 436
 437 int cmp_bed_sort_beg( const void *a, const void *b )
 438 {
 439     /* Martin A. Hansen, September 2008 */
 440
 441     /* Compare function for sorting a singly linked list of BED entries */
 442     /* according to begin position. */
 443
 444     node_sl *a_node = *( ( node_sl ** ) a );
 445     node_sl *b_node = *( ( node_sl ** ) b );
 446
 447     bed_entry *a_entry = ( bed_entry * ) a_node->val;
 448     bed_entry *b_entry = ( bed_entry * ) b_node->val;
 449
 450     if ( a_entry->chr_beg < b_entry->chr_beg ) {
 451         return -1;
 452     } else if ( a_entry->chr_beg > b_entry->chr_beg ) {
 453         return 1;
 454     } else {
 455         return 0;
 456     }
 457 }
 458
 459
 460 int cmp_bed_sort_strand_beg( const void *a, const void *b )
 461 {
 462     /* Martin A. Hansen, September 2008 */
 463
 464     /* Compare function for sorting a singly linked list of BED entries */
 465     /* according to strand AND begin position. */
 466
 467     node_sl *a_node = *( ( node_sl ** ) a );
 468     node_sl *b_node = *( ( node_sl ** ) b );
 469
 470     bed_entry *a_entry = ( bed_entry * ) a_node->val;
 471     bed_entry *b_entry = ( bed_entry * ) b_node->val;
 472
 473     if ( a_entry->strand < b_entry->strand ) {
 474         return -1;
 475     } else if ( a_entry->strand > b_entry->strand ) {
 476         return 1;
 477     } else if ( a_entry->chr_beg < b_entry->chr_beg ) {
 478         return -1;
 479     } else if ( a_entry->chr_beg > b_entry->chr_beg ) {
 480         return 1;
 481     } else {
 482         return 0;
 483     }
 484 }
 485
 486
 487 int cmp_bed_sort_chr_beg( const void *a, const void *b )
 488 {
 489     /* Martin A. Hansen, September 2008 */
 490
 491     /* Compare function for sorting a singly linked list of BED entries */
 492     /* according to chromosome name AND begin position. */
 493
 494     node_sl *a_node = *( ( node_sl ** ) a );
 495     node_sl *b_node = *( ( node_sl ** ) b );
 496
 497     bed_entry *a_entry = ( bed_entry * ) a_node->val;
 498     bed_entry *b_entry = ( bed_entry * ) b_node->val;
 499
 500     int diff = 0;
 501
 502     diff = strcmp( a_entry->chr, b_entry->chr );
 503
 504     if ( diff < 0 ) {
 505         return -1;
 506     } else if ( diff > 0 ) {
 507         return 1;
 508     } else if ( a_entry->chr_beg < b_entry->chr_beg ) {
 509         return -1;
 510     } else if ( a_entry->chr_beg > b_entry->chr_beg ) {
 511         return 1;
 512     } else {
 513         return 0;
 514     }
 515 }
 516
 517
 518 int cmp_bed_sort_chr_strand_beg( const void *a, const void *b )
 519 {
 520     /* Martin A. Hansen, September 2008 */
 521
 522     /* Compare function for sorting a singly linked list of BED entries */
 523     /* according to chromosome name AND strand AND begin position. */
 524
 525     node_sl *a_node = *( ( node_sl ** ) a );
 526     node_sl *b_node = *( ( node_sl ** ) b );
 527
 528     bed_entry *a_entry = ( bed_entry * ) a_node->val;
 529     bed_entry *b_entry = ( bed_entry * ) b_node->val;
 530
 531     int diff = 0;
 532
 533     diff = strcmp( a_entry->chr, b_entry->chr );
 534
 535     if ( diff < 0 ) {
 536         return -1;
 537     } else if ( diff > 0 ) {
 538         return 1;
 539     } else if ( a_entry->strand < b_entry->strand ) {
 540         return -1;
 541     } else if ( a_entry->strand > b_entry->strand ) {
 542         return 1;
 543     } else if ( a_entry->chr_beg < b_entry->chr_beg ) {
 544         return -1;
 545     } else if ( a_entry->chr_beg > b_entry->chr_beg ) {
 546         return 1;
 547     } else {
 548         return 0;
 549     }
 550 }
 551