Debug report option for assembleContigs.

src/main/java/com/milaboratory/mixcr/assembler/fullseq/FullSeqAssembler.java

@@ -1037,7 +1037,7 @@ OutputPort<int[]> createPort() {
      * Represents aggregated information about nucleotide states for all positions in all reads aggregated with
      * {@link #calculateRawData(Supplier)}.
      */
-    public static abstract class RawVariantsData {
+    public abstract class RawVariantsData {
         /**
          * Total number of reads
          */
@@ -1064,6 +1064,81 @@ public static abstract class RawVariantsData {
          */
         void destroy() {
         }
+
+        /**
+         * String representation of this state matrix
+         *
+         * @param qualityThreshold quality threshold (positions with quality lower then this value, wil be printed in lower case)
+         * @param readsFrom        range of read ids to print (beginning; inclusive)
+         * @param readsTo          range of read ids to print (end; exclusive)
+         */
+        public String toString(byte qualityThreshold, int readsFrom, int readsTo) {
+            int minPosition = Arrays.stream(points).min().getAsInt();
+            int maxPosition = Arrays.stream(points).max().getAsInt() + 1;
+
+            // Calculating maximal observed sequence length for each of the positions
+            int[] len = new int[maxPosition - minPosition];
+            OutputPort<int[]> port = createPort();
+            for (int position : points) {
+                int[] states = port.take();
+                for (int j = readsFrom; j < readsTo; j++) {
+                    int state = states[j];
+                    if (state != ABSENT_PACKED_VARIANT_INFO)
+                        len[position - minPosition] = Math.max(len[position - minPosition], variantIdToSequence.get(state >>> 8).size());
+                }
+            }
+            assert port.take() == null;
+
+            // Calculating position projection
+            int[] positionMap = new int[len.length];
+            for (int i = 1; i < len.length; i++)
+                positionMap[i] = positionMap[i - 1] + len[i - 1];
+            int maxLength = positionMap[len.length - 1] + len[len.length - 1];
+
+            // Allocating main array
+            char[][] result = new char[readsTo - readsFrom][maxLength];
+            for (char[] line : result)
+                Arrays.fill(line, ' ');
+
+            port = createPort();
+            for (int position : points) {
+                int[] states = port.take();
+                for (int j = readsFrom; j < readsTo; j++) {
+                    int state = states[j];
+                    if (state == ABSENT_PACKED_VARIANT_INFO)
+                        continue;
+                    String seq = variantIdToSequence.get(state >>> 8).toString();
+                    if ((state & 0x7F) < qualityThreshold)
+                        seq = seq.toLowerCase();
+                    for (int k = 0; k < len[position - minPosition]; k++)
+                        if (k < seq.length())
+                            result[j + readsFrom][positionMap[position - minPosition] + k] = seq.charAt(k);
+                        else
+                            result[j + readsFrom][positionMap[position - minPosition] + k] = '.';
+                }
+            }
+            assert port.take() == null;
+
+            return Arrays.stream(result)
+                    .map(String::new)
+                    .collect(Collectors.joining("\n"));
+        }
+
+        /**
+         * String representation of this state matrix
+         *
+         * @param qualityThreshold quality threshold (positions with quality lower then this value, wil be printed in lower case)
+         */
+        public String toString(byte qualityThreshold) {
+            return toString(qualityThreshold, 0, nReads);
+        }
+
+        /**
+         * String representation of this state matrix
+         */
+        public String toString() {
+            return toString((byte) 10);
+        }
     }
 
     private static final class VariantAggregator {

src/main/java/com/milaboratory/mixcr/cli/ActionAssembleContigs.java

@@ -58,7 +58,8 @@ public void go(ActionHelper helper) throws Exception {
         VDJCAlignerParameters alignerParameters;
         CloneAssemblerParameters cloneAssemblerParameters;
         try (ClnAReader reader = new ClnAReader(parameters.getInputFileName(), VDJCLibraryRegistry.getDefault());
-             PrimitivO tmpOut = new PrimitivO(new BufferedOutputStream(new FileOutputStream(parameters.getOutputFileName())))) {
+             PrimitivO tmpOut = new PrimitivO(new BufferedOutputStream(new FileOutputStream(parameters.getOutputFileName())));
+             BufferedWriter debugReport = parameters.debugReport == null ? null : new BufferedWriter(new OutputStreamWriter(new FileOutputStream(parameters.debugReport)))) {
 
             final CloneFactory cloneFactory = new CloneFactory(reader.getAssemblerParameters().getCloneFactoryParameters(),
                     reader.getAssemblingFeatures(), reader.getGenes(), reader.getAlignerParameters().getFeaturesToAlignMap());
@@ -82,6 +83,24 @@ public void go(ActionHelper helper) throws Exception {
                         throw new RuntimeException(e);
                     }
                 });
+
+                if (debugReport != null) {
+                    synchronized (debugReport) {
+                        try {
+                            debugReport.write("Clone: " + cloneAlignments.clone.getId());
+                            debugReport.newLine();
+                            debugReport.write(rawVariantsData.toString());
+                            debugReport.newLine();
+                            debugReport.newLine();
+                            debugReport.write("==========================================");
+                            debugReport.newLine();
+                            debugReport.newLine();
+                        } catch (IOException e) {
+                            throw new RuntimeException(e);
+                        }
+                    }
+                }
+
                 return fullSeqAssembler.callVariants(rawVariantsData);
             }, parameters.threads);
 
@@ -161,6 +180,10 @@ public static class FullSeqParameters extends ActionParametersWithOutput {
                 names = {"-r", "--report"})
         public String report;
 
+        @Parameter(description = "Report file.",
+                names = {"--debug-report"}, hidden = true)
+        public String debugReport;
+
         @Parameter(description = "JSON report file.",
                 names = {"--json-report"})
         public String jsonReport = null;

src/test/java/com/milaboratory/mixcr/assembler/fullseq/FullSeqAssemblerTest.java

@@ -37,7 +37,7 @@ public class FullSeqAssemblerTest {
     static final FullSeqAssemblerParameters DEFAULT_PARAMETERS =
             new FullSeqAssemblerParameters(0.1, 80, 120,
                     3, 7, 0.25, GeneFeature.VDJRegion,
-                    new QualityTrimmerParameters(20.0f, 8),20,  false);
+                    new QualityTrimmerParameters(20.0f, 8), 20, false);
 
     static final class MasterSequence {
         final int vPart, cdr3Part, jPart, cPart;
@@ -445,6 +445,7 @@ public void testLargeCloneNoMismatches() throws Exception {
                 align.usedGenes, align.parameters.alignerParameters.getFeaturesToAlignMap());
         FullSeqAssembler agg = new FullSeqAssembler(cloneFactory, DEFAULT_PARAMETERS, initialClone, align.parameters.alignerParameters);
         FullSeqAssembler.RawVariantsData prep = agg.calculateRawData(() -> CUtils.asOutputPort(alignments));
+        // System.out.println(prep);
         List<Clone> clones = new ArrayList<>(new CloneSet(Arrays.asList(agg.callVariants(prep))).getClones());
         clones.sort(Comparator.comparingDouble(Clone::getCount).reversed());
         for (Clone clone : clones) {