molecule_identifier.cpp
上传用户:yhdzpy8989
上传日期:2007-06-13
资源大小:13604k
文件大小:14k
源码类别:

生物技术

开发平台:

C/C++

  1. /*
  2.  * ===========================================================================
  3.  * PRODUCTION $Log: molecule_identifier.cpp,v $
  4.  * PRODUCTION Revision 1000.2  2004/06/01 18:28:46  gouriano
  5.  * PRODUCTION PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.13
  6.  * PRODUCTION
  7.  * ===========================================================================
  8.  */
  9. /*  $Id: molecule_identifier.cpp,v 1000.2 2004/06/01 18:28:46 gouriano Exp $
  10. * ===========================================================================
  11. *
  12. *                            PUBLIC DOMAIN NOTICE
  13. *               National Center for Biotechnology Information
  14. *
  15. *  This software/database is a "United States Government Work" under the
  16. *  terms of the United States Copyright Act.  It was written as part of
  17. *  the author's official duties as a United States Government employee and
  18. *  thus cannot be copyrighted.  This software/database is freely available
  19. *  to the public for use. The National Library of Medicine and the U.S.
  20. *  Government have not placed any restriction on its use or reproduction.
  21. *
  22. *  Although all reasonable efforts have been taken to ensure the accuracy
  23. *  and reliability of the software and data, the NLM and the U.S.
  24. *  Government do not and cannot warrant the performance or results that
  25. *  may be obtained by using this software or data. The NLM and the U.S.
  26. *  Government disclaim all warranties, express or implied, including
  27. *  warranties of performance, merchantability or fitness for any particular
  28. *  purpose.
  29. *
  30. *  Please cite the author in any work or product based on this material.
  31. *
  32. * ===========================================================================
  33. *
  34. * Authors:  Paul Thiessen
  35. *
  36. * File Description:
  37. *      Class to hold, and factory to generate, general
  38. *      (instance-independent) identifier for any molecule
  39. *
  40. * ===========================================================================
  41. */
  42. #ifdef _MSC_VER
  43. #pragma warning(disable:4018)   // disable signed/unsigned mismatch warning in MSVC
  44. #endif
  45. #include <ncbi_pch.hpp>
  46. #include <corelib/ncbistre.hpp>
  47. #include <corelib/ncbistl.hpp>
  48. #include <objects/seqloc/PDB_seq_id.hpp>
  49. #include <objects/general/Object_id.hpp>
  50. #include <objects/seqloc/Textseq_id.hpp>
  51. #include <objects/seqloc/PDB_mol_id.hpp>
  52. #include <memory>
  53. #include "molecule_identifier.hpp"
  54. #include "structure_set.hpp"
  55. #include "molecule.hpp"
  56. #include "sequence_set.hpp"
  57. #include "cn3d_tools.hpp"
  58. USING_NCBI_SCOPE;
  59. BEGIN_SCOPE(Cn3D)
  60. // there is one (global) list of molecule identifiers
  61. typedef list < MoleculeIdentifier > MoleculeIdentifierList;
  62. static MoleculeIdentifierList knownIdentifiers;
  63. const int MoleculeIdentifier::VALUE_NOT_SET = -1;
  64. const MoleculeIdentifier * MoleculeIdentifier::GetIdentifier(const Molecule *molecule,
  65.     string _pdbID, int _pdbChain, int _gi, string _accession)
  66. {
  67.     const StructureObject *object;
  68.     if (!molecule->GetParentOfType(&object)) return NULL;
  69.     // see if there's already an identifier that matches this molecule
  70.     MoleculeIdentifier *identifier = NULL;
  71.     MoleculeIdentifierList::iterator i, ie = knownIdentifiers.end();
  72.     for (i=knownIdentifiers.begin(); i!=ie; ++i) {
  73.         if ((object->mmdbID != StructureObject::NO_MMDB_ID && object->mmdbID == i->mmdbID &&
  74.              i->moleculeID != VALUE_NOT_SET && molecule->id == i->moleculeID) ||
  75.             (_pdbID.size() > 0 && _pdbID == i->pdbID &&
  76.              _pdbChain != VALUE_NOT_SET && _pdbChain == i->pdbChain) ||
  77.             (_gi != VALUE_NOT_SET && _gi == i->gi) ||
  78.             (_accession.size() > 0 && _accession == i->accession)) {
  79.             identifier = &(*i);
  80.             break;
  81.         }
  82.     }
  83.     // if no equivalent found, create a new one
  84.     if (!identifier) {
  85.         knownIdentifiers.resize(knownIdentifiers.size() + 1, MoleculeIdentifier());
  86.         identifier = &(knownIdentifiers.back());
  87.     }
  88.     // check consistency, and see if there's some more information we can fill out
  89.     if (object->mmdbID != StructureObject::NO_MMDB_ID) {
  90.         if (identifier->mmdbID != VALUE_NOT_SET &&
  91.             identifier->mmdbID != object->mmdbID && identifier->moleculeID != molecule->id) {
  92.             ERRORMSG("MoleculeIdentifier::GetIdentifier() - mmdb/molecule ID mismatch");
  93.         } else {
  94.             identifier->mmdbID = object->mmdbID;
  95.             identifier->moleculeID = molecule->id;
  96.         }
  97.     }
  98.     if (_pdbID.size() > 0) {
  99.         if (identifier->pdbID.size() > 0) {
  100.             if (identifier->pdbID != _pdbID || identifier->pdbChain != _pdbChain)
  101.                 ERRORMSG("MoleculeIdentifier::GetIdentifier() - pdbID/chain mismatch");
  102.         } else {
  103.             identifier->pdbID = _pdbID;
  104.             identifier->pdbChain = _pdbChain;
  105.         }
  106.     }
  107.     if (_gi != VALUE_NOT_SET) {
  108.         if (identifier->gi != VALUE_NOT_SET) {
  109.             if (identifier->gi != _gi)
  110.                 ERRORMSG("MoleculeIdentifier::GetIdentifier() - gi mismatch");
  111.         } else {
  112.             identifier->gi = _gi;
  113.         }
  114.     }
  115.     if (_accession.size() > 0) {
  116.         if (identifier->accession.size() > 0) {
  117.             if (identifier->accession != _accession)
  118.                 ERRORMSG("MoleculeIdentifier::GetIdentifier() - accession mismatch");
  119.         } else {
  120.             identifier->accession = _accession;
  121.         }
  122.     }
  123.     if (identifier->nResidues == 0)
  124.         identifier->nResidues = molecule->NResidues();
  125.     else if (identifier->nResidues != molecule->NResidues())
  126.         ERRORMSG("Length mismatch in molecule identifier " << identifier->ToString());
  127.     if (molecule->IsProtein() || molecule->IsNucleotide())
  128.         TRACEMSG("biopolymer molecule: identifier " << identifier << " (" << identifier->ToString() << ')');
  129.     return identifier;
  130. }
  131. const MoleculeIdentifier * MoleculeIdentifier::GetIdentifier(const Sequence *sequence,
  132.     string _pdbID, int _pdbChain, int _mmdbID, int _gi, string _accession)
  133. {
  134.     // see if there's already an identifier that matches this sequence
  135.     MoleculeIdentifier *identifier = NULL;
  136.     MoleculeIdentifierList::iterator i, ie = knownIdentifiers.end();
  137.     for (i=knownIdentifiers.begin(); i!=ie; ++i) {
  138.         if ((_pdbID.size() > 0 && _pdbID == i->pdbID &&
  139.              _pdbChain != VALUE_NOT_SET && _pdbChain == i->pdbChain) ||
  140.             (_gi != VALUE_NOT_SET && _gi == i->gi) ||
  141.             (_accession.size() > 0 && _accession == i->accession)) {
  142.             identifier = &(*i);
  143.             break;
  144.         }
  145.     }
  146.     // if no equivalent found, create a new one
  147.     if (!identifier) {
  148.         knownIdentifiers.resize(knownIdentifiers.size() + 1, MoleculeIdentifier());
  149.         identifier = &(knownIdentifiers.back());
  150.     }
  151.     // check consistency, and see if there's some more information we can fill out
  152.     if (_mmdbID != VALUE_NOT_SET) {
  153.         if (identifier->mmdbID != VALUE_NOT_SET) {
  154.             if (identifier->mmdbID != _mmdbID)
  155.                 ERRORMSG("MoleculeIdentifier::GetIdentifier() - mmdbID mismatch");
  156.         } else {
  157.             identifier->mmdbID = _mmdbID;
  158.         }
  159.     }
  160.     if (_pdbID.size() > 0) {
  161.         if (identifier->pdbID.size() > 0) {
  162.             if (identifier->pdbID != _pdbID || identifier->pdbChain != _pdbChain)
  163.                 ERRORMSG("MoleculeIdentifier::GetIdentifier() - pdbID/chain mismatch");
  164.         } else {
  165.             identifier->pdbID = _pdbID;
  166.             identifier->pdbChain = _pdbChain;
  167.         }
  168.     }
  169.     if (_gi != VALUE_NOT_SET) {
  170.         if (identifier->gi != VALUE_NOT_SET) {
  171.             if (identifier->gi != _gi) {
  172.                 if (identifier->accession.size() > 0 && identifier->accession == _accession) {
  173.                     // special case where accession is same, gi is different: two versions of same sequence
  174.                     ERRORMSG("The sequence being loaded (gi " << _gi << ") has the same accession code "
  175.                         "as a sequence already present (gi " << identifier->gi << "). Cn3D does not "
  176.                         "allow two versions of the same sequence to be present simultaneously. If you want "
  177.                         "to import gi " << _gi << ", you need to delete gi " << identifier->gi << " from "
  178.                         "the alignment and import windows, save & reload, then import again.");
  179.                     return NULL;
  180.                 } else {
  181.                     ERRORMSG("MoleculeIdentifier::GetIdentifier() - gi mismatch");
  182.                 }
  183.             }
  184.         } else {
  185.             identifier->gi = _gi;
  186.         }
  187.     }
  188.     if (_accession.size() > 0) {
  189.         if (identifier->accession.size() > 0) {
  190.             if (identifier->accession != _accession)
  191.                 ERRORMSG("MoleculeIdentifier::GetIdentifier() - accession mismatch");
  192.         } else {
  193.             identifier->accession = _accession;
  194.         }
  195.     }
  196.     if (identifier->nResidues == 0)
  197.         identifier->nResidues = sequence->Length();
  198.     else if (identifier->nResidues != sequence->Length())
  199.         ERRORMSG("Length mismatch in sequence identifier " << identifier->ToString());
  200. //    TESTMSG("sequence: identifier " << identifier << " (" << identifier->ToString() << ')');
  201.     return identifier;
  202. }
  203. string MoleculeIdentifier::ToString(void) const
  204. {
  205.     CNcbiOstrstream oss;
  206.     if (pdbID.size() > 0 && pdbChain != VALUE_NOT_SET) {
  207.         oss << pdbID;
  208.         if (pdbChain != ' ') {
  209.             oss <<  '_' << (char) pdbChain;
  210.         }
  211.     } else if (gi != VALUE_NOT_SET)
  212.         oss << "gi " << gi;
  213.     else if (accession.size() > 0)
  214.         oss << accession;
  215.     else if (mmdbID != VALUE_NOT_SET && moleculeID != VALUE_NOT_SET) {
  216.         oss << "mmdb " << mmdbID << " molecule " << moleculeID;
  217.     } else
  218.         oss << '?';
  219.     oss << '';
  220.     auto_ptr<char> chars(oss.str());    // frees memory upon function return
  221.     return string(oss.str());
  222. }
  223. const MoleculeIdentifier * MoleculeIdentifier::FindIdentifier(int mmdbID, int moleculeID)
  224. {
  225.     const MoleculeIdentifier *identifier = NULL;
  226.     MoleculeIdentifierList::const_iterator i, ie = knownIdentifiers.end();
  227.     for (i=knownIdentifiers.begin(); i!=ie; ++i) {
  228.         if (mmdbID == i->mmdbID && moleculeID == i->moleculeID) {
  229.             identifier = &(*i);
  230.             break;
  231.         }
  232.     }
  233.     return identifier;
  234. }
  235. void MoleculeIdentifier::ClearIdentifiers(void)
  236. {
  237.     knownIdentifiers.clear();
  238. }
  239. bool MoleculeIdentifier::MatchesSeqId(const ncbi::objects::CSeq_id& sid) const
  240. {
  241.     if (sid.IsGi())
  242.         return (sid.GetGi() == gi);
  243.     if (sid.IsPdb()) {
  244.         if (sid.GetPdb().GetMol().Get() == pdbID) {
  245.             if (sid.GetPdb().IsSetChain() && sid.GetPdb().GetChain() != pdbChain)
  246.                 return false;
  247.             else
  248.                 return true;
  249.         } else
  250.             return false;
  251.     }
  252.     if (sid.IsLocal() && sid.GetLocal().IsStr())
  253.         return (sid.GetLocal().GetStr() == accession || (accession.size() == 0 &&
  254.                     // special case where local accession is actually a PDB identifier + extra stuff
  255.                     sid.GetLocal().GetStr().size() >= 7 && sid.GetLocal().GetStr()[4] == ' ' &&
  256.                     sid.GetLocal().GetStr()[6] == ' ' && isalpha(sid.GetLocal().GetStr()[5]) &&
  257.                     sid.GetLocal().GetStr().substr(0, 4) == pdbID && sid.GetLocal().GetStr()[5] == pdbChain));
  258.     if (sid.IsGenbank() && sid.GetGenbank().IsSetAccession())
  259.         return (sid.GetGenbank().GetAccession() == accession);
  260.     if (sid.IsSwissprot() && sid.GetSwissprot().IsSetAccession())
  261.         return (sid.GetSwissprot().GetAccession() == accession);
  262.     if (sid.IsOther() && sid.GetOther().IsSetAccession())
  263.         return (sid.GetOther().GetAccession() == accession);
  264.     if (sid.IsEmbl() && sid.GetEmbl().IsSetAccession())
  265.         return (sid.GetEmbl().GetAccession() == accession);
  266.     if (sid.IsDdbj() && sid.GetDdbj().IsSetAccession())
  267.         return (sid.GetDdbj().GetAccession() == accession);
  268.     ERRORMSG("MoleculeIdentifier::MatchesSeqId() - can't match this type of Seq-id");
  269.     return false;
  270. }
  271. bool MoleculeIdentifier::CompareIdentifiers(const MoleculeIdentifier *a, const MoleculeIdentifier *b)
  272. {
  273.     // identifier sort - float sequences with PDB id's to the top, then gi's, then accessions
  274.     if (a->pdbID.size() > 0) {
  275.         if (b->pdbID.size() > 0) {
  276.             if (a->pdbID < b->pdbID)
  277.                 return true;
  278.             else if (a->pdbID > b->pdbID)
  279.                 return false;
  280.             else
  281.                 return (a->pdbChain < b->pdbChain);
  282.         } else
  283.             return true;
  284.     }
  285.     else if (a->gi != VALUE_NOT_SET) {
  286.         if (b->pdbID.size() > 0)
  287.             return false;
  288.         else if (b->gi != VALUE_NOT_SET)
  289.             return (a->gi < b->gi);
  290.         else
  291.             return true;
  292.     }
  293.     else if (a->accession.size() > 0) {
  294.         if (b->pdbID.size() > 0 || b->gi != VALUE_NOT_SET)
  295.             return false;
  296.         else if (b->accession.size() > 0)
  297.             return (a->accession < b->accession);
  298.     }
  299.     ERRORMSG("MoleculeIdentifier::CompareIdentifiers() - confused by identifier type");
  300.     return false;
  301. }
  302. END_SCOPE(Cn3D)
  303. /*
  304. * ---------------------------------------------------------------------------
  305. * $Log: molecule_identifier.cpp,v $
  306. * Revision 1000.2  2004/06/01 18:28:46  gouriano
  307. * PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.13
  308. *
  309. * Revision 1.13  2004/05/21 21:41:39  gorelenk
  310. * Added PCH ncbi_pch.hpp
  311. *
  312. * Revision 1.12  2004/03/15 18:25:36  thiessen
  313. * prefer prefix vs. postfix ++/-- operators
  314. *
  315. * Revision 1.11  2004/02/19 17:04:59  thiessen
  316. * remove cn3d/ from include paths; add pragma to disable annoying msvc warning
  317. *
  318. * Revision 1.10  2004/01/05 17:09:15  thiessen
  319. * abort import and warn if same accession different gi
  320. *
  321. * Revision 1.9  2003/09/03 18:14:01  thiessen
  322. * hopefully fix Seq-id issues
  323. *
  324. * Revision 1.8  2003/02/03 19:20:04  thiessen
  325. * format changes: move CVS Log to bottom of file, remove std:: from .cpp files, and use new diagnostic macros
  326. *
  327. * Revision 1.7  2002/02/22 14:24:00  thiessen
  328. * sort sequences in reject dialog ; general identifier comparison
  329. *
  330. * Revision 1.6  2002/01/24 20:08:16  thiessen
  331. * fix local id problem
  332. *
  333. * Revision 1.5  2001/12/12 14:04:14  thiessen
  334. * add missing object headers after object loader change
  335. *
  336. * Revision 1.4  2001/11/27 16:26:08  thiessen
  337. * major update to data management system
  338. *
  339. * Revision 1.3  2001/08/08 02:25:27  thiessen
  340. * add <memory>
  341. *
  342. * Revision 1.2  2001/07/04 19:39:17  thiessen
  343. * finish user annotation system
  344. *
  345. * Revision 1.1  2001/06/21 02:02:34  thiessen
  346. * major update to molecule identification and highlighting ; add toggle highlight (via alt)
  347. *
  348. */