Resolve "callPeaksMACS exits when data is bad"

Please fill in the appropriate checklist below (delete whatever is not relevant). These are the most common things requested on pull requests (PRs).

PR checklist

• This comment contains a description of changes (with reason)
• If you've fixed a bug or added code that should be tested, add tests!
• CHANGELOG.md is updated

Closes #26 (closed)

workflow/scripts/call_peaks_macs.py

@@ -138,8 +138,20 @@ def call_peaks_macs(experiment, xcor, control, prefix, genome_size, chrom_sizes)
     with open(xcor, 'r') as xcor_fh:
         firstline = xcor_fh.readline()
         frag_lengths = firstline.split()[2]  # third column
-        fragment_length = frag_lengths.split(',')[0]  # grab first value
-        logger.info("Fraglen %s", fragment_length)
+        frag_lengths_array = frag_lengths.split(',')
+        fragment_length = 0
+        fragment = False
+        # Loop through all values of fragment length
+        for f in frag_lengths.split(','):
+            fragment_length = f
+            logger.info("Fraglen %s", fragment_length)
+            if int(fragment_length) > 50:
+                fragment = True
+                break
+
+        if fragment == False:
+            logger.info('Error in cross-correlation analysis: %s', frag_lengths_array)
+            raise Exception("Error in cross-correlation analysis: %s" % frag_lengths_array)
 
     # Generate narrow peaks and preliminary signal tracks
 

workflow/tests/test_call_peaks_macs.py

@@ -3,10 +3,29 @@
 import pytest
 import os
 import utils
+from io import StringIO
+import call_peaks_macs
 
 test_output_path = os.path.dirname(os.path.abspath(__file__)) + \
                 '/../output/callPeaksMACS/'
 
+test_data_path = os.path.dirname(os.path.abspath(__file__)) + \
+		'/../../test_data/'
+        
+
+@pytest.mark.unit
+def test_fragment_length():
+    experiment = "experiment.tagAlign.gz"
+    control = "control.tagAlign.gz"
+    prefix = 'test'
+    genome_size = 'hs'
+    chrom_sizes = 'genomefile.txt'
+    cross_qc = os.path.join(test_data_path, 'test_cross.qc')
+    with pytest.raises(Exception) as excinfo:
+        call_peaks_macs.call_peaks_macs(experiment, cross_qc, control, prefix, genome_size, chrom_sizes)
+    assert str(excinfo.value) == "Error in cross-correlation analysis: ['0', '20', '33']"
+
+
 
 @pytest.mark.singleend
 def test_fc_signal_singleend():

.gitignore

@@ -108,5 +108,4 @@ report*.html*
 timeline*.html*
 /workflow/output/*
 /work/*
-/test_data/*
 /.nextflow/*

CHANGELOG.md

@@ -14,8 +14,10 @@ All notable changes to this project will be documented in this file.
 - Make gtf and geneName files as param inputs
 - Fix xcor to increase file size for --random-source
 - Fix skip diff test for paired-end data
+- Fix xcor to get lowest non zero value about 50
 - Fix references to display in Multiqc report
 
+
 ## [publish_1.0.6 ] - 2019-05-31
 ### Added
 - MIT LICENSE