Mercurial > repos > timpalpant > java_genomics_toolkit
diff java-genomics-toolkit/src/edu/unc/genomics/ngs/BaseAlignCounts.java @ 0:1daf3026d231
Upload alpha version
author | timpalpant |
---|---|
date | Mon, 13 Feb 2012 21:55:55 -0500 |
parents | |
children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/java-genomics-toolkit/src/edu/unc/genomics/ngs/BaseAlignCounts.java Mon Feb 13 21:55:55 2012 -0500 @@ -0,0 +1,94 @@ +package edu.unc.genomics.ngs; + +import java.io.BufferedWriter; +import java.io.IOException; +import java.nio.charset.Charset; +import java.nio.file.Files; +import java.nio.file.Path; +import java.util.Iterator; + +import org.apache.log4j.Logger; + +import com.beust.jcommander.Parameter; + +import edu.ucsc.genome.TrackHeader; +import edu.unc.genomics.Assembly; +import edu.unc.genomics.CommandLineTool; +import edu.unc.genomics.Interval; +import edu.unc.genomics.io.IntervalFile; + +public class BaseAlignCounts extends CommandLineTool { + + private static final Logger log = Logger.getLogger(BaseAlignCounts.class); + + @Parameter(names = {"-i", "--input"}, description = "Input file (reads)", required = true) + public IntervalFile<? extends Interval> intervalFile; + @Parameter(names = {"-a", "--assembly"}, description = "Genome assembly", required = true) + public Assembly assembly; + @Parameter(names = {"-x", "--extend"}, description = "Extend reads from 5' end (default = read length)") + public Integer extend; + @Parameter(names = {"-o", "--output"}, description = "Output file (Wig)", required = true) + public Path outputFile; + + @Override + public void run() throws IOException { + log.debug("Initializing output file"); + int mapped = 0; + try (BufferedWriter writer = Files.newBufferedWriter(outputFile, Charset.defaultCharset())) { + // Write the Wiggle track header to the output file + TrackHeader header = new TrackHeader("wiggle_0"); + header.setName("Converted " + intervalFile.getPath().getFileName()); + header.setDescription("Converted " + intervalFile.getPath().getFileName()); + writer.write(header.toString()); + writer.newLine(); + + // Process each chromosome in the assembly + for (String chr : assembly) { + log.debug("Processing chromosome " + chr); + // Write the contig header to the output file + writer.write("fixedStep chrom="+chr+" start=1 step=1 span=1"); + writer.newLine(); + + int start = 1; + while (start < assembly.getChrLength(chr)) { + int stop = start + DEFAULT_CHUNK_SIZE - 1; + int length = stop - start + 1; + int[] count = new int[length]; + + Iterator<? extends Interval> it = intervalFile.query(chr, start, stop); + while (it.hasNext()) { + Interval entry = it.next(); + int entryStop = entry.getStop(); + if (extend != null) { + if (entry.isWatson()) { + entryStop = entry.getStart() + extend; + } else { + entryStop = entry.getStart() - extend; + } + } + + for (int i = entry.getStart(); i <= entryStop; i++) { + count[i-start]++; + } + mapped++; + } + + // Write the count at each base pair to the output file + for (int i = 0; i < count.length; i++) { + writer.write(count[i]); + writer.newLine(); + } + + // Process the next chunk + start = stop + 1; + } + } + } + + log.info("Mapped "+mapped+" reads"); + } + + public static void main(String[] args) { + new BaseAlignCounts().instanceMain(args); + } +} \ No newline at end of file