Search in sources :

Example 11 with SamReaderFactory

use of htsjdk.samtools.SamReaderFactory in project ASCIIGenome by dariober.

the class Utils method countReadsInWindow.

/**
 * Count reads in interval using the given filters.
 * @param bam
 * @param gc
 * @param filters List of filters to apply
 * @return
 * @throws MalformedURLException
 */
public static long countReadsInWindow(String bam, GenomicCoords gc, List<SamRecordFilter> filters) throws MalformedURLException {
    /*  ------------------------------------------------------ */
    /* This chunk prepares SamReader from local bam or URL bam */
    UrlValidator urlValidator = new UrlValidator();
    SamReaderFactory srf = SamReaderFactory.make();
    srf.validationStringency(ValidationStringency.SILENT);
    SamReader samReader;
    if (urlValidator.isValid(bam)) {
        samReader = srf.open(SamInputResource.of(new URL(bam)).index(new URL(bam + ".bai")));
    } else {
        samReader = srf.open(new File(bam));
    }
    /*  ------------------------------------------------------ */
    long cnt = 0;
    Iterator<SAMRecord> sam = samReader.query(gc.getChrom(), gc.getFrom(), gc.getTo(), false);
    AggregateFilter aggregateFilter = new AggregateFilter(filters);
    while (sam.hasNext()) {
        SAMRecord rec = sam.next();
        if (!aggregateFilter.filterOut(rec)) {
            cnt++;
        }
    }
    try {
        samReader.close();
    } catch (IOException e) {
        e.printStackTrace();
    }
    return cnt;
}
Also used : SamReader(htsjdk.samtools.SamReader) AggregateFilter(htsjdk.samtools.filter.AggregateFilter) SamReaderFactory(htsjdk.samtools.SamReaderFactory) SAMRecord(htsjdk.samtools.SAMRecord) UrlValidator(org.apache.commons.validator.routines.UrlValidator) IOException(java.io.IOException) File(java.io.File) IndexedFastaSequenceFile(htsjdk.samtools.reference.IndexedFastaSequenceFile) URL(java.net.URL)

Example 12 with SamReaderFactory

use of htsjdk.samtools.SamReaderFactory in project ASCIIGenome by dariober.

the class Utils method sortAndIndexSamOrBam.

/**
 * Sort and index input sam or bam.
 * @throws IOException
 */
public static void sortAndIndexSamOrBam(String inSamOrBam, String sortedBam, boolean deleteOnExit) throws IOException {
    /*  ------------------------------------------------------ */
    /* This chunk prepares SamReader from local bam or URL bam */
    UrlValidator urlValidator = new UrlValidator();
    SamReaderFactory srf = SamReaderFactory.make();
    srf.validationStringency(ValidationStringency.SILENT);
    SamReader samReader;
    if (urlValidator.isValid(inSamOrBam)) {
        samReader = SamReaderFactory.makeDefault().open(SamInputResource.of(new URL(inSamOrBam)));
    } else {
        samReader = srf.open(new File(inSamOrBam));
    }
    /*  ------------------------------------------------------ */
    samReader.getFileHeader().setSortOrder(SortOrder.coordinate);
    File out = new File(sortedBam);
    if (deleteOnExit) {
        out.deleteOnExit();
        File idx = new File(out.getAbsolutePath().replaceAll("\\.bam$", "") + ".bai");
        idx.deleteOnExit();
    }
    SAMFileWriter outputSam = new SAMFileWriterFactory().setCreateIndex(true).makeSAMOrBAMWriter(samReader.getFileHeader(), false, out);
    for (final SAMRecord samRecord : samReader) {
        outputSam.addAlignment(samRecord);
    }
    samReader.close();
    outputSam.close();
}
Also used : SamReader(htsjdk.samtools.SamReader) SamReaderFactory(htsjdk.samtools.SamReaderFactory) SAMFileWriter(htsjdk.samtools.SAMFileWriter) SAMRecord(htsjdk.samtools.SAMRecord) UrlValidator(org.apache.commons.validator.routines.UrlValidator) SAMFileWriterFactory(htsjdk.samtools.SAMFileWriterFactory) File(java.io.File) IndexedFastaSequenceFile(htsjdk.samtools.reference.IndexedFastaSequenceFile) URL(java.net.URL)

Example 13 with SamReaderFactory

use of htsjdk.samtools.SamReaderFactory in project ASCIIGenome by dariober.

the class Utils method getSamReader.

/**
 * Prepare SamReader from local bam or URL bam
 */
public static SamReader getSamReader(String workFilename) throws MalformedURLException {
    UrlValidator urlValidator = new UrlValidator();
    SamReaderFactory srf = SamReaderFactory.make();
    srf.validationStringency(ValidationStringency.SILENT);
    SamReader samReader;
    if (urlValidator.isValid(workFilename)) {
        samReader = srf.open(SamInputResource.of(new URL(workFilename)).index(new URL(workFilename + ".bai")));
    } else {
        samReader = srf.open(new File(workFilename));
    }
    return samReader;
}
Also used : SamReader(htsjdk.samtools.SamReader) SamReaderFactory(htsjdk.samtools.SamReaderFactory) UrlValidator(org.apache.commons.validator.routines.UrlValidator) File(java.io.File) IndexedFastaSequenceFile(htsjdk.samtools.reference.IndexedFastaSequenceFile) URL(java.net.URL)

Example 14 with SamReaderFactory

use of htsjdk.samtools.SamReaderFactory in project ASCIIGenome by dariober.

the class FilterTest method STUBcanCallBS.

