etc/gatk-wdl/fm2/tasks/fgi.wdl

   1 version 1.0
   2
   3 task FilterSamReads {
   4         input {
   5                 File inputBam
   6                 File inputBamIndex
   7                 String outputBamPath
   8                 Int minMatchLen = 200
   9                 Boolean? createMd5File = false
  10
  11                 Int? compressionLevel
  12
  13                 Int? javaXmxMb = 1024
  14                 Int? memoryMb = javaXmxMb + 512
  15                 # One minute per input gigabyte.
  16                 Int? timeMinutes = 1 + ceil(size(inputBam, "G") * 1)
  17                 String? dockerImage = "quay.io/biocontainers/picard:2.23.8--0"
  18         }
  19         Array[String] FilterScriptContents = [
  20                 'function accept(rec) {',
  21                 '  if (rec.getReadUnmappedFlag()) return false;',
  22                 '  var cigar = rec.getCigar();',
  23                 '  if (cigar == null) return false;',
  24                 '  var readMatch = 0;',
  25                 '  for (var i=0;i < cigar.numCigarElements();++i) {',
  26                 '    var ce = cigar.getCigarElement(i);',
  27                 '    if (ce.getOperator().name() == "M") readMatch += ce.length;',
  28                 '  }',
  29                 '  if (readMatch > '+ minMatchLen +') return true;',
  30                 '}',
  31                 'accept(record);'
  32         ]
  33
  34         command {
  35                 set -e
  36                 mkdir -p "$(dirname ~{outputBamPath})"
  37                 JAVA_OPTS="-Xmx~{javaXmxMb}M -XX:ParallelGCThreads=1" picard \
  38                 FilterSamReads \
  39                 INPUT=~{sep=' INPUT=' inputBam} \
  40                 OUTPUT=~{outputBamPath} \
  41                 ~{"COMPRESSION_LEVEL=" + compressionLevel} \
  42                 JAVASCRIPT_FILE=${write_lines(FilterScriptContents)} \
  43                 FILTER=includeJavascript \
  44                 CREATE_INDEX=true \
  45                 CREATE_MD5_FILE=~{true="true" false="false" createMd5File}
  46         }
  47
  48         output {
  49                 File outputBam = outputBamPath
  50                 File outputBamIndex = sub(outputBamPath, "\.bam$", ".bai")
  51                 File? outputBamMd5 = outputBamPath + ".md5"
  52         }
  53
  54         runtime {
  55                 memory: "~{memoryMb}M"
  56                 time_minutes: timeMinutes
  57                 #docker: dockerImage
  58         }
  59 }
  60
  61 task callSNP {
  62         input {
  63                 File inputBam
  64                 File inputBamIndex
  65                 GatkIndex GatkIndex
  66                 String outputPath
  67                 File helperPl
  68         }
  69         File referenceFasta = GatkIndex.fastaFile
  70         String bcfFile = outputPath + "/mpileup.bcf"
  71         String snpFile = outputPath + "/snp.gz"
  72         command {
  73                 set -e
  74                 mkdir -p "~{outputPath}"
  75                 bcftools mpileup --threads 6 ~{inputBam} -d 30000 -Q 30 -f ~{referenceFasta} -p -Ob -o ~{bcfFile}
  76                 bcftools call -Oz -A -m ~{bcfFile} -o ~{snpFile}
  77                 bcftools index ~{snpFile}
  78                 bcftools query -f'%CHROM\t[%DP\t%QUAL\t%TGT\n]' -i 'POS==501' ~{snpFile} > ~{outputPath + "/snp0.txt"}
  79                 perl ~{helperPl} ~{outputPath + "/snp0.txt"} > ~{outputPath + "/../snp.txt"}
  80         }
  81
  82         output {
  83                 File outSNP0txt = outputPath + "/../snp0.txt"
  84                 File outSNPtxt = outputPath + "/../snp.txt"
  85         }
  86 }
  87
  88 task callSTR {
  89         input {
  90                 File inputBam
  91                 File inputBamIndex
  92                 String outputPath
  93                 File helperBED
  94                 File helperPl
  95         }
  96         String awkStr = "'\{print $1\"\\t\"$4\}'"
  97         command {
  98                 set -e
  99                 mkdir -p "~{outputPath}"
 100                 samtools mpileup -l ~{helperBED} ~{inputBam} |awk ~{awkStr} > ~{outputPath + "/str0.txt"}
 101                 perl ~{helperPl} ~{outputPath + "/str0.txt"} > ~{outputPath + "/../str.txt"}
 102         }
 103         output {
 104                 File outSTR0txt = outputPath + "/../str0.txt"
 105                 File outSTRtxt = outputPath + "/../str.txt"
 106         }
 107 }
 108
 109 struct GatkIndex {
 110     File fastaFile
 111     Array[File] indexFiles
 112 }