src/api/internal/SamHeaderValidator_p.cpp

   1 // ***************************************************************************
   2 // SamHeaderValidator.cpp (c) 2010 Derek Barnett
   3 // Marth Lab, Department of Biology, Boston College
   4 // ---------------------------------------------------------------------------
   5 // Last modified: 18 April 2011 (DB)
   6 // ---------------------------------------------------------------------------
   7 // Provides functionality for validating SamHeader data
   8 // ***************************************************************************
   9
  10 #include <api/SamConstants.h>
  11 #include <api/SamHeader.h>
  12 #include <api/internal/SamHeaderValidator_p.h>
  13 #include <api/internal/SamHeaderVersion_p.h>
  14 using namespace BamTools;
  15 using namespace BamTools::Internal;
  16
  17 #include <cctype>
  18 #include <iostream>
  19 #include <set>
  20 #include <sstream>
  21 using namespace std;
  22
  23 namespace BamTools {
  24 namespace Internal {
  25
  26 bool caseInsensitiveCompare(const string& lhs, const string& rhs) {
  27
  28     // can omit checking chars if lengths not equal
  29     const int lhsLength = lhs.length();
  30     const int rhsLength = rhs.length();
  31     if ( lhsLength != rhsLength )
  32         return false;
  33
  34     // do *basic* toupper checks on each string char's
  35     for ( int i = 0; i < lhsLength; ++i ) {
  36         if ( toupper( (int)lhs.at(i)) != toupper( (int)rhs.at(i)) )
  37             return false;
  38     }
  39
  40     // otherwise OK
  41     return true;
  42 }
  43
  44 } // namespace Internal
  45 } // namespace BamTools
  46
  47 // ------------------------------------------------------------------------
  48 // Allow validation rules to vary, as needed, between SAM header versions
  49 //
  50 // use SAM_VERSION_X_Y to tag important changes
  51 //
  52 // Together, they will allow for comparisons like:
  53 // if ( m_version < SAM_VERSION_2_0 ) {
  54 //     // use some older rule
  55 // else
  56 //     // use rule introduced with version 2.0
  57
  58 static const SamHeaderVersion SAM_VERSION_1_0 = SamHeaderVersion(1,0);
  59 static const SamHeaderVersion SAM_VERSION_1_1 = SamHeaderVersion(1,1);
  60 static const SamHeaderVersion SAM_VERSION_1_2 = SamHeaderVersion(1,2);
  61 static const SamHeaderVersion SAM_VERSION_1_3 = SamHeaderVersion(1,3);
  62 static const SamHeaderVersion SAM_VERSION_1_4 = SamHeaderVersion(1,4);
  63
  64 // TODO: This functionality is currently unused.
  65 //       Make validation "version-aware."
  66 //
  67 // ------------------------------------------------------------------------
  68
  69 const string SamHeaderValidator::ERROR_PREFIX = "ERROR: ";
  70 const string SamHeaderValidator::WARN_PREFIX  = "WARNING: ";
  71 const string SamHeaderValidator::NEWLINE      = "\n";
  72
  73 SamHeaderValidator::SamHeaderValidator(const SamHeader& header)
  74     : m_header(header)
  75 { }
  76
  77 SamHeaderValidator::~SamHeaderValidator(void) { }
  78
  79 bool SamHeaderValidator::Validate(bool verbose) {
  80
  81     // validate header components
  82     bool isValid = true;
  83     isValid &= ValidateMetadata();
  84     isValid &= ValidateSequenceDictionary();
  85     isValid &= ValidateReadGroupDictionary();
  86     isValid &= ValidateProgramChain();
  87
  88     // report errors if desired
  89     if ( verbose ) {
  90         PrintErrorMessages();
  91         PrintWarningMessages();
  92     }
  93
  94     // return validation status
  95     return isValid;
  96 }
  97
  98 bool SamHeaderValidator::ValidateMetadata(void) {
  99     bool isValid = true;
 100     isValid &= ValidateVersion();
 101     isValid &= ValidateSortOrder();
 102     isValid &= ValidateGroupOrder();
 103     return isValid;
 104 }
 105
 106 bool SamHeaderValidator::ValidateVersion(void) {
 107
 108     const string& version = m_header.Version;
 109
 110     // warn if version not present
 111     if ( version.empty() ) {
 112         AddWarning("Version (VN) missing. Not required, but strongly recommended");
 113         return true;
 114     }
 115
 116     // invalid if version does not contain a period
 117     const size_t periodFound = version.find(Constants::SAM_PERIOD);
 118     if ( periodFound == string::npos ) {
 119         AddError("Invalid version (VN) format: " + version);
 120         return false;
 121     }
 122
 123     // invalid if major version is empty or contains non-digits
 124     const string majorVersion = version.substr(0, periodFound);
 125     if ( majorVersion.empty() || !ContainsOnlyDigits(majorVersion) ) {
 126         AddError("Invalid version (VN) format: " + version);
 127         return false;
 128     }
 129
 130     // invalid if major version is empty or contains non-digits
 131     const string minorVersion = version.substr(periodFound + 1);
 132     if ( minorVersion.empty() || !ContainsOnlyDigits(minorVersion) ) {
 133         AddError("Invalid version (VN) format: " + version);
 134         return false;
 135     }
 136
 137     // TODO: check if version is not just syntactically OK,
 138     // but is also a valid SAM version ( 1.0 .. CURRENT )
 139
 140     // all checked out this far, then version is OK
 141     return true;
 142 }
 143
 144 // assumes non-empty input string
 145 bool SamHeaderValidator::ContainsOnlyDigits(const string& s) {
 146     const size_t nonDigitPosition = s.find_first_not_of(Constants::SAM_DIGITS);
 147     return ( nonDigitPosition == string::npos ) ;
 148 }
 149
 150 bool SamHeaderValidator::ValidateSortOrder(void) {
 151
 152     const string& sortOrder = m_header.SortOrder;
 153
 154     // warn if sort order not present
 155     if ( sortOrder.empty() ) {
 156         AddWarning("Sort order (SO) missing. Not required, but strongly recommended");
 157         return true;
 158     }
 159
 160     // if sort order is valid keyword
 161     if ( sortOrder == Constants::SAM_HD_SORTORDER_COORDINATE ||
 162          sortOrder == Constants::SAM_HD_SORTORDER_QUERYNAME  ||
 163          sortOrder == Constants::SAM_HD_SORTORDER_UNSORTED
 164        )
 165     {
 166         return true;
 167     }
 168
 169     // otherwise
 170     AddError("Invalid sort order (SO): " + sortOrder);
 171     return false;
 172 }
 173
 174 bool SamHeaderValidator::ValidateGroupOrder(void) {
 175
 176     const string& groupOrder = m_header.GroupOrder;
 177
 178     // if no group order, no problem, just return OK
 179     if ( groupOrder.empty() )
 180         return true;
 181
 182     // if group order is valid keyword
 183     if ( groupOrder == Constants::SAM_HD_GROUPORDER_NONE  ||
 184          groupOrder == Constants::SAM_HD_GROUPORDER_QUERY ||
 185          groupOrder == Constants::SAM_HD_GROUPORDER_REFERENCE
 186        )
 187     {
 188         return true;
 189     }
 190
 191     // otherwise
 192     AddError("Invalid group order (GO): " + groupOrder);
 193     return false;
 194 }
 195
 196 bool SamHeaderValidator::ValidateSequenceDictionary(void) {
 197
 198     bool isValid = true;
 199
 200     // check for unique sequence names
 201     isValid &= ContainsUniqueSequenceNames();
 202
 203     // iterate over sequences
 204     const SamSequenceDictionary& sequences = m_header.Sequences;
 205     SamSequenceConstIterator seqIter = sequences.ConstBegin();
 206     SamSequenceConstIterator seqEnd  = sequences.ConstEnd();
 207     for ( ; seqIter != seqEnd; ++seqIter ) {
 208         const SamSequence& seq = (*seqIter);
 209         isValid &= ValidateSequence(seq);
 210     }
 211
 212     // return validation state
 213     return isValid;
 214 }
 215
 216 bool SamHeaderValidator::ContainsUniqueSequenceNames(void) {
 217
 218     bool isValid = true;
 219     set<string> sequenceNames;
 220     set<string>::iterator nameIter;
 221
 222     // iterate over sequences
 223     const SamSequenceDictionary& sequences = m_header.Sequences;
 224     SamSequenceConstIterator seqIter = sequences.ConstBegin();
 225     SamSequenceConstIterator seqEnd  = sequences.ConstEnd();
 226     for ( ; seqIter != seqEnd; ++seqIter ) {
 227         const SamSequence& seq = (*seqIter);
 228
 229         // lookup sequence name
 230         const string& name = seq.Name;
 231         nameIter = sequenceNames.find(name);
 232
 233         // error if found (duplicate entry)
 234         if ( nameIter != sequenceNames.end() ) {
 235             AddError("Sequence name (SN): " + name + " is not unique");
 236             isValid = false;
 237         }
 238
 239         // otherwise ok, store name
 240         sequenceNames.insert(name);
 241     }
 242
 243     // return validation state
 244     return isValid;
 245 }
 246
 247 bool SamHeaderValidator::ValidateSequence(const SamSequence& seq) {
 248     bool isValid = true;
 249     isValid &= CheckNameFormat(seq.Name);
 250     isValid &= CheckLengthInRange(seq.Length);
 251     return isValid;
 252 }
 253
 254 bool SamHeaderValidator::CheckNameFormat(const string& name) {
 255
 256     // invalid if name is empty
 257     if ( name.empty() ) {
 258         AddError("Sequence entry (@SQ) is missing SN tag");
 259         return false;
 260     }
 261
 262     // invalid if first character is a reserved char
 263     const char firstChar = name.at(0);
 264     if ( firstChar == Constants::SAM_EQUAL || firstChar == Constants::SAM_STAR ) {
 265         AddError("Invalid sequence name (SN): " + name);
 266         return false;
 267     }
 268     // otherwise OK
 269     return true;
 270 }
 271
 272 bool SamHeaderValidator::CheckLengthInRange(const string& length) {
 273
 274     // invalid if empty
 275     if ( length.empty() ) {
 276         AddError("Sequence entry (@SQ) is missing LN tag");
 277         return false;
 278     }
 279
 280     // convert string length to numeric
 281     stringstream lengthStream(length);
 282     unsigned int sequenceLength;
 283     lengthStream >> sequenceLength;
 284
 285     // invalid if length outside accepted range
 286     if ( sequenceLength < Constants::SAM_SQ_LENGTH_MIN || sequenceLength > Constants::SAM_SQ_LENGTH_MAX ) {
 287         AddError("Sequence length (LN): " + length + " out of range");
 288         return false;
 289     }
 290
 291     // otherwise OK
 292     return true;
 293 }
 294
 295 bool SamHeaderValidator::ValidateReadGroupDictionary(void) {
 296
 297     bool isValid = true;
 298
 299     // check for unique read group IDs & platform units
 300     isValid &= ContainsUniqueIDsAndPlatformUnits();
 301
 302     // iterate over read groups
 303     const SamReadGroupDictionary& readGroups = m_header.ReadGroups;
 304     SamReadGroupConstIterator rgIter = readGroups.ConstBegin();
 305     SamReadGroupConstIterator rgEnd  = readGroups.ConstEnd();
 306     for ( ; rgIter != rgEnd; ++rgIter ) {
 307         const SamReadGroup& rg = (*rgIter);
 308         isValid &= ValidateReadGroup(rg);
 309     }
 310
 311     // return validation state
 312     return isValid;
 313 }
 314
 315 bool SamHeaderValidator::ContainsUniqueIDsAndPlatformUnits(void) {
 316
 317     bool isValid = true;
 318     set<string> readGroupIds;
 319     set<string> platformUnits;
 320     set<string>::iterator idIter;
 321     set<string>::iterator puIter;
 322
 323     // iterate over sequences
 324     const SamReadGroupDictionary& readGroups = m_header.ReadGroups;
 325     SamReadGroupConstIterator rgIter = readGroups.ConstBegin();
 326     SamReadGroupConstIterator rgEnd  = readGroups.ConstEnd();
 327     for ( ; rgIter != rgEnd; ++rgIter ) {
 328         const SamReadGroup& rg = (*rgIter);
 329
 330         // --------------------------------
 331         // check for unique ID
 332
 333         // lookup read group ID
 334         const string& id = rg.ID;
 335         idIter = readGroupIds.find(id);
 336
 337         // error if found (duplicate entry)
 338         if ( idIter != readGroupIds.end() ) {
 339             AddError("Read group ID (ID): " + id + " is not unique");
 340             isValid = false;
 341         }
 342
 343         // otherwise ok, store id
 344         readGroupIds.insert(id);
 345
 346         // --------------------------------
 347         // check for unique platform unit
 348
 349         // lookup platform unit
 350         const string& pu = rg.PlatformUnit;
 351         puIter = platformUnits.find(pu);
 352
 353         // error if found (duplicate entry)
 354         if ( puIter != platformUnits.end() ) {
 355             AddError("Platform unit (PU): " + pu + " is not unique");
 356             isValid = false;
 357         }
 358
 359         // otherwise ok, store platform unit
 360         platformUnits.insert(pu);
 361     }
 362
 363     // return validation state
 364     return isValid;
 365 }
 366
 367 bool SamHeaderValidator::ValidateReadGroup(const SamReadGroup& rg) {
 368     bool isValid = true;
 369     isValid &= CheckReadGroupID(rg.ID);
 370     isValid &= CheckSequencingTechnology(rg.SequencingTechnology);
 371     return isValid;
 372 }
 373
 374 bool SamHeaderValidator::CheckReadGroupID(const string& id) {
 375
 376     // invalid if empty
 377     if ( id.empty() ) {
 378         AddError("Read group entry (@RG) is missing ID tag");
 379         return false;
 380     }
 381
 382     // otherwise OK
 383     return true;
 384 }
 385
 386 bool SamHeaderValidator::CheckSequencingTechnology(const string& technology) {
 387
 388     // if no technology provided, no problem, just return OK
 389     if ( technology.empty() )
 390         return true;
 391
 392     // if technology is valid keyword
 393     if ( caseInsensitiveCompare(technology, Constants::SAM_RG_SEQTECHNOLOGY_CAPILLARY)  ||
 394          caseInsensitiveCompare(technology, Constants::SAM_RG_SEQTECHNOLOGY_HELICOS)    ||
 395          caseInsensitiveCompare(technology, Constants::SAM_RG_SEQTECHNOLOGY_ILLUMINA)   ||
 396          caseInsensitiveCompare(technology, Constants::SAM_RG_SEQTECHNOLOGY_IONTORRENT) ||
 397          caseInsensitiveCompare(technology, Constants::SAM_RG_SEQTECHNOLOGY_LS454)      ||
 398          caseInsensitiveCompare(technology, Constants::SAM_RG_SEQTECHNOLOGY_PACBIO)     ||
 399          caseInsensitiveCompare(technology, Constants::SAM_RG_SEQTECHNOLOGY_SOLID)
 400        )
 401     {
 402         return true;
 403     }
 404
 405     // otherwise
 406     AddError("Invalid read group sequencing platform (PL): " + technology);
 407     return false;
 408 }
 409
 410 bool SamHeaderValidator::ValidateProgramChain(void) {
 411     bool isValid = true;
 412     isValid &= ContainsUniqueProgramIds();
 413     isValid &= ValidatePreviousProgramIds();
 414     return isValid;
 415 }
 416
 417 bool SamHeaderValidator::ContainsUniqueProgramIds(void) {
 418
 419     bool isValid = true;
 420     set<string> programIds;
 421     set<string>::iterator pgIdIter;
 422
 423     // iterate over program records
 424     const SamProgramChain& programs = m_header.Programs;
 425     SamProgramConstIterator pgIter = programs.ConstBegin();
 426     SamProgramConstIterator pgEnd  = programs.ConstEnd();
 427     for ( ; pgIter != pgEnd; ++pgIter ) {
 428         const SamProgram& pg = (*pgIter);
 429
 430         // lookup program ID
 431         const string& pgId = pg.ID;
 432         pgIdIter = programIds.find(pgId);
 433
 434         // error if found (duplicate entry)
 435         if ( pgIdIter != programIds.end() ) {
 436             AddError("Program ID (ID): " + pgId + " is not unique");
 437             isValid = false;
 438         }
 439
 440         // otherwise ok, store ID
 441         programIds.insert(pgId);
 442     }
 443
 444     // return validation state
 445     return isValid;
 446 }
 447
 448 bool SamHeaderValidator::ValidatePreviousProgramIds(void) {
 449
 450     bool isValid = true;
 451
 452     // iterate over program records
 453     const SamProgramChain& programs = m_header.Programs;
 454     SamProgramConstIterator pgIter = programs.ConstBegin();
 455     SamProgramConstIterator pgEnd  = programs.ConstEnd();
 456     for ( ; pgIter != pgEnd; ++pgIter ) {
 457         const SamProgram& pg = (*pgIter);
 458
 459         // ignore record for validation if PreviousProgramID is empty
 460         const string& ppId = pg.PreviousProgramID;
 461         if ( ppId.empty() )
 462             continue;
 463
 464         // see if program "chain" contains an entry for ppId
 465         if ( !programs.Contains(ppId) ) {
 466             AddError("PreviousProgramID (PP): " + ppId + " is not a known ID");
 467             isValid = false;
 468         }
 469     }
 470
 471     // return validation state
 472     return isValid;
 473 }
 474 void SamHeaderValidator::AddError(const string& message) {
 475     m_errorMessages.push_back(ERROR_PREFIX + message + NEWLINE);
 476 }
 477
 478 void SamHeaderValidator::AddWarning(const string& message) {
 479     m_warningMessages.push_back(WARN_PREFIX + message + NEWLINE);
 480 }
 481
 482 void SamHeaderValidator::PrintErrorMessages(void) {
 483
 484     // skip if no error messages
 485     if ( m_errorMessages.empty() ) return;
 486
 487     // print error header line
 488     cerr << "* SAM header has " << m_errorMessages.size() << " errors:" << endl;
 489
 490     // print each error message
 491     vector<string>::const_iterator errorIter = m_errorMessages.begin();
 492     vector<string>::const_iterator errorEnd  = m_errorMessages.end();
 493     for ( ; errorIter != errorEnd; ++errorIter )
 494         cerr << (*errorIter);
 495 }
 496
 497 void SamHeaderValidator::PrintWarningMessages(void) {
 498
 499     // skip if no warning messages
 500     if ( m_warningMessages.empty() ) return;
 501
 502     // print warning header line
 503     cerr << "* SAM header has " << m_warningMessages.size() << " warnings:" << endl;
 504
 505     // print each warning message
 506     vector<string>::const_iterator warnIter = m_warningMessages.begin();
 507     vector<string>::const_iterator warnEnd  = m_warningMessages.end();
 508     for ( ; warnIter != warnEnd; ++warnIter )
 509         cerr << (*warnIter);
 510 }