// @Test
public void STUBcanCallBS() {
    System.out.println("\u203E");
    int f_incl = 0;
    int F_excl = 0;
    String chrom = "chrY";
    int from = 1;
    int to = 1;
    SamReaderFactory srf = SamReaderFactory.make();
    SamReader samReader = srf.open(new File("test_data/mjb050_oxBS.bam"));
    SAMFileHeader fh = samReader.getFileHeader();
    IntervalList il = new IntervalList(fh);
    Interval interval = new Interval(chrom, from, to);
    il.add(interval);
    List<SamRecordFilter> filters = FlagToFilter.flagToFilterList(f_incl, F_excl);
    SamLocusIterator samLocIter = new SamLocusIterator(samReader, il, true);
    samLocIter.setSamFilters(filters);
    while (samLocIter.hasNext()) {
        LocusInfo locus = samLocIter.next();
        int M = 0;
        int U = 0;
        int mism = 0;
        for (RecordAndOffset recOff : locus.getRecordAndPositions()) {
            int pos = locus.getPosition();
            // Code to get ref sequence at pos
            // If ref sequence is C, count read bases if:
            char refbase = Character.toUpperCase('\0');
            char rb = Character.toUpperCase((char) recOff.getReadBase());
            boolean isTopStrand = (new ReadFromTopStrandFilter(true)).filterOut(recOff.getRecord());
            if (refbase == 'C') {
                if (isTopStrand) {
                    // -ve 2nd pair
                    if (rb == 'C') {
                        M++;
                    } else if (rb == 'T') {
                        U++;
                    } else {
                        mism++;
                    }
                }
            } else if (refbase == 'G') {
                if (!isTopStrand) {
                    if (rb == 'G') {
                        M++;
                    } else if (rb == 'A') {
                        U++;
                    } else {
                        mism++;
                    }
                // System.out.println(locus.getPosition() + " ");
                }
            } else {
            // Not a C on ref
            }
        }
    // if(locus.getRecordAndPositions().size() > 0){
    // System.out.println(locus.getPosition() + " " + locus.getRecordAndPositions().size() + " " +
    // locus.getRecordAndPositions().get(0).getRecord().getFlags());
    // }
    }
}
Also used : SamReaderFactory(htsjdk.samtools.SamReaderFactory) SamRecordFilter(htsjdk.samtools.filter.SamRecordFilter) SamReader(htsjdk.samtools.SamReader) SamLocusIterator(htsjdk.samtools.util.SamLocusIterator) IntervalList(htsjdk.samtools.util.IntervalList) LocusInfo(htsjdk.samtools.util.SamLocusIterator.LocusInfo) RecordAndOffset(htsjdk.samtools.util.SamLocusIterator.RecordAndOffset) SAMFileHeader(htsjdk.samtools.SAMFileHeader) File(java.io.File) Interval(htsjdk.samtools.util.Interval)

Example 15 with SamReaderFactory

use of htsjdk.samtools.SamReaderFactory in project ASCIIGenome by dariober.

the class FilterTest method canFilterFromIntFlag.

// @Test
public void canFilterFromIntFlag() {
    int f_incl = 131;
    int F_excl = 72;
    String chrom = "chrY";
    int from = 1;
    int to = 100;
    SamReaderFactory srf = SamReaderFactory.make();
    SamReader samReader = srf.open(new File("test_data/mjb050_oxBS.bam"));
    SAMFileHeader fh = samReader.getFileHeader();
    IntervalList il = new IntervalList(fh);
    Interval interval = new Interval(chrom, from, to);
    il.add(interval);
    List<SamRecordFilter> filters = FlagToFilter.flagToFilterList(f_incl, F_excl);
    SamLocusIterator samLocIter = new SamLocusIterator(samReader, il, true);
    samLocIter.setSamFilters(filters);
    while (samLocIter.hasNext()) {
        LocusInfo locus = samLocIter.next();
        if (locus.getRecordAndPositions().size() > 0) {
        // System.out.println(locus.getPosition() + " " + locus.getRecordAndPositions().size() + " " +
        // locus.getRecordAndPositions().get(0).getRecord().getFlags());
        }
    }
}
Also used : SamReader(htsjdk.samtools.SamReader) SamLocusIterator(htsjdk.samtools.util.SamLocusIterator) SamReaderFactory(htsjdk.samtools.SamReaderFactory) SamRecordFilter(htsjdk.samtools.filter.SamRecordFilter) IntervalList(htsjdk.samtools.util.IntervalList) LocusInfo(htsjdk.samtools.util.SamLocusIterator.LocusInfo) SAMFileHeader(htsjdk.samtools.SAMFileHeader) File(java.io.File) Interval(htsjdk.samtools.util.Interval)

Aggregations

SamReaderFactory (htsjdk.samtools.SamReaderFactory)57 SamReader (htsjdk.samtools.SamReader)51 File (java.io.File)43 SAMRecord (htsjdk.samtools.SAMRecord)27 IOException (java.io.IOException)26 SAMFileHeader (htsjdk.samtools.SAMFileHeader)18 SAMRecordIterator (htsjdk.samtools.SAMRecordIterator)17 SAMSequenceDictionary (htsjdk.samtools.SAMSequenceDictionary)17 IndexedFastaSequenceFile (htsjdk.samtools.reference.IndexedFastaSequenceFile)14 ArrayList (java.util.ArrayList)13 List (java.util.List)11 SAMFileWriterFactory (htsjdk.samtools.SAMFileWriterFactory)10 HashSet (java.util.HashSet)10 SAMReadGroupRecord (htsjdk.samtools.SAMReadGroupRecord)9 SAMFileWriter (htsjdk.samtools.SAMFileWriter)8 SAMSequenceRecord (htsjdk.samtools.SAMSequenceRecord)8 URL (java.net.URL)8 HashMap (java.util.HashMap)8 BufferedReader (java.io.BufferedReader)7 PrintWriter (java.io.PrintWriter)7