2
|
1 package edu.unc.genomics.nucleosomes;
|
|
2
|
|
3 import java.io.BufferedWriter;
|
|
4 import java.io.IOException;
|
|
5 import java.nio.charset.Charset;
|
|
6 import java.nio.file.Files;
|
|
7 import java.nio.file.Path;
|
|
8
|
|
9 import org.apache.log4j.Logger;
|
|
10
|
|
11 import com.beust.jcommander.Parameter;
|
|
12
|
|
13 import edu.unc.genomics.CommandLineTool;
|
|
14 import edu.unc.genomics.CommandLineToolException;
|
|
15 import edu.unc.genomics.io.WigFile;
|
|
16 import edu.unc.genomics.io.WigFileException;
|
|
17
|
|
18 public class Phasogram extends CommandLineTool {
|
|
19
|
|
20 private static final Logger log = Logger.getLogger(Phasogram.class);
|
|
21
|
|
22 @Parameter(names = {"-i", "--input"}, description = "Input wig file (read counts)", required = true)
|
|
23 public WigFile inputFile;
|
|
24 @Parameter(names = {"-m", "--max"}, description = "Maximum phase shift", required = true)
|
|
25 public int maxPhase;
|
|
26 @Parameter(names = {"-o", "--output"}, description = "Output file (histogram)", required = true)
|
|
27 public Path outputFile;
|
|
28
|
|
29 public void run() throws IOException {
|
|
30 long[] phaseCounts = new long[maxPhase+1];
|
|
31
|
|
32 // Process each chromosome in the assembly
|
|
33 for (String chr : inputFile.chromosomes()) {
|
|
34 log.debug("Processing chromosome " + chr);
|
|
35
|
|
36 int start = inputFile.getChrStart(chr);
|
|
37 while (start < inputFile.getChrStop(chr)) {
|
|
38 int stop = start + DEFAULT_CHUNK_SIZE - 1;
|
|
39
|
|
40 try {
|
|
41 float[] data = WigFile.flattenData(inputFile.query(chr, start, stop), start, stop);
|
|
42 for (int i = 0; i < data.length-maxPhase; i++) {
|
|
43 for (int j = 0; j <= maxPhase; j++) {
|
|
44 phaseCounts[j] += data[i];
|
|
45 }
|
|
46 }
|
|
47
|
|
48 } catch (WigFileException e) {
|
|
49 log.fatal("Error querying data from Wig file!");
|
|
50 e.printStackTrace();
|
|
51 throw new CommandLineToolException("Error querying data from Wig file!");
|
|
52 }
|
|
53
|
|
54 // Process the next chunk
|
|
55 start = stop - maxPhase + 1;
|
|
56 }
|
|
57 }
|
|
58
|
|
59 log.debug("Writing output to disk");
|
|
60 try (BufferedWriter writer = Files.newBufferedWriter(outputFile, Charset.defaultCharset())) {
|
|
61 for (int i = 0; i < phaseCounts.length; i++) {
|
|
62 writer.write(i+"\t"+phaseCounts[i]);
|
|
63 writer.newLine();
|
|
64 }
|
|
65 }
|
|
66 }
|
|
67
|
|
68 public static void main(String[] args) {
|
|
69 new Phasogram().instanceMain(args);
|
|
70 }
|
|
71 } |