sam_header.c

   1 #include "sam_header.h"
   2 #include <stdio.h>
   3 #include <string.h>
   4 #include <ctype.h>
   5 #include <stdlib.h>
   6 #include <stdarg.h>
   7
   8 const char *o_hd_tags[] = {"SO","GO",NULL};
   9 const char *r_hd_tags[] = {"VN",NULL};
  10 const char *o_sq_tags[] = {"AS","M5","UR","SP",NULL};
  11 const char *r_sq_tags[] = {"SN","LN",NULL};
  12 const char *o_rg_tags[] = {"LB","DS","PU","PI","CN","DT","PL",NULL};
  13 const char *r_rg_tags[] = {"ID","SM",NULL};
  14 const char *o_pg_tags[] = {"VN","CL",NULL};
  15 const char *r_pg_tags[] = {"ID",NULL};
  16 const char *types[]          = {"HD","SQ","RG","PG","CO",NULL};
  17 const char **optional_tags[] = {o_hd_tags,o_sq_tags,o_rg_tags,o_pg_tags,NULL,NULL};
  18 const char **required_tags[] = {r_hd_tags,r_sq_tags,r_rg_tags,r_pg_tags,NULL,NULL};
  19
  20 void debug(const char *format, ...)
  21 {
  22     va_list ap;
  23     va_start(ap, format);
  24     vfprintf(stderr, format, ap);
  25     va_end(ap);
  26 }
  27
  28 void error(const char *format, ...)
  29 {
  30     va_list ap;
  31     va_start(ap, format);
  32     vfprintf(stderr, format, ap);
  33     va_end(ap);
  34     exit(-1);
  35 }
  36
  37 list_t *list_append(list_t *root, void *data)
  38 {
  39     list_t *l = root;
  40     while (l && l->next)
  41         l = l->next;
  42     if ( l )
  43     {
  44         l->next = malloc(sizeof(list_t));
  45         l = l->next;
  46     }
  47     else
  48     {
  49         l = malloc(sizeof(list_t));
  50         root = l;
  51     }
  52     l->data = data;
  53     l->next = NULL;
  54     return root;
  55 }
  56
  57 void list_free(list_t *root)
  58 {
  59     list_t *l = root;
  60     while (root)
  61     {
  62         l = root;
  63         root = root->next;
  64         free(l);
  65     }
  66 }
  67
  68
  69
  70 // Look for a tag "XY" in a predefined const char *[] array.
  71 int tag_exists(const char *tag, const char **tags)
  72 {
  73     int itag=0;
  74     if ( !tags ) return -1;
  75     while ( tags[itag] )
  76     {
  77         if ( tags[itag][0]==tag[0] && tags[itag][1]==tag[1] ) return itag;
  78         itag++;
  79     }
  80     return -1;
  81 }
  82
  83
  84
  85 // Mimics the behaviour of getline, except it returns pointer to the next chunk of the text
  86 //  or NULL if everything has been read. The lineptr should be freed by the caller. The
  87 //  newline character is stripped.
  88 const char *nextline(char **lineptr, size_t *n, const char *text)
  89 {
  90     int len;
  91     const char *to = text;
  92
  93     if ( !*to ) return NULL;
  94
  95     while ( *to && *to!='\n' && *to!='\r' ) to++;
  96     len = to - text + 1;
  97
  98     if ( *to )
  99     {
 100         // Advance the pointer for the next call
 101         if ( *to=='\n' ) to++;
 102         else if ( *to=='\r' && *(to+1)=='\n' ) to+=2;
 103     }
 104     if ( !len )
 105         return to;
 106
 107     if ( !*lineptr )
 108     {
 109         *lineptr = malloc(len);
 110         *n = len;
 111     }
 112     else if ( *n<len )
 113     {
 114         *lineptr = realloc(*lineptr, len);
 115         *n = len;
 116     }
 117     if ( !*lineptr )
 118             error("FIXME\n");
 119
 120     memcpy(*lineptr,text,len);
 121     (*lineptr)[len-1] = 0;
 122
 123     return to;
 124 }
 125
 126 // name points to "XY", value_from points to the first character of the value string and
 127 //  value_to points to the last character of the value string.
 128 HeaderTag *new_tag(const char *name, const char *value_from, const char *value_to)
 129 {
 130     HeaderTag *tag = malloc(sizeof(HeaderTag));
 131     int len = value_to-value_from+1;
 132
 133     tag->key[0] = name[0];
 134     tag->key[1] = name[1];
 135     tag->value = malloc(len+1);
 136     memcpy(tag->value,value_from,len+1);
 137     tag->value[len] = 0;
 138     return tag;
 139 }
 140
 141 HeaderTag *header_line_has_tag(HeaderLine *hline, const char *key)
 142 {
 143     list_t *tags = hline->tags;
 144     while (tags)
 145     {
 146         HeaderTag *tag = tags->data;
 147         if ( tag->key[0]==key[0] && tag->key[1]==key[1] ) return tag;
 148         tags = tags->next;
 149     }
 150     return NULL;
 151 }
 152
 153 #if 0
 154 // Is there a HeaderLine with all required fields identical to those given in the hline?
 155 HeaderLine *sam_header_has_line(HeaderDict *dict, HeaderLine *hline)
 156 {
 157     HeaderLine *found=NULL;
 158
 159     while (dict)
 160     {
 161         HeaderLine *dline = dict->data;
 162
 163         if ( hline->type[0]!=dline->type[0] || hline->type[1]!=dline->type[1] )
 164         {
 165             dict = dict->next;
 166             continue;
 167         }
 168
 169         int itype = tag_exists(hline->type,types);
 170         if ( itype==-1 ) error("[sam_header_has_line] Unknown type [%c%c]\n", hline->type[0],hline->type[1]);
 171
 172         int ireq=0, differ=0;
 173         while ( required_tags[itype] && required_tags[itype][ireq] )
 174         {
 175             HeaderTag *t1, *t2;
 176             t1 = header_line_has_tag(hline,required_tags[itype][ireq]);
 177             t2 = header_line_has_tag(dline,required_tags[itype][ireq]);
 178             if ( !t1 || !t2 ) error("[sam_header_has_line] Missing a required tag [%c%c]\n",
 179                 required_tags[itype][ireq][0],required_tags[itype][ireq][1]);
 180             if ( strcmp(t1->value,t2->value) )
 181             ireq++;
 182         }
 183         dict = dict->next;
 184     }
 185     return found;
 186 }
 187 #endif
 188
 189 HeaderLine *sam_header_line_parse(const char *headerLine)
 190 {
 191     HeaderLine *hline;
 192     HeaderTag *tag;
 193     const char *from, *to;
 194     from = headerLine;
 195
 196     if ( *from != '@' ) error("[sam_header_line_parse] expected '@', got [%s]\n", headerLine);
 197     to = ++from;
 198
 199     while (*to && *to!='\t') to++;
 200     if ( to-from != 2 ) error("[sam_header_line_parse] expected '@XY', got [%s]\n", headerLine);
 201
 202     hline = malloc(sizeof(HeaderLine));
 203     hline->type[0] = from[0];
 204     hline->type[1] = from[1];
 205     hline->tags = NULL;
 206
 207     int itype = tag_exists(hline->type, types);
 208
 209     from = to;
 210     while (*to && *to=='\t') to++;
 211     if ( to-from != 1 )
 212         error("[sam_header_line_parse] multiple tabs on line [%s] (%d)\n", headerLine,(int)(to-from));
 213     from = to;
 214     while (*from)
 215     {
 216         while (*to && *to!='\t') to++;
 217
 218         if ( !required_tags[itype] && !optional_tags[itype] )
 219             tag = new_tag("  ",from,to-1);
 220         else
 221             tag = new_tag(from,from+3,to-1);
 222
 223         if ( header_line_has_tag(hline,tag->key) )
 224                 debug("The tag '%c%c' present (at least) twice on line [%s]\n", tag->key[0],tag->key[1], headerLine);
 225         hline->tags = list_append(hline->tags, tag);
 226
 227         from = to;
 228         while (*to && *to=='\t') to++;
 229         if ( *to && to-from != 1 )
 230                 error("[sam_header_line_parse] multiple tabs on line [%s] (%d)\n", headerLine,(int)(to-from));
 231
 232         from = to;
 233     }
 234     return hline;
 235 }
 236
 237
 238 // Must be of an existing type, all tags must be recognised and all required tags must be present
 239 int sam_header_line_validate(HeaderLine *hline)
 240 {
 241     list_t *tags;
 242     HeaderTag *tag;
 243     int itype, itag;
 244
 245     // Is the type correct?
 246     itype = tag_exists(hline->type, types);
 247     if ( itype==-1 )
 248     {
 249         debug("The type [%c%c] not recognised.\n", hline->type[0],hline->type[1]);
 250         return 0;
 251     }
 252
 253     // Has all required tags?
 254     itag = 0;
 255     while ( required_tags[itype] && required_tags[itype][itag] )
 256     {
 257         if ( !header_line_has_tag(hline,required_tags[itype][itag]) )
 258         {
 259             debug("The tag [%c%c] required for [%c%c] not present.\n", required_tags[itype][itag][0],required_tags[itype][itag][1],
 260                 hline->type[0],hline->type[1]);
 261             return 0;
 262         }
 263         itag++;
 264     }
 265
 266     // Are all tags recognised?
 267     tags = hline->tags;
 268     while ( tags )
 269     {
 270         tag = tags->data;
 271         if ( !tag_exists(tag->key,required_tags[itype]) && !tag_exists(tag->key,optional_tags[itype]) )
 272         {
 273             debug("Unknown tag [%c%c] for [%c%c].\n", tag->key[0],tag->key[1], hline->type[0],hline->type[1]);
 274             return 0;
 275         }
 276         tags = tags->next;
 277     }
 278
 279     return 1;
 280 }
 281
 282 void print_header_line(HeaderLine *hline)
 283 {
 284     list_t *tags = hline->tags;
 285     HeaderTag *tag;
 286
 287     printf("@%c%c", hline->type[0],hline->type[1]);
 288     while (tags)
 289     {
 290         tag = tags->data;
 291         printf("\t%c%c:%s", tag->key[0],tag->key[1],tag->value);
 292         tags = tags->next;
 293     }
 294     printf("\n");
 295 }
 296
 297
 298 void sam_header_free(HeaderDict *header)
 299 {
 300     list_t *hlines = header;
 301     while (hlines)
 302     {
 303         HeaderLine *hline = hlines->data;
 304         list_t *tags = hline->tags;
 305         while (tags)
 306         {
 307             HeaderTag *tag = tags->data;
 308             free(tag->value);
 309             free(tag);
 310             tags = tags->next;
 311         }
 312         list_free(hline->tags);
 313         free(hline);
 314         hlines = hlines->next;
 315     }
 316     list_free(header);
 317 }
 318
 319 // Returns a newly allocated string
 320 char *sam_header_write(const HeaderDict *header)
 321 {
 322     char *out = NULL;
 323     int len=0, nout=0;
 324     const list_t *hlines;
 325
 326     // Calculate the length of the string to allocate
 327     hlines = header;
 328     while (hlines)
 329     {
 330         len += 4;   // @XY and \n
 331
 332         HeaderLine *hline = hlines->data;
 333         list_t *tags = hline->tags;
 334         while (tags)
 335         {
 336             HeaderTag *tag = tags->data;
 337             len += strlen(tag->value) + 1;                  // \t
 338             if ( tag->key[0]!=' ' || tag->key[1]!=' ' )
 339                 len += strlen(tag->value) + 3;              // XY:
 340             tags = tags->next;
 341         }
 342         hlines = hlines->next;
 343     }
 344
 345     nout = 0;
 346     out  = malloc(len+1);
 347     hlines = header;
 348     while (hlines)
 349     {
 350         HeaderLine *hline = hlines->data;
 351
 352         nout += sprintf(out+nout,"@%c%c",hline->type[0],hline->type[1]);
 353
 354         list_t *tags = hline->tags;
 355         while (tags)
 356         {
 357             HeaderTag *tag = tags->data;
 358             nout += sprintf(out+nout,"\t");
 359             if ( tag->key[0]!=' ' || tag->key[1]!=' ' )
 360                 nout += sprintf(out+nout,"%c%c:", tag->key[0],tag->key[1]);
 361             nout += sprintf(out+nout,"%s", tag->value);
 362             tags = tags->next;
 363         }
 364         hlines = hlines->next;
 365         nout += sprintf(out+nout,"\n");
 366     }
 367     out[len] = 0;
 368     return out;
 369 }
 370
 371 HeaderDict *sam_header_parse(const char *headerText)
 372 {
 373     list_t *hlines = NULL;
 374     HeaderLine *hline;
 375     const char *text;
 376     char *buf=NULL;
 377     size_t nbuf = 0;
 378
 379     if ( !headerText )
 380         error("FIXME");
 381
 382     text = headerText;
 383     while ( (text=nextline(&buf, &nbuf, text)) )
 384     {
 385         hline = sam_header_line_parse(buf);
 386         if ( sam_header_line_validate(hline) )
 387             hlines = list_append(hlines, hline);
 388         else
 389         {
 390             sam_header_free(hlines);
 391             return NULL;
 392         }
 393     }
 394     if ( buf ) free(buf);
 395
 396     return hlines;
 397 }
 398
 399 khash_t(str) *sam_header_lookup_table(const HeaderDict *dict, char type[2], char key_tag[2], char value_tag[2])
 400 {
 401     const list_t *l   = dict;
 402     khash_t(str) *tbl = kh_init(str);
 403     khiter_t k;
 404     int ret;
 405
 406     while (l)
 407     {
 408         HeaderLine *hline = l->data;
 409         if ( hline->type[0]!=type[0] || hline->type[1]!=type[1] )
 410         {
 411             l = l->next;
 412             continue;
 413         }
 414
 415         HeaderTag *key, *value;
 416         key   = header_line_has_tag(hline,key_tag);
 417         value = header_line_has_tag(hline,value_tag);
 418         if ( !key || !value )
 419         {
 420             l = l->next;
 421             continue;
 422         }
 423
 424         k = kh_get(str, tbl, key->value);
 425         if ( k != kh_end(tbl) )
 426             debug("[sam_header_lookup_table] They key %s not unique.\n", key->value);
 427         k = kh_put(str, tbl, key->value, &ret);
 428         kh_value(tbl, k) = value->value;
 429
 430         l = l->next;
 431     }
 432     return tbl;
 433 }
 434
 435
 436 #if 0
 437 TODO
 438 HeaderDict *sam_header_merge(int n, const HeaderDict **dicts)
 439 {
 440     HeaderDict *out=NULL;
 441     int idict;
 442
 443     for (idict=0; idict<n; idict++)
 444     {
 445         list_t *hlines = dicts[idict];
 446         while (hlines)
 447         {
 448             HeaderLine *hline = sam_header_has_line(out, hlines->data);
 449             sam_header_line_merge(hline,hlines->data);
 450             hlines = hlines->next;
 451         }
 452     }
 453 }
 454 #endif
 455