Package htsjdk.samtools

Examples of htsjdk.samtools.SAMFileReader


        return 0;
    }

    private void standardReheader(final SAMFileHeader replacementHeader) {
        final SAMFileReader recordReader = new SAMFileReader(INPUT);
        if (replacementHeader.getSortOrder() != recordReader.getFileHeader().getSortOrder()) {
            throw new PicardException("Sort orders of INPUT (" + recordReader.getFileHeader().getSortOrder().name() +
            ") and HEADER (" + replacementHeader.getSortOrder().name() + ") do not agree.");
        }
        final SAMFileWriter writer = new SAMFileWriterFactory().makeSAMOrBAMWriter(replacementHeader, true, OUTPUT);

        final ProgressLogger progress = new ProgressLogger(Log.getInstance(ReplaceSamHeader.class));
        for (final SAMRecord rec : recordReader) {
            rec.setHeader(replacementHeader);
            writer.addAlignment(rec);
            progress.record(rec);
        }
        writer.close();
        recordReader.close();
    }
View Full Code Here


                                   final long stopAfter,
                                   final Collection<SinglePassSamProgram> programs) {

        // Setup the standard inputs
        IOUtil.assertFileIsReadable(input);
        final SAMFileReader in = new SAMFileReader(input);

        // Optionally load up the reference sequence and double check sequence dictionaries
        final ReferenceSequenceFileWalker walker;
        if (referenceSequence == null) {
            walker = null;
        }
        else {
            IOUtil.assertFileIsReadable(referenceSequence);
            walker = new ReferenceSequenceFileWalker(referenceSequence);

            if (!in.getFileHeader().getSequenceDictionary().isEmpty()) {
                SequenceUtil.assertSequenceDictionariesEqual(in.getFileHeader().getSequenceDictionary(),
                                                             walker.getSequenceDictionary());
            }
        }

        // Check on the sort order of the BAM file
        {
            final SortOrder sort = in.getFileHeader().getSortOrder();
            if (sort != SortOrder.coordinate) {
                if (assumeSorted) {
                    log.warn("File reports sort order '" + sort + "', assuming it's coordinate sorted anyway.");
                }
                else {
                    throw new PicardException("File " + input.getAbsolutePath() + " should be coordinate sorted but " +
                                              "the header says the sort order is " + sort + ". If you believe the file " +
                                              "to be coordinate sorted you may pass ASSUME_SORTED=true");
                }
            }
        }

        // Call the abstract setup method!
        boolean anyUseNoRefReads = false;
        for (final SinglePassSamProgram program : programs) {
            program.setup(in.getFileHeader(), input);
            anyUseNoRefReads = anyUseNoRefReads || program.usesNoRefReads();
        }


        final ProgressLogger progress = new ProgressLogger(log);

        for (final SAMRecord rec : in) {
            final ReferenceSequence ref;
            if (walker == null || rec.getReferenceIndex() == SAMRecord.NO_ALIGNMENT_REFERENCE_INDEX) {
                ref = null;
            }
            else {
                ref = walker.get(rec.getReferenceIndex());
            }

            for (final SinglePassSamProgram program : programs) {
                program.acceptRead(rec, ref);
            }

            progress.record(rec);

            // See if we need to terminate early?
            if (stopAfter > 0 && progress.getCount() >= stopAfter) {
                break;
            }

            // And see if we're into the unmapped reads at the end
            if (!anyUseNoRefReads && rec.getReferenceIndex() == SAMRecord.NO_ALIGNMENT_REFERENCE_INDEX) {
                break;
            }
        }

        in.close();

        for (final SinglePassSamProgram program : programs) {
            program.finish();
        }
    }
