Java SamReaderFactory类代码示例

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
/**
 * A method that creates a SamReader object that's passed on to the HTSJDK.
 * Each time someone tries to open a SamReader on an URL,
 * HTSJDK checks if there's a custom reader factory and if it's there, this method is called.
 *
 * @param url target file URL
 * @return A SamReader object on a specified file URL
 */
@Override
public SamReader open(URL url) {
    PerformanceMonitor.start();
    AmazonS3URI amazonURI = new AmazonS3URI(url.toString());
    S3Client client = new S3Client();
    S3InputStreamFactory streamFactory = new S3InputStreamFactory(client);

    //download index file if is possible, and then start download .bam file
    final Optional<SeekableStream> indexStream;
    try {
        IndexLoader loader = new IndexLoader(client);
        indexStream = loader.loadIndex(amazonURI);
    } catch (IOException e) {
        throw new RuntimeIOException(e.getMessage() + " failed to download index", e);
    }

    SeekableStream stream = new S3SeekableStream(amazonURI, client, streamFactory);
    SamReaderFactory factory = SamReaderFactory.makeDefault();
    SamInputResource inputResource = SamInputResource.of(stream);

    indexStream.ifPresent(inputResource::index);

    return factory.open(inputResource);
}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
public static void main(String args[]) throws IOException, ParseException {
      Options options = new Options();
      options.addOption("u","uniquehits",false,"only output hits with a single mapping");
      options.addOption("s","nosuboptimal",false,"do not include hits whose score is not equal to the best score for the read");
      options.addOption("p","pairedend",false,"output paired-end hits");
      options.addOption("j","junctions",false,"output junction mapping reads (reads with gaps)");
      CommandLineParser parser = new GnuParser();
      CommandLine cl = parser.parse( options, args, false );            
  	uniqueOnly = cl.hasOption("uniquehits");
  	filterSubOpt = cl.hasOption("nosuboptimal");
  	inclPairedEnd = cl.hasOption("pairedend");
  	inclJunction = cl.hasOption("junctions");
  	SamReaderFactory factory =
          SamReaderFactory.makeDefault()
              .enable(SamReaderFactory.Option.INCLUDE_SOURCE_IN_RECORDS, SamReaderFactory.Option.VALIDATE_CRC_CHECKSUMS)
              .validationStringency(ValidationStringency.SILENT);
SamReader reader = factory.open(SamInputResource.of(System.in));
      CloseableIterator<SAMRecord> iter = reader.iterator();
      while (iter.hasNext()) {
          SAMRecord record = iter.next();
          if (record.getReadUnmappedFlag()) {continue; }
          processRecord(record);
      }
      iter.close();
      reader.close();
  }
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
private static void sliceFromS3(@NotNull final CommandLine cmd) throws IOException, EmptyFileException {
    final OkHttpClient httpClient = SbpS3Client.create(Integer.parseInt(cmd.getOptionValue(MAX_CONCURRENT_REQUESTS)));
    final URL bamUrl = SbpS3UrlGenerator.generateUrl(cmd.getOptionValue(BUCKET), cmd.getOptionValue(INPUT));
    final URL indexUrl = SbpS3UrlGenerator.generateUrl(cmd.getOptionValue(BUCKET), cmd.getOptionValue(INDEX));
    final String outputPath = cmd.getOptionValue(OUTPUT);
    final String bedPath = cmd.getOptionValue(BED);
    final int maxBufferSize = readMaxBufferSize(cmd);
    final File indexFile = downloadIndex(indexUrl);
    final SamReader reader = SamReaderFactory.makeDefault().open(SamInputResource.of(bamUrl).index(indexFile));
    LOGGER.info("Generating query intervals from BED file: {}", bedPath);
    final QueryInterval[] intervals = getIntervalsFromBED(bedPath, reader.getFileHeader());
    final BAMFileSpan span = bamSpanForIntervals(indexFile, reader.getFileHeader(), intervals);
    final List<Chunk> expandedChunks = expandChunks(span.getChunks());
    LOGGER.info("Generated {} query intervals which map to {} bam chunks", intervals.length, expandedChunks.size());
    final SamInputResource bamResource =
            SamInputResource.of(new CachingSeekableHTTPStream(httpClient, bamUrl, expandedChunks, maxBufferSize)).index(indexFile);
    final SamReader cachingReader = SamReaderFactory.makeDefault().open(bamResource);

    LOGGER.info("Slicing bam...");
    final CloseableIterator<SAMRecord> iterator = getIterator(cachingReader, intervals, span.toCoordinateArray());
    writeToSlice(outputPath, cachingReader.getFileHeader(), iterator);
    cachingReader.close();
    reader.close();
    indexFile.deleteOnExit();
}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
@NotNull
public Multimap<Chromosome, CobaltCount> fromBam(@NotNull final String referenceBam, @NotNull final String tumorBam)
        throws IOException, ExecutionException, InterruptedException {
    final File tumorFile = new File(tumorBam);
    final File referenceFile = new File(referenceBam);

    final SamReaderFactory readerFactory = SamReaderFactory.make();
    final String chromosomeLengthFileName = ChromosomeLengthFile.generateFilename(outputDirectory, tumor);
    final List<ChromosomeLength> lengths;
    try (SamReader reader = readerFactory.open(tumorFile)) {
        lengths = ChromosomeLengthFactory.create(reader.getFileHeader());
    }
    ChromosomeLengthFile.write(chromosomeLengthFileName, lengths);

    LOGGER.info("Calculating Read Count from {}", tumorFile.toString());
    final List<Future<ChromosomeReadCount>> tumorFutures = createFutures(readerFactory, tumorFile, lengths);

    LOGGER.info("Calculating Read Count from {}", referenceFile.toString());
    final List<Future<ChromosomeReadCount>> referenceFutures = createFutures(readerFactory, referenceFile, lengths);

    final Multimap<String, ReadCount> tumorCounts = fromFutures(tumorFutures);
    final Multimap<String, ReadCount> referenceCounts = fromFutures(referenceFutures);

    LOGGER.info("Read Count Complete");
    return CobaltCountFactory.merge(referenceCounts, tumorCounts);
}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
public static void main(String[] args) throws IOException {
	File bamFile = new File(
			"./data/miniCaviar_IDT_NEB.runA.NA12878.bwa.bam");
	File newBamFile = new File(
			"./data/miniCaviar_IDT_NEB.runA.NA12878.bwa.chrom1.bam");
	SAMFileWriter outputSam = null;
	final SamReader reader = SamReaderFactory.makeDefault()
			.validationStringency(ValidationStringency.SILENT)
			.open(bamFile);

	outputSam = new SAMFileWriterFactory().makeBAMWriter(
			reader.getFileHeader(), true, newBamFile);

	int currentReads = 0;
	for (final SAMRecord samRecord : reader) {
		if (samRecord.getReferenceIndex() == 0)
			writeBam(samRecord, outputSam);

	}

	reader.close();
	outputSam.close();
}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
public static void main(String[] args) {
	String bam = "./data/HG00115.chrom11.ILLUMINA.bwa.GBR.exome.20130415.bam";
	SamReader reader = SamReaderFactory.makeDefault()
			.validationStringency(ValidationStringency.SILENT)
			.open(new File(bam));
	int length = 0;
	for (final SAMRecord record : reader) {
		if (record.getReadUnmappedFlag())
			continue;
		length += record.getReadLength();
	}

	ReferenceSequenceFile rsf = ReferenceSequenceFileFactory
			.getReferenceSequenceFile(new File("data/hs37d5.fa"));
	ReferenceSequence rs = rsf.getSequence("11");

	System.out.println(length * 1.0 / rs.length());
}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
public static void main(String[] args) throws IOException {
	SamReader reader = SamReaderFactory
			.makeDefault()
			.validationStringency(ValidationStringency.SILENT)
			.open(new File(
					"./data/miniCaviar_IDT_NEB.runA.NA12878.bwa.chrom1.bam"));

	SAMFileHeader samFileHeader = reader.getFileHeader();

	System.out.println(samFileHeader.getTextHeader());

	int i = 0;
	for (SAMRecord record : reader) {
		if (i < 100)
			System.out.println(record.getAlignmentStart());
		i++;
	}
	reader.close();
}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
private boolean setSamSeqDictFromBam(String bamfile) {

		/*  ------------------------------------------------------ */
		/* This chunk prepares SamReader from local bam            */
		SamReaderFactory srf=SamReaderFactory.make();
		srf.validationStringency(ValidationStringency.SILENT);
		SamReader samReader;
		samReader= srf.open(new File(bamfile));
		/*  ------------------------------------------------------ */
		
		SAMSequenceDictionary seqDict = samReader.getFileHeader().getSequenceDictionary();
		if(seqDict != null && !seqDict.isEmpty()){
			this.setSamSeqDictSource(new File(bamfile).getAbsolutePath());
			this.setSamSeqDict(seqDict);
			return true;
		}
		return false;
	}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
