SHOGUN  v3.2.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
IntronList.cpp
Go to the documentation of this file.
00001 
00002 #include <stdio.h>
00003 #include <string.h>
00004 
00005 #include <shogun/mathematics/Math.h>
00006 #include <shogun/lib/config.h>
00007 #include <shogun/io/SGIO.h>
00008 #include <shogun/structure/IntronList.h>
00009 
00010 using namespace shogun;
00011 
00012 CIntronList::CIntronList()
00013 :CSGObject()
00014 {
00015     m_length = 0;
00016     m_all_pos = NULL;
00017     m_intron_list = NULL;
00018     m_quality_list = NULL;
00019 }
00020 CIntronList::~CIntronList()
00021 {
00022     for (int i=0; i<m_length; i++)
00023     {
00024         SG_FREE(m_intron_list[i]);
00025         SG_FREE(m_quality_list[i]);
00026     }
00027     SG_FREE(m_intron_list);
00028     SG_FREE(m_quality_list);
00029     SG_FREE(m_all_pos);
00030 }
00031 void CIntronList::init_list(int32_t* all_pos, int32_t len)
00032 {
00033     m_length = len;
00034     m_all_pos = SG_MALLOC(int32_t, len);
00035     memcpy(m_all_pos, all_pos, len*sizeof(int32_t));
00036     m_intron_list = SG_MALLOC(int32_t*, len);
00037     m_quality_list = SG_MALLOC(int32_t*, len);
00038 
00039     //initialize all elements with an array of length one
00040     int32_t* one;
00041     for (int i=0;i<m_length;i++)
00042     {
00043         one = SG_MALLOC(int32_t, 1);//use malloc here because mem can be increased efficiently with realloc later
00044         m_intron_list[i] = one;
00045         m_intron_list[i][0] = 1;
00046         one = SG_MALLOC(int32_t, 1);
00047         m_quality_list[i] = one;
00048         m_quality_list[i][0] = 1;
00049     }
00050 }
00051 void CIntronList::read_introns(int32_t* start_pos, int32_t* end_pos, int32_t* quality, int32_t len)
00052 {
00053     int k=0;
00054     for(int i=0;i<m_length;i++)//iterate over candidate positions
00055     {
00056         while (k<len)
00057         {
00058             //SG_PRINT("i:%i, m_all_pos[i]:%i, k:%i, end_pos[k]: %i\n", i, m_all_pos[i], k, end_pos[k])
00059             if (k>0)
00060                 if (end_pos[k]<end_pos[k-1])
00061                     SG_ERROR("end pos array is not sorted: end_pos[k-1]:%i end_pos[k]:%i\n", end_pos[k-1], end_pos[k])
00062             if (end_pos[k]>=m_all_pos[i])
00063                 break;
00064             else
00065                 k++;
00066 
00067         }
00068         while (k<len && end_pos[k]==m_all_pos[i])
00069         {
00070             //SG_PRINT("\tk:%i, end_pos[k]: %i, start_pos[k]:%i\n", k, end_pos[k], start_pos[k])
00071             ASSERT(start_pos[k]<end_pos[k])
00072             ASSERT(end_pos[k]<=m_all_pos[m_length-1])
00073             // intron list
00074             //------------
00075             int32_t from_list_len = m_intron_list[i][0];
00076             int32_t* new_list = SG_REALLOC(int32_t, m_intron_list[i], from_list_len, (from_list_len+1));
00077             if (new_list == NULL)
00078                 SG_ERROR("IntronList: Out of mem 4")
00079             new_list[from_list_len]= start_pos[k];
00080             new_list[0]++;
00081             m_intron_list[i] = new_list;
00082             // quality list
00083             //--------------
00084             int32_t q_list_len = m_quality_list[i][0];
00085             //SG_PRINT("\t q_list_len:%i, from_list_len:%i \n",q_list_len, from_list_len)
00086             ASSERT(q_list_len==from_list_len)
00087             new_list = SG_REALLOC(int32_t, m_quality_list[i], q_list_len, (q_list_len+1));
00088             if (new_list == NULL)
00089                 SG_ERROR("IntronList: Out of mem 5")
00090             new_list[q_list_len]= quality[k];
00091             new_list[0]++;
00092             m_quality_list[i] = new_list;
00093 
00094             k++;
00095         }
00096     }
00097 }
00102 void CIntronList::get_intron_support(int32_t* values, int32_t from_pos, int32_t to_pos)
00103 {
00104     if (from_pos>=m_length)
00105         SG_ERROR("from_pos (%i) is not < m_length (%i)\n",to_pos, m_length)
00106     if (to_pos>=m_length)
00107         SG_ERROR("to_pos (%i) is not < m_length (%i)\n",to_pos, m_length)
00108     int32_t* from_list = m_intron_list[to_pos];
00109     int32_t* q_list = m_quality_list[to_pos];
00110 
00111     //SG_PRINT("from_list[0]: %i\n", from_list[0])
00112 
00113     int32_t coverage = 0;
00114     int32_t quality = 0;
00115     for (int i=1;i<from_list[0]; i++)
00116     {
00117         //SG_PRINT("from_list[%i]: %i, m_all_pos[from_pos]:%i\n", i,  from_list[i], m_all_pos[from_pos])
00118         if (from_list[i]==m_all_pos[from_pos])
00119         {
00120             //SG_PRINT("found intron: %i->%i\n", from_pos, to_pos )
00121             coverage = coverage+1;
00122             quality = CMath::max(quality, q_list[i]);
00123         }
00124     }
00125     values[0] = coverage;
00126     values[1] = quality;
00127 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

SHOGUN Machine Learning Toolbox - Documentation