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

生物技术

开发平台:

C/C++

  1. /*
  2.  * ===========================================================================
  3.  * PRODUCTION $Log: sage_dload.cpp,v $
  4.  * PRODUCTION Revision 1000.4  2004/06/01 19:42:46  gouriano
  5.  * PRODUCTION PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.5
  6.  * PRODUCTION
  7.  * ===========================================================================
  8.  */
  9. /*  $Id: sage_dload.cpp,v 1000.4 2004/06/01 19:42: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:  Mike DiCuccio
  35.  *
  36.  * File Description:
  37.  *
  38.  */
  39. #include <ncbi_pch.hpp>
  40. #include <objtools/data_loaders/table/sage_dload.hpp>
  41. #include <sqlite/sqlite.hpp>
  42. #include <objects/general/User_object.hpp>
  43. #include <objects/general/User_field.hpp>
  44. #include <objtools/manip/sage_manip.hpp>
  45. #include <objects/seq/Seq_annot.hpp>
  46. #include <objects/seqfeat/Seq_feat.hpp>
  47. #include <objects/seqloc/Seq_loc.hpp>
  48. #include <objects/seqloc/Seq_interval.hpp>
  49. #include <objects/seqset/Seq_entry.hpp>
  50. #include <objects/seqset/Bioseq_set.hpp>
  51. #include <objmgr/scope.hpp>
  52. #include <objmgr/impl/data_source.hpp>
  53. #include <objmgr/impl/synonyms.hpp>
  54. #include <objmgr/impl/handle_range_map.hpp>
  55. BEGIN_NCBI_SCOPE
  56. USING_SCOPE(objects);
  57. bool CSageDataLoader::SIdHandleByContent::operator()(const CSeq_id_Handle& h1,
  58.                                                      const CSeq_id_Handle& h2) const
  59. {
  60.     CConstRef<CSeq_id> id1 = h1.GetSeqId();
  61.     CConstRef<CSeq_id> id2 = h2.GetSeqId();
  62.     return (*id1 < *id2);
  63. }
  64. CSageDataLoader::CSageDataLoader(const string& input_file,
  65.                                  const string& temp_file,
  66.                                  bool delete_file)
  67.     : CDataLoader(input_file)
  68. {
  69.     //
  70.     // create our SQLite DB
  71.     //
  72.     m_Table.Reset(new CSQLiteTable(input_file, temp_file, delete_file));
  73.     //
  74.     // now, store some precalculated info about our table
  75.     //
  76.     {{
  77.         // extract the column names
  78.         list<string> cols;
  79.         m_Table->GetColumnTitles(cols);
  80.         m_Cols.reserve(cols.size());
  81.         std::copy(cols.begin(), cols.end(), back_inserter(m_Cols));
  82.     }}
  83.     // determine our column mapping
  84.     int i = 0;
  85.     m_ColAssign.resize(m_Cols.size(), eUnknown);
  86.     fill(m_ColIdx, m_ColIdx + eMaxKnownCol, -1);
  87.     ITERATE(vector<string>, iter, m_Cols) {
  88.         string str(*iter);
  89.         NStr::ToLower(str);
  90.         if (str == "contig"  ||
  91.             str == "contig_accession"  ||
  92.             str == "accession") {
  93.             m_ColAssign[i] = eAccession;
  94.             m_ColIdx[eAccession] = i;
  95.         } else if (str == "tag_position_on_contig"  ||
  96.                    str == "pos"  ||
  97.                    str == "contig_pos") {
  98.             m_ColAssign[i] = eFrom;
  99.             m_ColIdx[eFrom] = i;
  100.         } else if (str == "tag"  ||
  101.                    str == "catgtag") {
  102.             m_ColAssign[i] = eTag;
  103.             m_ColIdx[eTag] = i;
  104.         } else if (str == "method") {
  105.             m_ColAssign[i] = eMethod;
  106.             m_ColIdx[eMethod] = i;
  107.         } else if (str == "tag_orientation"  ||
  108.                    str == "strand"  ||
  109.                    str == "orientation") {
  110.             m_ColAssign[i] = eStrand;
  111.             m_ColIdx[eStrand] = i;
  112.         }
  113.         ++i;
  114.     }
  115.     string acc_col = "col" + NStr::IntToString(m_ColIdx[eAccession]);
  116.     CSQLite& sqlite = m_Table->SetDB();
  117.     // create an index on accession
  118.     try {
  119.         sqlite.Execute("create index IDX_accession "
  120.                        "on TableData (" + acc_col + ")");
  121.     }
  122.     catch (...) {
  123.         // index already exists - ignored
  124.     }
  125.     // extract a list of the accessions we have
  126.     CRef<CSQLiteQuery> q
  127.         (sqlite.Compile("select distinct " + acc_col +
  128.                         " from TableData order by " + acc_col));
  129.     int count;
  130.     const char** data = NULL;
  131.     const char** cols = NULL;
  132.     while (q->NextRow(count, data, cols)) {
  133.         CRef<CSeq_id> id(new CSeq_id(data[0]));
  134.         if (id->Which() == CSeq_id::e_not_set) {
  135.             LOG_POST(Error << "failed to index id = " << data[0]);
  136.             continue;
  137.         }
  138.         CSeq_id_Handle handle = CSeq_id_Handle::GetHandle(*id);
  139.         m_Ids.insert(TIdMap::value_type(handle, data[0]));
  140.         _TRACE("  id = " << data[0]);
  141.     }
  142.     LOG_POST(Info << "CSageDataLoader: " << m_Ids.size() << " distinct ids");
  143. }
  144. // Request from a datasource using handles and ranges instead of seq-loc
  145. // The TSEs loaded in this call will be added to the tse_set.
  146. void CSageDataLoader::GetRecords(const CSeq_id_Handle& idh,
  147.                                  EChoice choice)
  148. {
  149.     string acc_col = "col" + NStr::IntToString(m_ColIdx[eAccession]);
  150.     //
  151.     // find out if we've already loaded annotations for this seq-id
  152.     //
  153.     TEntries::iterator iter = m_Entries.find(idh);
  154.     if (iter != m_Entries.end()) {
  155.         return;
  156.     }
  157.     //
  158.     // now, find out if this ID is in our list of ids
  159.     //
  160.     pair<TIdMap::iterator, TIdMap::iterator> id_iter = m_Ids.equal_range(idh);
  161.     if (id_iter.first == id_iter.second) {
  162.         return;
  163.     }
  164.     //
  165.     // okay, this ID is correct, so load all features with this id
  166.     //
  167.     string sql("select * from TableData where " + acc_col +
  168.                " in (");
  169.     string tmp;
  170.     for ( ;  id_iter.first != id_iter.second;  ++id_iter.first) {
  171.         TIdMap::iterator iter = id_iter.first;
  172.         if ( !tmp.empty() ) {
  173.             tmp += ", ";
  174.         }
  175.         tmp += "'" + iter->second + "'";
  176.     }
  177.     sql += tmp + ")";
  178.     CSQLiteTable::TIterator row_iter = m_Table->Begin(sql);
  179.     CRef<CSeq_annot> annot(new CSeq_annot());
  180.     vector<string> data;
  181.     for ( ;  *row_iter;  ++(*row_iter)) {
  182.         list<string> temp;
  183.         (*row_iter).GetRow(temp);
  184.         data.resize(temp.size());
  185.         std::copy(temp.begin(), temp.end(), data.begin());
  186.         // create a new feature
  187.         CRef<CSeq_feat> feat(new CSeq_feat());
  188.         CSeq_loc& loc = feat->SetLocation();
  189.         loc.SetInt().SetId().Assign(*idh.GetSeqId());
  190.         CUser_object& user = feat->SetData().SetUser()
  191.             .SetCategory(CUser_object::eCategory_Experiment);
  192.         CUser_object& experiment =
  193.             user.SetExperiment(CUser_object::eExperiment_Sage);
  194.         CSageData manip(experiment);
  195.         // fill in our columns
  196.         bool neg_strand = false;
  197.         TSeqPos from = 0;
  198.         TSeqPos len  = 0;
  199.         for (int i = 0;  i < data.size();  ++i) {
  200.             switch (m_ColAssign[i]) {
  201.             case eAccession:
  202.                 // already handled as ID...
  203.                 break;
  204.             case eTag:
  205.                 manip.SetTag(data[i]);
  206.                 len = data[i].length();
  207.                 break;
  208.             case eStrand:
  209.                 if (NStr::Compare(data[i], "-") == 0) {
  210.                     neg_strand = true;
  211.                 }
  212.                 break;
  213.             case eFrom:
  214.                 from = NStr::StringToInt(data[i]);
  215.                 break;
  216.             case eMethod:
  217.                 manip.SetMethod(data[i]);
  218.                 break;
  219.             case eUnknown:
  220.             default:
  221.                 manip.SetField(m_Cols[i], data[i]);
  222.                 break;
  223.             }
  224.         }
  225.         size_t to = from + len;
  226.         if (neg_strand) {
  227.             loc.SetInt().SetStrand(eNa_strand_minus);
  228.             size_t offs = to - from;
  229.             from -= offs;
  230.             to   -= offs;
  231.         } else {
  232.             loc.SetInt().SetStrand(eNa_strand_plus);
  233.         }
  234.         loc.SetInt().SetFrom(from);
  235.         loc.SetInt().SetTo  (to);
  236.         annot->SetData().SetFtable().push_back(feat);
  237.     }
  238.     CRef<CSeq_entry> entry;
  239.     if (annot) {
  240.         // we then add the object to the data loader
  241.         // we need to create a dummy TSE for it first
  242.         entry.Reset(new CSeq_entry());
  243.         entry->SetSet().SetSeq_set();
  244.         entry->SetSet().SetAnnot().push_back(annot);
  245.         GetDataSource()->AddTSE(*entry);
  246.         _TRACE("CSageDataLoader(): loaded "
  247.             << annot->GetData().GetFtable().size()
  248.             << " features for " << idh.AsString());
  249.     }
  250.     // we always save an entry here.  If the entry is empty,
  251.     // we have no information about this sequence, but we at
  252.     // least don't need to repeat an expensive search
  253.     m_Entries[idh] = entry;
  254. }
  255. END_NCBI_SCOPE
  256. /*
  257.  * ===========================================================================
  258.  * $Log: sage_dload.cpp,v $
  259.  * Revision 1000.4  2004/06/01 19:42:46  gouriano
  260.  * PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.5
  261.  *
  262.  * Revision 1.5  2004/05/21 21:42:53  gorelenk
  263.  * Added PCH ncbi_pch.hpp
  264.  *
  265.  * Revision 1.4  2003/11/28 13:41:10  dicuccio
  266.  * Fixed to match new API in CDataLoader
  267.  *
  268.  * Revision 1.3  2003/11/13 21:01:03  dicuccio
  269.  * Altered to support finding multiple accession strings matching a single handle
  270.  *
  271.  * Revision 1.2  2003/10/30 21:43:08  jcherry
  272.  * Fixed typo
  273.  *
  274.  * Revision 1.1  2003/10/02 17:40:17  dicuccio
  275.  * Initial revision
  276.  *
  277.  * ===========================================================================
  278.  */