QcCommand.scala 4.79 KB
Newer Older
Peter van 't Hof's avatar
Peter van 't Hof committed
1
2
3
4
package nl.lumc.sasc.biopet.pipelines.flexiprep

import java.io.File

Peter van 't Hof's avatar
Peter van 't Hof committed
5
import nl.lumc.sasc.biopet.core.summary.{ SummaryQScript, Summarizable }
Peter van 't Hof's avatar
Peter van 't Hof committed
6
import nl.lumc.sasc.biopet.core.{ BiopetFifoPipe, BiopetCommandLineFunction, BiopetPipe }
Peter van 't Hof's avatar
Peter van 't Hof committed
7
import nl.lumc.sasc.biopet.extensions.{ Cat, Gzip, Sickle, Cutadapt }
Peter van 't Hof's avatar
Peter van 't Hof committed
8
import nl.lumc.sasc.biopet.extensions.seqtk.SeqtkSeq
Peter van 't Hof's avatar
Peter van 't Hof committed
9
10
11
12
13
14
import nl.lumc.sasc.biopet.utils.config.Configurable
import org.broadinstitute.gatk.utils.commandline.{ Output, Input }

/**
 * Created by pjvan_thof on 9/22/15.
 */
Peter van 't Hof's avatar
Peter van 't Hof committed
15
class QcCommand(val root: Configurable, val fastqc: Fastqc) extends BiopetCommandLineFunction with Summarizable {
Peter van 't Hof's avatar
Peter van 't Hof committed
16
17
18
19
20
21
22
23
24
25
26
27

  val flexiprep = root match {
    case f: Flexiprep => f
    case _            => throw new IllegalArgumentException("This class may only be used inside Flexiprep")
  }

  @Input(required = true)
  var input: File = _

  @Output(required = true)
  var output: File = _

Peter van 't Hof's avatar
Peter van 't Hof committed
28
29
30
31
  var compress = true

  var read: String = _

Peter van 't Hof's avatar
Peter van 't Hof committed
32
33
34
35
36
37
  override def defaultCoreMemory = 2.0
  override def defaultThreads = 3

  val seqtk = new SeqtkSeq(root)
  var clip: Option[Cutadapt] = None
  var trim: Option[Sickle] = None
Peter van 't Hof's avatar
Peter van 't Hof committed
38
39
40
  var outputCommand: BiopetCommandLineFunction = null

  def jobs = (Some(seqtk) :: clip :: trim :: Some(outputCommand) :: Nil).flatten
Peter van 't Hof's avatar
Peter van 't Hof committed
41

Peter van 't Hof's avatar
Peter van 't Hof committed
42
43
44
45
46
47
48
  def summaryFiles = Map()

  def summaryStats = Map()