View Full Code Here

    @Override
    protected int doWork() {
        IOUtil.assertFileIsReadable(INPUT);
        IOUtil.assertFileIsWritable(METRICS);

        final SAMFileReader in = new SAMFileReader(INPUT);
        final SAMFileHeader.SortOrder order = in.getFileHeader().getSortOrder();
        SAMFileWriter out = null;
        if (OUTPUT != null) {
            IOUtil.assertFileIsWritable(OUTPUT);
            out = new SAMFileWriterFactory().makeSAMOrBAMWriter(in.getFileHeader(), true, OUTPUT);
        }

        final Histogram<Integer> histo = new Histogram<Integer>("clipped_bases", "read_count");

        // Combine any adapters and custom adapter pairs from the command line into an array for use in clipping
        final AdapterPair[] adapters;
        {
            final List<AdapterPair> tmp = new ArrayList<AdapterPair>();
            tmp.addAll(ADAPTERS);
            if (FIVE_PRIME_ADAPTER != null && THREE_PRIME_ADAPTER != null) {
                tmp.add(new CustomAdapterPair(FIVE_PRIME_ADAPTER, THREE_PRIME_ADAPTER));
            }
            adapters = tmp.toArray(new AdapterPair[tmp.size()]);
        }

        ////////////////////////////////////////////////////////////////////////
        // Main loop that consumes reads, clips them and writes them to the output
        ////////////////////////////////////////////////////////////////////////
        final ProgressLogger progress = new ProgressLogger(log, 1000000, "Read");
        final SAMRecordIterator iterator = in.iterator();

        final AdapterMarker adapterMarker = new AdapterMarker(ADAPTER_TRUNCATION_LENGTH, adapters).
                setMaxPairErrorRate(MAX_ERROR_RATE_PE).setMinPairMatchBases(MIN_MATCH_BASES_PE).
                setMaxSingleEndErrorRate(MAX_ERROR_RATE_SE).setMinSingleEndMatchBases(MIN_MATCH_BASES_SE).
                setNumAdaptersToKeep(NUM_ADAPTERS_TO_KEEP).
View Full Code Here

     * This is factored out of doWork only for unit testing.
     */
    int writeSamText(PrintStream printStream) {
        try {
            IOUtil.assertFileIsReadable(INPUT);
            final SAMFileReader in = new SAMFileReader(INPUT);
            final AsciiWriter writer = new AsciiWriter(printStream);
            final SAMFileHeader header = in.getFileHeader();
            if (header.getTextHeader() != null) {
                writer.write(header.getTextHeader());
            } else {
                // Headers that are too large are not retained as text, so need to regenerate text
                new SAMTextHeaderCodec().encode(writer, header, true);
View Full Code Here

        public SeparateEndAlignmentIterator(final List<File> read1Alignments, final List<File> read2Alignments) {
            final List<SAMFileHeader> headers = new ArrayList<SAMFileHeader>();
            final List<SAMFileReader> read1 = new ArrayList<SAMFileReader>(read1Alignments.size());
            final List<SAMFileReader> read2 = new ArrayList<SAMFileReader>(read2Alignments.size());
            for (final File f : read1Alignments) {
                final SAMFileReader r = new SAMFileReader(f);
                headers.add(r.getFileHeader());
                read1.add(r);
            }
            for (final File f : read2Alignments) {
                final SAMFileReader r = new SAMFileReader(f);
                headers.add(r.getFileHeader());
                read2.add(r);
            }

            final SamFileHeaderMerger headerMerger = new SamFileHeaderMerger(SAMFileHeader.SortOrder.coordinate, headers, false);
            read1Iterator = new PeekableIterator<SAMRecord>(
View Full Code Here

    final File SUMMARY_OUT = new File(METRICS_FILE_PREFIX + SUMMARY_FILE_EXTENSION);
    final File DETAILS_OUT = new File(METRICS_FILE_PREFIX + DETAIL_FILE_EXTENSION);
    final File PLOTS_OUT = new File(METRICS_FILE_PREFIX + PDF_FILE_EXTENSION);
    assertIoFiles(SUMMARY_OUT, DETAILS_OUT, PLOTS_OUT);

    final SAMFileReader samReader = new SAMFileReader(INPUT);
    if (!ASSUME_SORTED && samReader.getFileHeader().getSortOrder() != SAMFileHeader.SortOrder.coordinate) {
      throw new PicardException("The input file " + INPUT.getAbsolutePath() + " does not appear to be coordinate sorted");
    }

    final ReferenceSequenceFileWalker refWalker = new ReferenceSequenceFileWalker(REFERENCE);
    final ProgressLogger progressLogger = new ProgressLogger(log);

    final RrbsMetricsCollector metricsCollector = new RrbsMetricsCollector(METRIC_ACCUMULATION_LEVEL, samReader.getFileHeader().getReadGroups(),
        C_QUALITY_THRESHOLD, NEXT_BASE_QUALITY_THRESHOLD, MINIMUM_READ_LENGTH, MAX_MISMATCH_RATE);

    for (final SAMRecord samRecord : samReader) {
      progressLogger.record(samRecord);
      if (!samRecord.getReadUnmappedFlag() && !isSequenceFiltered(samRecord.getReferenceName())) {
View Full Code Here

                OUTPUT = new File(baseFileName + BAMIndex.BAMIndexSuffix);
            }
        }

        IOUtil.assertFileIsWritable(OUTPUT);
        final SAMFileReader bam;

        if (inputUrl != null) {
            // remote input
            bam = new SAMFileReader(inputUrl, null, false);
        } else {
            // input from a normal file
            IOUtil.assertFileIsReadable(inputFile);
            bam = new SAMFileReader(inputFile);
        }

        if (!bam.isBinary()) {
            throw new SAMException("Input file must be bam file, not sam file.");
        }

        if (!bam.getFileHeader().getSortOrder().equals(SAMFileHeader.SortOrder.coordinate)) {
            throw new SAMException("Input bam file must be sorted by coordinates");
        }

        BAMIndexer.createIndex(bam, OUTPUT);
View Full Code Here

    protected int doWork() {
        IOUtil.assertFileIsReadable(INPUT);
        IOUtil.assertFileIsWritable(OUTPUT);

        final SAMFileReader in = new SAMFileReader(INPUT);

        // create the read group we'll be using
        final SAMReadGroupRecord rg = new SAMReadGroupRecord(RGID);
        rg.setLibrary(RGLB);
        rg.setPlatform(RGPL);
        rg.setSample(RGSM);
        rg.setPlatformUnit(RGPU);
        if (RGCN != null) rg.setSequencingCenter(RGCN);
        if (RGDS != null) rg.setDescription(RGDS);
        if (RGDT != null) rg.setRunDate(RGDT);
        if (RGPI != null) rg.setPredictedMedianInsertSize(RGPI);

        log.info(String.format("Created read group ID=%s PL=%s LB=%s SM=%s%n", rg.getId(), rg.getPlatform(), rg.getLibrary(), rg.getSample()));

        // create the new header and output file
        final SAMFileHeader inHeader = in.getFileHeader();
        final SAMFileHeader outHeader = inHeader.clone();
        outHeader.setReadGroups(Arrays.asList(rg));
        if (SORT_ORDER != null) outHeader.setSortOrder(SORT_ORDER);

        final SAMFileWriter outWriter = new SAMFileWriterFactory().makeSAMOrBAMWriter(outHeader,
                                                                                      outHeader.getSortOrder() == inHeader.getSortOrder(),
                                                                                      OUTPUT);

        final ProgressLogger progress = new ProgressLogger(log);
        for (final SAMRecord read : in) {
            read.setAttribute(SAMTag.RG.name(), RGID);
            outWriter.addAlignment(read);
            progress.record(read);
        }

        // cleanup
        in.close();
        outWriter.close();
        return 0;
    }
View Full Code Here

        int superSized = 0;
        int tandemPairs = 0;
        double chimeraSizeMinimum  = Math.max(getOutieMode(), (double)CHIMERA_KB_MIN);

        for (File f : INPUT) {
            SAMFileReader reader = new SAMFileReader(f);

            if (reader.getFileHeader().getSortOrder() != SAMFileHeader.SortOrder.coordinate) {
                throw new PicardException("SAM file must " + f.getName() + " must be sorted in coordintate order");
            }

            for (SAMRecord sam : reader) {

                // We're getting all our info from the first of each pair.
                if (!sam.getFirstOfPairFlag()) {
                    continue;
                }

                // Ignore unmapped read pairs
                if (sam.getReadUnmappedFlag()) {
                    if (!sam.getMateUnmappedFlag()) {
                        fragments++;
                        continue;
                    }

                    // If both ends are unmapped and we've hit unaligned reads we're done
                    if (sam.getReferenceIndex() == SAMRecord.NO_ALIGNMENT_REFERENCE_INDEX)  {
                        break;
                    }
                    continue;
                }

                if (sam.getMateUnmappedFlag()) {
                    fragments++;
                    continue;
                }

                // Ignore low-quality reads.  If we don't have the mate mapping quality, assume it's OK
                if ((sam.getAttribute(SAMTag.MQ.name()) != null &&
                     sam.getIntegerAttribute(SAMTag.MQ.name()) < MINIMUM_MAPPING_QUALITY) ||
                    sam.getMappingQuality() < MINIMUM_MAPPING_QUALITY) {
                    continue;
                }

                final int absInsertSize = Math.abs(sam.getInferredInsertSize());
                if (absInsertSize > chimeraSizeMinimum) {
                    superSized++;
                }
                else if (sam.getMateNegativeStrandFlag() == sam.getReadNegativeStrandFlag()) {
                    tandemPairs++;
                }
                else if (!sam.getMateReferenceIndex().equals(sam.getReferenceIndex())) {
                    crossChromPairs++;
                }
                else {
                    final PairOrientation pairOrientation = SamPairUtil.getPairOrientation(sam);
                    if (pairOrientation == PairOrientation.RF) {
                        outieHistogram.increment(absInsertSize);
                        outies++;
                        if (sam.getDuplicateReadFlag()) {
                            outieDupes++;
                        }
                    }
                    else if(pairOrientation == PairOrientation.FR) {
                        innieHistogram.increment(absInsertSize);
                        innies++;
                        if (sam.getDuplicateReadFlag()) {
                            innieDupes++;
                        }
                    }
                    else {
                        throw new IllegalStateException("This should never happen");
                    }
                }
            }
            reader.close();
        }

        MetricsFile<JumpingLibraryMetrics,Integer> metricsFile = getMetricsFile();
        JumpingLibraryMetrics metrics = new JumpingLibraryMetrics();
        metrics.JUMP_PAIRS = outies;
View Full Code Here

        int samplePerFile = SAMPLE_FOR_MODE/INPUT.size();

        Histogram<Integer> histo = new Histogram<Integer>();

        for (File f : INPUT) {
            SAMFileReader reader = new SAMFileReader(f);
            int sampled = 0;
            for (Iterator<SAMRecord> it = reader.iterator(); it.hasNext() && sampled < samplePerFile;) {
                SAMRecord sam = it.next();
                if (!sam.getFirstOfPairFlag()) {
                    continue;
                }
                // If we get here we've hit the end of the aligned reads
                if (sam.getReadUnmappedFlag() && sam.getReferenceIndex() == SAMRecord.NO_ALIGNMENT_REFERENCE_INDEX) {
                    break;
                }
                else if (sam.getReadUnmappedFlag() || sam.getMateUnmappedFlag()) {
                    continue;
                }
                else
                {
                    if ( (sam.getAttribute(SAMTag.MQ.name()) == null ||
                           sam.getIntegerAttribute(SAMTag.MQ.name()) >= MINIMUM_MAPPING_QUALITY) &&
                        sam.getMappingQuality() >= MINIMUM_MAPPING_QUALITY &&
                        sam.getMateNegativeStrandFlag() != sam.getReadNegativeStrandFlag() &&
                        sam.getMateReferenceIndex().equals(sam.getReferenceIndex())) {
                        if (SamPairUtil.getPairOrientation(sam) == PairOrientation.RF) {
                            histo.increment(Math.abs(sam.getInferredInsertSize()));
                            sampled++;
                        }
                    }
                }

            }
            reader.close();
        }

        return histo.size() > 0 ?  histo.getMode() : 0;
    }
View Full Code Here

TOP

Related Classes of htsjdk.samtools.SAMFileReader

Copyright © 2018 www.massapicom. All rights reserved.
All source code are property of their respective owners. Java is a trademark of Sun Microsystems, Inc and owned by ORACLE Inc. Contact coftware#gmail.com.