public static long getAlignedReadCount(String bam) 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(bam)){
			samReader = SamReaderFactory.makeDefault().open(
					SamInputResource.of(new URL(bam)).index(new URL(bam + ".bai"))
			);
		} else {
			samReader= srf.open(new File(bam));
		}
		/*  ------------------------------------------------------ */

		List<SAMSequenceRecord> sequences = samReader.getFileHeader().getSequenceDictionary().getSequences();
		long alnCount= 0;
		for(SAMSequenceRecord x : sequences){
			alnCount += samReader.indexing().getIndex().getMetaData(x.getSequenceIndex()).getAlignedRecordCount();
		}
		samReader.close();
		return alnCount;
    }
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
public static boolean bamHasIndex(String bam) 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(bam)){
			samReader = SamReaderFactory.makeDefault().open(
					SamInputResource.of(new URL(bam)).index(new URL(bam + ".bai"))
			);
		} else {
			samReader= srf.open(new File(bam));
		}
		/*  ------------------------------------------------------ */

		// SamReaderFactory srf=SamReaderFactory.make();
		// srf.validationStringency(ValidationStringency.SILENT);
		// SamReader samReader = srf.open(new File(bam));
		boolean hasIndex= samReader.hasIndex();
		samReader.close();
		return hasIndex;
		
	}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
/**
 * Reads a map of sample name/bam filename from fileIn. For each bam file reads its indices, using samtools and gets a number of
 * aligned records from index information. Number of aligned records for each file is recorded into output file, Format of
 * output file: one line per bam file each line: "sample name" \t "number of aligned records"
 *
 * @param fileIn
 *            contains list of sample names and bam files to process, each line in file represent one file
 * @throws IOException
 *
 */

//perl version: 12-21 str
public static Map<String, Long> printStatsToFile(String fileIn) throws IOException {
    Map<String, String> files = parseFile(fileIn);
    Map<String, Long> result = new LinkedHashMap<>();

    for (Map.Entry<String, String> entry : files.entrySet()) {
        try (SamReader sam = SamReaderFactory.makeDefault().open(new File(entry.getValue()))) {

            SamReader.Indexing ind = sam.indexing();
            AbstractBAMFileIndex index = (AbstractBAMFileIndex)ind.getIndex();
            long count = 0;
            for (int i = 0; i < index.getNumberOfReferences(); i++) {
                BAMIndexMetaData meta = index.getMetaData(i);
                count += meta.getAlignedRecordCount();
            }
            result.put(entry.getKey(), count);
        }
    }
    return result;
}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
@BeforeClass
public void initHeaders() throws IOException {
    try (final SamReader normalBamReader = SamReaderFactory.makeDefault().open(NORMAL_BAM_FILE);
         final SamReader tumorBamReader = SamReaderFactory.makeDefault().open(TUMOR_BAM_FILE)) {
        normalHeader = normalBamReader.getFileHeader();
        tumorHeader = tumorBamReader.getFileHeader();

        normalHetPulldownExpected = new Pulldown(normalHeader);
        normalHetPulldownExpected.add(new AllelicCount(new SimpleInterval("1", 11522, 11522), 7, 4));
        normalHetPulldownExpected.add(new AllelicCount(new SimpleInterval("1", 12098, 12098), 8, 6));
        normalHetPulldownExpected.add(new AllelicCount(new SimpleInterval("1", 14630, 14630), 9, 8));
        normalHetPulldownExpected.add(new AllelicCount(new SimpleInterval("2", 14689, 14689), 6, 9));
        normalHetPulldownExpected.add(new AllelicCount(new SimpleInterval("2", 14982, 14982), 6, 5));

        tumorHetPulldownExpected = new Pulldown(tumorHeader);
        tumorHetPulldownExpected.add(new AllelicCount(new SimpleInterval("1", 11522, 11522), 7, 4));
        tumorHetPulldownExpected.add(new AllelicCount(new SimpleInterval("1", 12098, 12098), 8, 6));
        tumorHetPulldownExpected.add(new AllelicCount(new SimpleInterval("1", 14630, 14630), 9, 8));
        tumorHetPulldownExpected.add(new AllelicCount(new SimpleInterval("2", 14689, 14689), 6, 9));
        tumorHetPulldownExpected.add(new AllelicCount(new SimpleInterval("2", 14982, 14982), 6, 5));
    }
}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
@Test(dataProvider = "outputWriterProvider")
public void testWritingHeader(final File outputFile, final SAMProgramRecord record,
        final boolean addProgramGroup) throws Exception {
    Assert.assertFalse(outputFile.exists(),
            "broken test: test output file exists " + outputFile);
    final RTOutputBamArgumentCollection args = new RTOutputBamArgumentCollection();
    args.outputName = outputFile.getAbsolutePath();
    args.addOutputSAMProgramRecord = addProgramGroup;
    final GATKReadWriter writer =
            args.outputWriter(new SAMFileHeader(), (record == null) ? null : () -> record, true,
                    null
            );
    writer.close();
    Assert.assertTrue(outputFile.exists(), "not output written");
    final SAMFileHeader writtenHeader =
            SamReaderFactory.makeDefault().getFileHeader(outputFile);
    final SAMFileHeader expectedHeader = new SAMFileHeader();
    expectedHeader.setSortOrder(SAMFileHeader.SortOrder.unsorted);
    if (addProgramGroup && record != null) {
        expectedHeader.addProgramRecord(record);
    }
    Assert.assertEquals(writtenHeader, expectedHeader);
}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
/**
 * Simple implementation of a gather operations that uses SAMFileReaders and Writers in order to concatenate
 * multiple BAM files.
 */
private static void gatherNormally(final List<File> inputs, final File output, final boolean createIndex, final boolean createMd5,
                                   final File referenceFasta) {
    final SAMFileHeader header;
    {
        header = SamReaderFactory.makeDefault().referenceSequence(referenceFasta).getFileHeader(inputs.get(0));
    }

    final SAMFileWriter out = new SAMFileWriterFactory().setCreateIndex(createIndex).setCreateMd5File(createMd5).makeSAMOrBAMWriter(header, true, output);

    for (final File f : inputs) {
        log.info("Gathering " + f.getAbsolutePath());
        final SamReader in = SamReaderFactory.makeDefault().referenceSequence(referenceFasta).open(f);
        for (final SAMRecord rec : in) out.addAlignment(rec);
        CloserUtil.close(in);
    }

    out.close();
}
 

import htsjdk.samtools.SamReaderFactory; //导入依赖的package包/类
protected int doWork() {
    IOUtil.assertFileIsReadable(INPUT);
    IOUtil.assertFileIsWritable(OUTPUT);
    final SamReader reader = SamReaderFactory.makeDefault().referenceSequence(REFERENCE_SEQUENCE).open(INPUT);
    final SAMFileWriter writer = new SAMFileWriterFactory().makeWriter(reader.getFileHeader(), true, OUTPUT, REFERENCE_SEQUENCE);

    if (CREATE_INDEX && writer.getFileHeader().getSortOrder() != SAMFileHeader.SortOrder.coordinate) {
        throw new PicardException("Can't CREATE_INDEX unless sort order is coordinate");
    }

    final ProgressLogger progress = new ProgressLogger(Log.getInstance(SamFormatConverter.class));
    for (final SAMRecord rec : reader) {
        writer.addAlignment(rec);
        progress.record(rec);
    }
    CloserUtil.close(reader);
    writer.close();
    return 0;
}