开通博客赚积分
发布资源赚积分
分类

源码开发语言/平台
当前位置：首页 > 源码/资料 > 数值算法/人工智能 > 生物技术 > 查看源码
seq_vector.cpp
资源名称：ncbi_cxx_win--Jul_08_2004.zip　[点击查看]
上传用户：yhdzpy8989
上传日期：2007-06-13
资源大小：13604k
文件大小：17k
源码类别：

生物技术

开发平台：

C/C++

seq_vector.cpp：源码内容
  1. /*
  2.  * ===========================================================================
  3.  * PRODUCTION $Log: seq_vector.cpp,v $
  4.  * PRODUCTION Revision 1000.3  2004/06/01 19:24:24  gouriano
  5.  * PRODUCTION PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.68
  6.  * PRODUCTION
  7.  * ===========================================================================
  8.  */
  10. /*  $Id: seq_vector.cpp,v 1000.3 2004/06/01 19:24:24 gouriano Exp $
  11. * ===========================================================================
  12. *
  13. *                            PUBLIC DOMAIN NOTICE
  14. *               National Center for Biotechnology Information
  15. *
  16. *  This software/database is a "United States Government Work" under the
  17. *  terms of the United States Copyright Act.  It was written as part of
  18. *  the author's official duties as a United States Government employee and
  19. *  thus cannot be copyrighted.  This software/database is freely available
  20. *  to the public for use. The National Library of Medicine and the U.S.
  21. *  Government have not placed any restriction on its use or reproduction.
  22. *
  23. *  Although all reasonable efforts have been taken to ensure the accuracy
  24. *  and reliability of the software and data, the NLM and the U.S.
  25. *  Government do not and cannot warrant the performance or results that
  26. *  may be obtained by using this software or data. The NLM and the U.S.
  27. *  Government disclaim all warranties, express or implied, including
  28. *  warranties of performance, merchantability or fitness for any particular
  29. *  purpose.
  30. *
  31. *  Please cite the author in any work or product based on this material.
  32. *
  33. * ===========================================================================
  34. *
  35. * Author: Aleksey Grichenko, Eugene Vasilchenko
  36. *
  37. * File Description:
  38. *   Sequence data container for object manager
  39. *
  40. */
  43. #include <ncbi_pch.hpp>
  44. #include <objmgr/seq_vector.hpp>
  45. #include <objmgr/seq_vector_ci.hpp>
  46. #include <corelib/ncbimtx.hpp>
  47. #include <objmgr/impl/data_source.hpp>
  48. #include <objects/seq/seqport_util.hpp>
  49. #include <objects/seqloc/Seq_loc.hpp>
  50. #include <objmgr/seq_map.hpp>
  51. #include <objmgr/objmgr_exception.hpp>
  52. #include <algorithm>
  53. #include <map>
  54. #include <vector>
  56. BEGIN_NCBI_SCOPE
  57. BEGIN_SCOPE(objects)
  60. ////////////////////////////////////////////////////////////////////
  61. //
  62. //  CSeqVector::
  63. //
  66. CSeqVector::CSeqVector(void)
  67.     : m_Size(0)
  68. {
  69.     m_Iterator.x_SetVector(*this);
  70. }
  73. CSeqVector::CSeqVector(const CSeqVector& vec)
  74.     : m_SeqMap(vec.m_SeqMap),
  75.       m_Scope(vec.m_Scope),
  76.       m_Size(vec.m_Size),
  77.       m_Mol(vec.m_Mol),
  78.       m_Strand(vec.m_Strand),
  79.       m_Coding(vec.m_Coding)
  80. {
  81.     m_Iterator.x_SetVector(*this);
  82. }
  85. CSeqVector::CSeqVector(CConstRef<CSeqMap> seqMap, CScope& scope,
  86.                        EVectorCoding coding, ENa_strand strand)
  87.     : m_SeqMap(seqMap),
  88.       m_Scope(&scope),
  89.       m_Size(seqMap->GetLength(&scope)),
  90.       m_Mol(seqMap->GetMol()),
  91.       m_Strand(strand),
  92.       m_Coding(CSeq_data::e_not_set)
  93. {
  94.     m_Iterator.x_SetVector(*this);
  95.     SetCoding(coding);
  96. }
  98. CSeqVector::CSeqVector(const CSeqMap& seqMap, CScope& scope,
  99.                        EVectorCoding coding, ENa_strand strand)
  100.     : m_SeqMap(&seqMap),
  101.       m_Scope(&scope),
  102.       m_Size(seqMap.GetLength(&scope)),
  103.       m_Mol(seqMap.GetMol()),
  104.       m_Strand(strand),
  105.       m_Coding(CSeq_data::e_not_set)
  106. {
  107.     m_Iterator.x_SetVector(*this);
  108.     SetCoding(coding);
  109. }
  112. CSeqVector::CSeqVector(const CSeq_loc& loc, CScope& scope,
  113.                        EVectorCoding coding, ENa_strand strand)
  114.     : m_SeqMap(CSeqMap::CreateSeqMapForSeq_loc(loc, &scope)),
  115.       m_Scope(&scope),
  116.       m_Size(m_SeqMap->GetLength(&scope)),
  117.       m_Mol(m_SeqMap->GetMol()),
  118.       m_Strand(strand),
  119.       m_Coding(CSeq_data::e_not_set)
  120. {
  121.     m_Iterator.x_SetVector(*this);
  122.     SetCoding(coding);
  123. }
  126. CSeqVector::~CSeqVector(void)
  127. {
  128. }
  131. CSeqVector& CSeqVector::operator= (const CSeqVector& vec)
  132. {
  133.     if ( &vec != this ) {
  134.         m_SeqMap = vec.m_SeqMap;
  135.         m_Scope  = vec.m_Scope;
  136.         m_Size   = vec.m_Size;
  137.         m_Mol    = vec.m_Mol;
  138.         m_Strand = vec.m_Strand;
  139.         m_Coding = vec.m_Coding;
  140.         m_Iterator.x_SetVector(*this);
  141.     }
  142.     return *this;
  143. }
  146. bool CSeqVector::CanGetRange(TSeqPos from, TSeqPos to) const
  147. {
  148.     return m_SeqMap->CanResolveRange(m_Scope.GetScopeOrNull(),
  149.                                      from, to, m_Strand);
  150. }
  153. CSeqVector::TResidue CSeqVector::x_GetGapChar(TCoding coding)
  154. {
  155.     switch (coding) {
  156.     case CSeq_data::e_Iupacna: // DNA - N
  157.         return 'N';
  158.     
  159.     case CSeq_data::e_Ncbi8na: // DNA - bit representation
  160.     case CSeq_data::e_Ncbi4na:
  161.         return 0x0f;  // all bits set == any base
  163.     case CSeq_data::e_Ncbieaa: // Proteins - X
  164.     case CSeq_data::e_Iupacaa:
  165.         return 'X';
  166.     
  167.     case CSeq_data::e_Ncbi8aa: // Protein - numeric representation
  168.     case CSeq_data::e_Ncbistdaa:
  169.         return 21;
  171.     case CSeq_data::e_not_set:
  172.         return 0;     // It's not good to throw an exception here
  174.     case CSeq_data::e_Ncbi2na: // Codings without gap symbols
  175.     case CSeq_data::e_Ncbipaa: //### Not sure about this
  176.     case CSeq_data::e_Ncbipna: //### Not sure about this
  177.     default:
  178.         NCBI_THROW(CSeqVectorException, eCodingError,
  179.                    "Can not indicate gap using the selected coding");
  180.     }
  181. }
  184. DEFINE_STATIC_FAST_MUTEX(s_ConvertTableMutex2);
  186. const char* CSeqVector::sx_GetConvertTable(TCoding src,
  187.                                            TCoding dst,
  188.                                            bool reverse)
  189. {
  190.     CFastMutexGuard guard(s_ConvertTableMutex2);
  191.     typedef pair<pair<TCoding, TCoding>, bool> TKey;
  192.     typedef map<TKey, vector<char> > TTables;
  193.     static TTables tables;
  195.     TKey key(pair<TCoding, TCoding>(src, dst), reverse);
  196.     TTables::iterator it = tables.find(key);
  197.     if ( it != tables.end() ) {
  198.         // already created
  199.         return it->second.empty()? 0: &it->second[0];
  200.     }
  201.     it = tables.insert(TTables::value_type(key, vector<char>())).first;
  202.     if ( !CSeqportUtil::IsCodeAvailable(src) ||
  203.          !CSeqportUtil::IsCodeAvailable(dst) ) {
  204.         // invalid types
  205.         return 0;
  206.     }
  208.     const size_t COUNT = kMax_UChar+1;
  209.     const unsigned kInvalidCode = kMax_UChar;
  211.     pair<unsigned, unsigned> srcIndex = CSeqportUtil::GetCodeIndexFromTo(src);
  212.     if ( srcIndex.second >= COUNT ) {
  213.         // too large range
  214.         return 0;
  215.     }
  217.     if ( reverse ) {
  218.         // check if src needs complement conversion
  219.         try {
  220.             CSeqportUtil::GetIndexComplement(src, srcIndex.first);
  221.         }
  222.         catch ( exception& /*noComplement*/ ) {
  223.             reverse = false;
  224.         }
  225.     }
  227.     if ( dst != src ) {
  228.         pair<unsigned, unsigned> dstIndex =
  229.             CSeqportUtil::GetCodeIndexFromTo(dst);
  230.         if ( dstIndex.second >= COUNT ) {
  231.             // too large range
  232.             return 0;
  233.         }
  235.         try {
  236.             // check for types compatibility
  237.             CSeqportUtil::GetMapToIndex(src, dst, srcIndex.first);
  238.         }
  239.         catch ( exception& /*badType*/ ) {
  240.             // incompatible types
  241.             return 0;
  242.         }
  243.     }
  244.     else if ( !reverse ) {
  245.         // no need to convert at all
  246.         return 0;
  247.     }
  249.     it->second.resize(COUNT, char(kInvalidCode));
  250.     for ( unsigned i = srcIndex.first; i <= srcIndex.second; ++i ) {
  251.         try {
  252.             unsigned code = i;
  253.             if ( reverse ) {
  254.                 code = CSeqportUtil::GetIndexComplement(src, code);
  255.             }
  256.             if ( dst != src ) {
  257.                 code = CSeqportUtil::GetMapToIndex(src, dst, code);
  258.             }
  259.             code = min(kInvalidCode, code);
  260.             it->second[i] = char(code);
  261.         }
  262.         catch ( exception& /*noConversion or noComplement*/ ) {
  263.         }
  264.     }
  265.     return &it->second[0];
  266. }
  269. void CSeqVector::SetCoding(TCoding coding)
  270. {
  271.     if (m_Coding != coding) {
  272.         m_Coding = coding;
  273.     }
  274.     m_Iterator.SetCoding(coding);
  275. }
  278. void CSeqVector::SetIupacCoding(void)
  279. {
  280.     SetCoding(IsProtein()? CSeq_data::e_Iupacaa: CSeq_data::e_Iupacna);
  281. }
  284. void CSeqVector::SetNcbiCoding(void)
  285. {
  286.     SetCoding(IsProtein()? CSeq_data::e_Ncbistdaa: CSeq_data::e_Ncbi4na);
  287. }
  290. void CSeqVector::SetCoding(EVectorCoding coding)
  291. {
  292.     switch ( coding ) {
  293.     case CBioseq_Handle::eCoding_Iupac:
  294.         SetIupacCoding();
  295.         break;
  296.     case CBioseq_Handle::eCoding_Ncbi:
  297.         SetNcbiCoding();
  298.         break;
  299.     default:
  300.         SetCoding(CSeq_data::e_not_set);
  301.         break;
  302.     }
  303. }
  306. END_SCOPE(objects)
  307. END_NCBI_SCOPE
  309. /*
  310. * ---------------------------------------------------------------------------
  311. * $Log: seq_vector.cpp,v $
  312. * Revision 1000.3  2004/06/01 19:24:24  gouriano
  313. * PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.68
  314. *
  315. * Revision 1.68  2004/05/21 21:42:13  gorelenk
  316. * Added PCH ncbi_pch.hpp
  317. *
  318. * Revision 1.67  2004/04/12 16:49:16  vasilche
  319. * Allow null scope in CSeqMap_CI and CSeqVector.
  320. *
  321. * Revision 1.66  2004/02/25 18:58:47  shomrat
  322. * Added a new construtor based on Seq-loc in scope
  323. *
  324. * Revision 1.65  2003/11/19 22:18:04  grichenk
  325. * All exceptions are now CException-derived. Catch "exception" rather
  326. * than "runtime_error".
  327. *
  328. * Revision 1.64  2003/10/08 14:16:55  vasilche
  329. * Removed circular reference CSeqVector <-> CSeqVector_CI.
  330. *
  331. * Revision 1.63  2003/09/30 16:22:04  vasilche
  332. * Updated internal object manager classes to be able to load ID2 data.
  333. * SNP blobs are loaded as ID2 split blobs - readers convert them automatically.
  334. * Scope caches results of requests for data to data loaders.
  335. * Optimized CSeq_id_Handle for gis.
  336. * Optimized bioseq lookup in scope.
  337. * Reduced object allocations in annotation iterators.
  338. * CScope is allowed to be destroyed before other objects using this scope are
  339. * deleted (feature iterators, bioseq handles etc).
  340. * Optimized lookup for matching Seq-ids in CSeq_id_Mapper.
  341. * Added 'adaptive' option to objmgr_demo application.
  342. *
  343. * Revision 1.62  2003/09/05 17:29:40  grichenk
  344. * Structurized Object Manager exceptions
  345. *
  346. * Revision 1.61  2003/08/29 13:34:47  vasilche
  347. * Rewrote CSeqVector/CSeqVector_CI code to allow better inlining.
  348. * CSeqVector::operator[] made significantly faster.
  349. * Added possibility to have platform dependent byte unpacking functions.
  350. *
  351. * Revision 1.60  2003/08/21 18:43:29  vasilche
  352. * Added CSeqVector::IsProtein() and CSeqVector::IsNucleotide() methods.
  353. *
  354. * Revision 1.59  2003/08/21 17:04:10  vasilche
  355. * Fixed bug in making conversion tables.
  356. *
  357. * Revision 1.58  2003/08/21 13:32:04  vasilche
  358. * Optimized CSeqVector iteration.
  359. * Set some CSeqVector values (mol type, coding) in constructor instead of detecting them while iteration.
  360. * Remove unsafe bit manipulations with coding.
  361. *
  362. * Revision 1.57  2003/06/24 19:46:43  grichenk
  363. * Changed cache from vector<char> to char*. Made
  364. * CSeqVector::operator[] inline.
  365. *
  366. * Revision 1.56  2003/06/17 20:35:39  grichenk
  367. * CSeqVector_CI-related improvements
  368. *
  369. * Revision 1.55  2003/06/13 17:22:28  grichenk
  370. * Check if seq-map is not null before using it
  371. *
  372. * Revision 1.54  2003/06/12 18:39:21  vasilche
  373. * Added default constructor of CSeqVector.
  374. * Cleared cache iterator on CSeqVector assignment.
  375. *
  376. * Revision 1.53  2003/06/11 19:32:55  grichenk
  377. * Added molecule type caching to CSeqMap, simplified
  378. * coding and sequence type calculations in CSeqVector.
  379. *
  380. * Revision 1.52  2003/06/04 13:48:56  grichenk
  381. * Improved double-caching, fixed problem with strands.
  382. *
  383. * Revision 1.51  2003/06/02 16:06:38  dicuccio
  384. * Rearranged src/objects/ subtree.  This includes the following shifts:
  385. *     - src/objects/asn2asn --> arc/app/asn2asn
  386. *     - src/objects/testmedline --> src/objects/ncbimime/test
  387. *     - src/objects/objmgr --> src/objmgr
  388. *     - src/objects/util --> src/objmgr/util
  389. *     - src/objects/alnmgr --> src/objtools/alnmgr
  390. *     - src/objects/flat --> src/objtools/flat
  391. *     - src/objects/validator --> src/objtools/validator
  392. *     - src/objects/cddalignview --> src/objtools/cddalignview
  393. * In addition, libseq now includes six of the objects/seq... libs, and libmmdb
  394. * replaces the three libmmdb? libs.
  395. *
  396. * Revision 1.50  2003/05/27 19:44:06  grichenk
  397. * Added CSeqVector_CI class
  398. *
  399. * Revision 1.49  2003/05/20 15:44:38  vasilche
  400. * Fixed interaction of CDataSource and CDataLoader in multithreaded app.
  401. * Fixed some warnings on WorkShop.
  402. * Added workaround for memory leak on WorkShop.
  403. *
  404. * Revision 1.48  2003/05/05 21:00:26  vasilche
  405. * Fix assignment of empty CSeqVector.
  406. *
  407. * Revision 1.47  2003/04/29 19:51:13  vasilche
  408. * Fixed interaction of Data Loader garbage collector and TSE locking mechanism.
  409. * Made some typedefs more consistent.
  410. *
  411. * Revision 1.46  2003/04/24 16:12:38  vasilche
  412. * Object manager internal structures are splitted more straightforward.
  413. * Removed excessive header dependencies.
  414. *
  415. * Revision 1.45  2003/02/07 16:28:05  vasilche
  416. * Fixed delayed seq vector coding setting.
  417. *
  418. * Revision 1.44  2003/02/06 19:05:28  vasilche
  419. * Fixed old cache data copying.
  420. * Delayed sequence type (protein/dna) resolution.
  421. *
  422. * Revision 1.43  2003/02/05 17:59:17  dicuccio
  423. * Moved formerly private headers into include/objects/objmgr/impl
  424. *
  425. * Revision 1.42  2003/02/05 15:05:28  vasilche
  426. * Added string::reserve()
  427. *
  428. * Revision 1.41  2003/01/28 17:17:06  vasilche
  429. * Fixed bug processing minus strands.
  430. * Used CSeqMap::ResolvedRangeIterator with strand coordinate translation.
  431. *
  432. * Revision 1.40  2003/01/23 19:33:57  vasilche
  433. * Commented out obsolete methods.
  434. * Use strand argument of CSeqVector instead of creation reversed seqmap.
  435. * Fixed ordering operators of CBioseqHandle to be consistent.
  436. *
  437. * Revision 1.39  2003/01/23 18:26:02  ucko
  438. * Explicitly #include <algorithm>
  439. *
  440. * Revision 1.38  2003/01/22 20:11:54  vasilche
  441. * Merged functionality of CSeqMapResolved_CI to CSeqMap_CI.
  442. * CSeqMap_CI now supports resolution and iteration over sequence range.
  443. * Added several caches to CScope.
  444. * Optimized CSeqVector().
  445. * Added serveral variants of CBioseqHandle::GetSeqVector().
  446. * Tried to optimize annotations iterator (not much success).
  447. * Rewritten CHandleRange and CHandleRangeMap classes to avoid sorting of list.
  448. *
  449. * Revision 1.37  2003/01/03 19:45:12  dicuccio
  450. * Replaced kPosUnknwon with kInvalidSeqPos (non-static variable; worka-round for
  451. * MSVC)
  452. *
  453. * Revision 1.36  2002/12/26 16:39:24  vasilche
  454. * Object manager class CSeqMap rewritten.
  455. *
  456. * Revision 1.35  2002/12/06 15:36:00  grichenk
  457. * Added overlap type for annot-iterators
  458. *
  459. * Revision 1.34  2002/10/03 13:45:39  grichenk
  460. * CSeqVector::size() made const
  461. *
  462. * Revision 1.33  2002/09/26 20:15:11  grichenk
  463. * Fixed cache lengths checks
  464. *
  465. * Revision 1.32  2002/09/16 20:13:01  grichenk
  466. * Fixed Iupac coding setter
  467. *
  468. * Revision 1.31  2002/09/16 13:52:49  dicuccio
  469. * Fixed bug in calculating total range of cached interval to retrieve -
  470. * must clamp the cached range to the desired range.
  471. *
  472. * Revision 1.30  2002/09/12 19:59:25  grichenk
  473. * Fixed bugs in calculating cached intervals
  474. *
  475. * Revision 1.29  2002/09/10 19:55:52  grichenk
  476. * Fixed reverse-complement position
  477. *
  478. * Revision 1.28  2002/09/09 21:38:37  grichenk
  479. * Fixed problem with GetSeqData() for minus strand
  480. *
  481. * Revision 1.27  2002/09/03 21:27:01  grichenk
  482. * Replaced bool arguments in CSeqVector constructor and getters
  483. * with enums.
  484. *
  485. * Revision 1.26  2002/07/08 20:51:02  grichenk
  486. * Moved log to the end of file
  487. * Replaced static mutex (in CScope, CDataSource) with the mutex
  488. * pool. Redesigned CDataSource data locking.
  489. *
  490. * Revision 1.25  2002/05/31 20:58:19  grichenk
  491. * Fixed GetSeqData() bug
  492. *
  493. * Revision 1.24  2002/05/31 17:53:00  grichenk
  494. * Optimized for better performance (CTSE_Info uses atomic counter,
  495. * delayed annotations indexing, no location convertions in
  496. * CAnnot_Types_CI if no references resolution is required etc.)
  497. *
  498. * Revision 1.23  2002/05/17 17:14:53  grichenk
  499. * +GetSeqData() for getting a range of characters from a seq-vector
  500. *
  501. * Revision 1.22  2002/05/09 14:16:31  grichenk
  502. * sm_SizeUnknown -> kPosUnknown, minor fixes for unsigned positions
  503. *
  504. * Revision 1.21  2002/05/06 03:28:47  vakatov
  505. * OM/OM1 renaming
  506. *
  507. * Revision 1.20  2002/05/03 21:28:10  ucko
  508. * Introduce T(Signed)SeqPos.
  509. *
  510. * Revision 1.19  2002/05/03 18:36:19  grichenk
  511. * Fixed members initialization
  512. *
  513. * Revision 1.18  2002/05/02 20:42:38  grichenk
  514. * throw -> THROW1_TRACE
  515. *
  516. * Revision 1.17  2002/04/29 16:23:28  grichenk
  517. * GetSequenceView() reimplemented in CSeqVector.
  518. * CSeqVector optimized for better performance.
  519. *
  520. * Revision 1.16  2002/04/26 14:37:21  grichenk
  521. * Limited CSeqVector cache size
  522. *
  523. * Revision 1.15  2002/04/25 18:14:47  grichenk
  524. * e_not_set coding gap symbol set to 0
  525. *
  526. * Revision 1.14  2002/04/25 16:37:21  grichenk
  527. * Fixed gap coding, added GetGapChar() function
  528. *
  529. * Revision 1.13  2002/04/23 19:01:07  grichenk
  530. * Added optional flag to GetSeqVector() and GetSequenceView()
  531. * for switching to IUPAC encoding.
  532. *
  533. * Revision 1.12  2002/04/22 20:03:08  grichenk
  534. * Updated comments
  535. *
  536. * Revision 1.11  2002/04/17 21:07:59  grichenk
  537. * String pre-allocation added
  538. *
  539. * Revision 1.10  2002/04/03 18:06:48  grichenk
  540. * Fixed segmented sequence bugs (invalid positioning of literals
  541. * and gaps). Improved CSeqVector performance.
  542. *
  543. * Revision 1.8  2002/03/28 18:34:58  grichenk
  544. * Fixed convertions bug
  545. *
  546. * Revision 1.7  2002/03/08 21:24:35  gouriano
  547. * fixed errors with unresolvable references
  548. *
  549. * Revision 1.6  2002/02/21 19:27:06  grichenk
  550. * Rearranged includes. Added scope history. Added searching for the
  551. * best seq-id match in data sources and scopes. Updated tests.
  552. *
  553. * Revision 1.5  2002/02/15 20:35:38  gouriano
  554. * changed implementation of HandleRangeMap
  555. *
  556. * Revision 1.4  2002/01/30 22:09:28  gouriano
  557. * changed CSeqMap interface
  558. *
  559. * Revision 1.3  2002/01/23 21:59:32  grichenk
  560. * Redesigned seq-id handles and mapper
  561. *
  562. * Revision 1.2  2002/01/16 16:25:58  gouriano
  563. * restructured objmgr
  564. *
  565. * Revision 1.1  2002/01/11 19:06:24  gouriano
  566. * restructured objmgr
  567. *
  568. *
  569. * ===========================================================================
  570. */