update to allow per-contig trim setting in utgcns

marbl · Sep 24, 2024 · 39ef0f2 · 39ef0f2
1 parent 0d51e7d
commit 39ef0f2
Show file tree

Hide file tree

Showing 4 changed files with 17 additions and 5 deletions.
diff --git a/src/stores/tgTig.C b/src/stores/tgTig.C
@@ -37,6 +37,7 @@ tgTig::saveToRecord(tgTigRecord &tr) {
   tr._suggestBubble       = _suggestBubble;
   tr._suggestCircular     = _suggestCircular;
   tr._circularLength      = _circularLength;
+  tr._suggestNoTrim       = _suggestNoTrim;
   tr._spare               = _spare;
 
   tr._trimBgn             = _trimBgn;
@@ -59,6 +60,7 @@ tgTig::restoreFromRecord(tgTigRecord &tr) {
   _suggestBubble       = tr._suggestBubble;
   _suggestCircular     = tr._suggestCircular;
   _circularLength      = tr._circularLength;
+  _suggestNoTrim       = tr._suggestNoTrim;
   _spare               = tr._spare;
 
   _trimBgn             = tr._trimBgn;
@@ -116,6 +118,7 @@ tgTig::operator=(tgTig & tg) {
   _suggestBubble       = tg._suggestBubble;
   _suggestCircular     = tg._suggestCircular;
   _circularLength      = tg._circularLength;
+  _suggestNoTrim       = tg._suggestNoTrim;
   _spare               = tg._spare;
 
   _layoutLen           = tg._layoutLen;
@@ -178,6 +181,7 @@ tgTig::clear(void) {
   _suggestBubble        = 0;
   _suggestCircular      = 0;
   _circularLength       = 0;
+  _suggestNoTrim        = 0;
   _spare                = 0;
 
   _layoutLen            = 0;

diff --git a/src/stores/tgTig.H b/src/stores/tgTig.H
@@ -226,7 +226,8 @@ struct tgTigRecordV1 {
   uint32              _suggestCircular : 1;
   uint32              _suggestBubble   : 1;
   uint32              _circularLength  : 20;
-  uint32              _spare           : 7;
+  uint32              _suggestNoTrim   : 1;
+  uint32              _spare           : 6;
 
   uint32              _layoutLen;
   uint32              _basesLen;
@@ -247,7 +248,8 @@ struct tgTigRecord {
   uint32              _suggestCircular : 1;
   uint32              _suggestBubble   : 1;
   uint32              _circularLength  : 20;
-  uint32              _spare           : 7;
+  uint32              _suggestNoTrim   : 1;
+  uint32              _spare           : 6;
 
   uint32              _trimBgn;
   uint32              _trimEnd;
@@ -267,6 +269,7 @@ struct tgTigRecord {
     _suggestCircular   = that._suggestCircular;
     _suggestBubble     = that._suggestBubble;
     _circularLength    = that._circularLength;
+    _suggestNoTrim     = that._suggestNoTrim;
     _spare             = that._spare;
     _trimBgn           = that._trimBgn;
     _trimEnd           = that._trimEnd;
@@ -414,7 +417,8 @@ public:
   uint32       _suggestCircular : 1  = false;          //  Bogart found overlaps making a circle.
   uint32       _suggestBubble   : 1  = false;          //  Bogart thinks this is a bubble.
   uint32       _circularLength  : 20 = 0;              //  Length of the end-to-end circularizing overlap.
-  uint32       _spare           : 7  = 0;
+  uint32       _suggestNoTrim   : 1 = false;           //  If set to true, ignore any trim value specified on the command line to consensus and do no trimming
+  uint32       _spare           : 6  = 0;
 
   uint32       _layoutLen = 0;           //  The max coord in the layout.  Same as _basesLen if _bases exists.
 

diff --git a/src/utgcns/layoutToPackage.C b/src/utgcns/layoutToPackage.C
@@ -142,6 +142,10 @@ loadVerkkoLayouts(sqCache              *reads,
       nReads = strtouint32(W[1]);
     }
 
+	else if (strcmp(W[0], "trm") == 0) {
+	   tig->_suggestNoTrim = (strtouint32(W[1]) != 0);
+    }
+
     else if (strcmp(W[0], "end") == 0) {
       if (nReads != 0)
         fprintf(stderr, "ERROR: Tig '%d' reads doesn't match number expected\n", tig->tigID()), err++;

diff --git a/src/utgcns/unitigConsensus.C b/src/utgcns/unitigConsensus.C
@@ -474,7 +474,7 @@ unitigConsensus::generateTemplateStitch(void) {
     bool   moreToExtend  = (readEnd < readLen);
 
 
-    int32 maxDifference = std::min(MAX_PADDING*10, (int32)ceil(0.30*olapLen));
+    int32 maxDifference = std::min(MAX_PADDING, (int32)ceil(0.30*olapLen));
     //  Reset if the edit distance is waay more than our error rate allows or it's very short and we haven't topped out on error.  This seems to be a quirk with
     //  edlib when aligning to N's - I got startLocation = endLocation = 0 and editDistance = alignmentLength.
     if ((double)result.editDistance / result.alignmentLength > bandErrRate || 
@@ -940,7 +940,7 @@ unitigConsensus::generatePBDAG(char aligner_, u32toRead &reads_) {
   _templateToCNS  = new uint32 [tiglen + 1];
   _templateLength = tiglen;
 
-  std::string cns = ag.consensusNoSplit(_minCoverage, _templateToCNS, _templateLength);
+  std::string cns = ag.consensusNoSplit((_tig->_suggestNoTrim == 0 ? _minCoverage : 0), _templateToCNS, _templateLength);
 
   delete [] tigseq;