Commit 4a4fa9be authored by Peter van 't Hof's avatar Peter van 't Hof
Browse files

Added Qual filter option

parent a5030119
...@@ -51,6 +51,7 @@ class VcfFilter(val root: Configurable) extends BiopetJavaCommandLineFunction { ...@@ -51,6 +51,7 @@ class VcfFilter(val root: Configurable) extends BiopetJavaCommandLineFunction {
object VcfFilter extends ToolCommand { object VcfFilter extends ToolCommand {
case class Args(inputVcf: File = null, case class Args(inputVcf: File = null,
outputVcf: File = null, outputVcf: File = null,
minQualscore: Option[Double] = None,
minSampleDepth: Int = -1, minSampleDepth: Int = -1,
minTotalDepth: Int = -1, minTotalDepth: Int = -1,
minAlternateDepth: Int = -1, minAlternateDepth: Int = -1,
...@@ -105,6 +106,9 @@ object VcfFilter extends ToolCommand { ...@@ -105,6 +106,9 @@ object VcfFilter extends ToolCommand {
opt[Unit]("filterNoCalls") unbounded () action { (x, c) => opt[Unit]("filterNoCalls") unbounded () action { (x, c) =>
c.copy(filterNoCalls = true) c.copy(filterNoCalls = true)
} text ("Filter when there are only no calls") } text ("Filter when there are only no calls")
opt[Double]("minQualscore") unbounded () action { (x, c) =>
c.copy(minQualscore = Some(x))
} text ("Min qual score")
} }
var commandArgs: Args = _ var commandArgs: Args = _
...@@ -122,7 +126,8 @@ object VcfFilter extends ToolCommand { ...@@ -122,7 +126,8 @@ object VcfFilter extends ToolCommand {
writer.writeHeader(header) writer.writeHeader(header)
for (record <- reader) { for (record <- reader) {
if (filterRefCalls(record) && if (minQualscore(record) &&
filterRefCalls(record) &&
filterNoCalls(record) && filterNoCalls(record) &&
minTotalDepth(record) && minTotalDepth(record) &&
minSampleDepth(record) && minSampleDepth(record) &&
...@@ -139,6 +144,11 @@ object VcfFilter extends ToolCommand { ...@@ -139,6 +144,11 @@ object VcfFilter extends ToolCommand {
writer.close writer.close
} }
def minQualscore(record: VariantContext): Boolean = {
if (commandArgs.minQualscore.isEmpty) return true
record.getPhredScaledQual >= commandArgs.minQualscore.get
}
def filterRefCalls(record: VariantContext): Boolean = { def filterRefCalls(record: VariantContext): Boolean = {
if (commandArgs.filterNoCalls) record.getGenotypes.exists(g => !g.isHomRef) if (commandArgs.filterNoCalls) record.getGenotypes.exists(g => !g.isHomRef)
else true else true
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment