.get_chromosome_rle

Author: zhengxwen

NAMESPACE

@@ -9,6 +9,7 @@ useDynLib(SeqArray,
     SEQ_SetSpaceChrom, SEQ_SetSpaceAnnotID,
     SEQ_SplitSelection, SEQ_SplitSelectionX,
     SEQ_GetSpace, SEQ_GetSpaceSample, SEQ_GetSpaceVariant,
+    SEQ_Get_ChromRLE,
     SEQ_Summary, SEQ_System,
     SEQ_VCF_NumLines, SEQ_VCF_Split, SEQ_VCF_Parse,
     SEQ_ToVCF_Init, SEQ_ToVCF_Done,

R/Internal.R

@@ -255,6 +255,17 @@ process_count <- 1L
 
 
 
+#######################################################################
+# Get RLE-coded chromosome
+#
+.get_chromosome_rle <- function(gdsfile)
+{
+    v <- .Call(SEQ_Get_ChromRLE, gdsfile)
+    structure(list(lengths=v[[1L]], values=v[[2L]]), class="rle")
+}
+
+
+
 #######################################################################
 # Open and close a connection,
 # Please always call '.close_conn' after '.open_bin' and '.open_text'

src/Index.h

@@ -76,6 +76,8 @@ class ErrSeqArray;
 template<typename TYPE> class COREARRAY_DLL_LOCAL C_RLE
 {
 public:
+	friend class CChromIndex;
+
 	/// constructor
 	C_RLE()
 	{
@@ -249,25 +251,26 @@ class COREARRAY_DLL_LOCAL CChromIndex
 
 	typedef vector<TRange> TRangeList;
 
+	/// map to TRangeList from chromosome coding
+	map<string, TRangeList> Map;
+
 	/// constructor
 	CChromIndex();
-
 	/// clear
 	void Clear();
-
 	/// represent chromosome codes as a RLE object in Map
 	void AddChrom(PdGDSFolder Root);
-
 	/// the total length of a TRangeList object
 	size_t RangeTotalLength(const TRangeList &RngList);
 
 	/// whether it is empty
 	inline bool Empty() const { return Map.empty(); }
 	/// return chromosome coding with the index 
 	inline const string &operator [](size_t idx) { return RleChr[idx]; }
-
-	/// map to TRangeList from chromosome coding
-	map<string, TRangeList> Map;
+	/// return the RLE representation of chromosome: value
+	inline vector<string> &RLE_Values() { return RleChr.Values; }
+	/// return the RLE representation of chromosome: length
+	inline vector<C_UInt32> &RLE_Lengths() { return RleChr.Lengths; }
 
 protected:
 	/// indexing chromosome

src/SeqArray.cpp

@@ -1363,7 +1363,7 @@ COREARRAY_DLL_EXPORT SEXP SEQ_ClearVarMap(SEXP gdsfile)
 
 
 // ===========================================================
-// Clear VarMap in a GDS file
+// Get or clear the memory buffer storing variant positions
 // ===========================================================
 
 COREARRAY_DLL_EXPORT SEXP SEQ_BufferPosition(SEXP gdsfile, SEXP clear)
@@ -1385,6 +1385,36 @@ COREARRAY_DLL_EXPORT SEXP SEQ_BufferPosition(SEXP gdsfile, SEXP clear)
 
 
 
+// ===========================================================
+// Get RLE-coded chromosome
+// ===========================================================
+
+COREARRAY_DLL_EXPORT SEXP SEQ_Get_ChromRLE(SEXP gdsfile)
+{
+	COREARRAY_TRY
+		CFileInfo &File = GetFileInfo(gdsfile);
+		CChromIndex &Chrom = File.Chromosome();
+		vector<string> &V = Chrom.RLE_Values();
+		vector<C_UInt32> &L = Chrom.RLE_Lengths();
+		if (V.size() != L.size())
+			throw ErrSeqArray("Invalid RLE of chromosome!");
+		rv_ans = PROTECT(NEW_LIST(2));
+		SEXP l = NEW_INTEGER(L.size());
+		SET_ELEMENT(rv_ans, 0, l);
+		memcpy(INTEGER(l), &L[0], sizeof(int)*L.size());
+		SEXP v = NEW_CHARACTER(V.size());
+		SET_ELEMENT(rv_ans, 1, v);
+		for (size_t i=0; i < V.size(); i++)
+		{
+			string &s = V[i];
+			SET_STRING_ELT(v, i, Rf_mkCharLen(&s[0], s.size()));
+		}
+		UNPROTECT(1);
+	COREARRAY_CATCH
+}
+
+
+
 // ===========================================================
 // Get system configuration
 // ===========================================================
@@ -1677,6 +1707,7 @@ COREARRAY_DLL_EXPORT void R_init_SeqArray(DllInfo *info)
 
 		CALL(SEQ_IntAssign, 2),             CALL(SEQ_AppendFill, 3),
 		CALL(SEQ_ClearVarMap, 1),           CALL(SEQ_BufferPosition, 2),
+		CALL(SEQ_Get_ChromRLE, 1),
 
 		CALL(SEQ_bgzip_create, 1),
 		CALL(SEQ_ToVCF_Init, 6),            CALL(SEQ_ToVCF_Done, 0),