  override def addToQscriptSummary(qscript: SummaryQScript, name: String): Unit = {
    clip match {
      case Some(job) => qscript.addSummarizable(job, s"clipping_$read")
Peter van 't Hof's avatar
Peter van 't Hof committed
49
      case _         =>
Peter van 't Hof's avatar
Peter van 't Hof committed
50
51
52
    }
    trim match {
      case Some(job) => qscript.addSummarizable(job, s"trimming_$read")
Peter van 't Hof's avatar
Peter van 't Hof committed
53
      case _         =>
Peter van 't Hof's avatar
Peter van 't Hof committed
54
55
56
    }
  }

Peter van 't Hof's avatar
Peter van 't Hof committed
57
58
  override def beforeGraph(): Unit = {
    super.beforeGraph()
Peter van 't Hof's avatar
Peter van 't Hof committed
59
    require(read != null)
Peter van 't Hof's avatar
Peter van 't Hof committed
60
    deps ::= input
Peter van 't Hof's avatar
Peter van 't Hof committed
61
    outputFiles :+= output
Peter van 't Hof's avatar
Peter van 't Hof committed
62
63
  }

Peter van 't Hof's avatar
Peter van 't Hof committed
64
  override def beforeCmd(): Unit = {
Peter van 't Hof's avatar
Peter van 't Hof committed
65
    seqtk.input = input
Peter van 't Hof's avatar
Peter van 't Hof committed
66
    seqtk.output = new File(output.getParentFile, input.getName + ".seqtk.fq")
Peter van 't Hof's avatar
Peter van 't Hof committed
67
68
69
70
71
72
73
74
75
76
    seqtk.Q = fastqc.encoding match {
      case null => None
      case enc if enc.contains("Sanger / Illumina 1.9") => None
      case enc if enc.contains("Illumina <1.3") => Option(64)
      case enc if enc.contains("Illumina 1.3") => Option(64)
      case enc if enc.contains("Illumina 1.5") => Option(64)
      case _ => None
    }
    if (seqtk.Q.isDefined) seqtk.V = true

Peter van 't Hof's avatar
Peter van 't Hof committed
77
    clip = if (!flexiprep.skipClip) {
Peter van 't Hof's avatar
Peter van 't Hof committed
78
79
      val foundAdapters = fastqc.foundAdapters.map(_.seq)
      if (foundAdapters.nonEmpty) {
Peter van 't Hof's avatar
Peter van 't Hof committed
80
81
82
        val cutadept = new Cutadapt(root)
        cutadept.fastq_input = seqtk.output
        cutadept.fastq_output = new File(output.getParentFile, input.getName + ".cutadept.fq")
Peter van 't Hof's avatar
Peter van 't Hof committed
83
        cutadept.stats_output = new File(flexiprep.outputDir, s"${flexiprep.sampleId.getOrElse("x")}-${flexiprep.libId.getOrElse("x")}.$read.clip.stats")
Peter van 't Hof's avatar
Peter van 't Hof committed
84
85
86
        if (cutadept.default_clip_mode == "3") cutadept.opt_adapter ++= foundAdapters
        else if (cutadept.default_clip_mode == "5") cutadept.opt_front ++= foundAdapters
        else if (cutadept.default_clip_mode == "both") cutadept.opt_anywhere ++= foundAdapters
Peter van 't Hof's avatar
Peter van 't Hof committed
87
88
89
90
        Some(cutadept)
      } else None
    } else None

Peter van 't Hof's avatar
Peter van 't Hof committed
91
    trim = if (!flexiprep.skipTrim) {
Peter van 't Hof's avatar
Peter van 't Hof committed
92
      val sickle = new Sickle(root)
Peter van 't Hof's avatar
Peter van 't Hof committed
93
      sickle.output_stats = new File(flexiprep.outputDir, s"${flexiprep.sampleId.getOrElse("x")}-${flexiprep.libId.getOrElse("x")}.$read.trim.stats")
Peter van 't Hof's avatar
Peter van 't Hof committed
94
95
96
97
98
      sickle.input_R1 = clip match {
        case Some(clip) => clip.fastq_output
        case _          => seqtk.output
      }
      sickle.output_R1 = new File(output.getParentFile, input.getName + ".sickle.fq")
Peter van 't Hof's avatar
Peter van 't Hof committed
99
100
      Some(sickle)
    } else None
Peter van 't Hof's avatar
Peter van 't Hof committed
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129

    val outputFile = (clip, trim) match {
      case (_, Some(trim)) => trim.output_R1
      case (Some(clip), _) => clip.fastq_output
      case _               => seqtk.output
    }

    if (compress) outputCommand = {
      val gzip = new Gzip(root)
      gzip.input = outputFile :: Nil
      gzip.output = output
      gzip
    }
    else outputCommand = {
      val cat = new Cat(root)
      cat.input = outputFile :: Nil
      cat.output = output
      cat
    }

    seqtk.beforeGraph()
    clip.foreach(_.beforeGraph())
    trim.foreach(_.beforeGraph())
    outputCommand.beforeGraph()

    seqtk.beforeCmd()
    clip.foreach(_.beforeCmd())
    trim.foreach(_.beforeCmd())
    outputCommand.beforeCmd()
Peter van 't Hof's avatar
Peter van 't Hof committed
130
131
132
  }

  def cmdLine = {
Peter van 't Hof's avatar
Peter van 't Hof committed
133
134

    val cmd = (clip, trim) match {
Peter van 't Hof's avatar
Peter van 't Hof committed
135
136
137
138
      case (Some(clip), Some(trim)) => new BiopetFifoPipe(root, seqtk :: clip :: trim :: outputCommand :: Nil)
      case (Some(clip), _)          => new BiopetFifoPipe(root, seqtk :: clip :: outputCommand :: Nil)
      case (_, Some(trim))          => new BiopetFifoPipe(root, seqtk :: trim :: outputCommand :: Nil)
      case _                        => new BiopetFifoPipe(root, seqtk :: outputCommand :: Nil)
Peter van 't Hof's avatar
Peter van 't Hof committed
139
    }
Peter van 't Hof's avatar
Peter van 't Hof committed
140
141

    //val cmds = (Some(seqtk) :: clip :: trim :: Some(new Gzip(root)) :: Nil).flatten
Peter van 't Hof's avatar
Peter van 't Hof committed
142
143
144
145
    cmd.beforeGraph()
    cmd.commandLine
  }
}