From 2648b41398c5405de080d0f3193efd3d1d5a9b22 Mon Sep 17 00:00:00 2001 From: Karthik Gururaj Date: Tue, 25 Feb 2014 21:44:20 -0800 Subject: [PATCH] Added vectorized PairHMM implementation by Mohammad and Mustafa into the Maven build of GATK. C++ code has PAPI calls for reading hardware counters Followed Khalid's suggestion for packing libVectorLoglessCaching into the jar file with Maven Native library part of git repo 1. Renamed directory structure from public/c++/VectorPairHMM to public/VectorPairHMM/src/main/c++ as per Khalid's suggestion 2. Use java.home in public/VectorPairHMM/pom.xml to pass environment variable JRE_HOME to the make process. This is needed because the Makefile needs to compile JNI code with the flag -I/../include (among others). Assuming that the Maven build process uses a JDK (and not just a JRE), the variable java.home points to the JRE inside maven. 3. Dropped all pretense at cross-platform compatibility. Removed Mac profile from pom.xml for VectorPairHMM Moved JNI_README 1. Added the catch UnsatisfiedLinkError exception in PairHMMLikelihoodCalculationEngine.java to fall back to LOGLESS_CACHING in case the native library could not be loaded. Made VECTOR_LOGLESS_CACHING as the default implementation. 2. Updated the README with Mauricio's comments 3. baseline.cc is used within the library - if the machine supports neither AVX nor SSE4.1, the native library falls back to un-vectorized C++ in baseline.cc. 4. pairhmm-1-base.cc: This is not part of the library, but is being heavily used for debugging/profiling. Can I request that we keep it there for now? In the next release, we can delete it from the repository. 5. I agree with Mauricio about the ifdefs. I am sure you already know, but just to reassure you the debug code is not compiled into the library (because of the ifdefs) and will not affect performance. 1. Changed logger.info to logger.warn in PairHMMLikelihoodCalculationEngine.java 2. Committing the right set of files after rebase Added public license text to all C++ files Added license to Makefile Add package info to Sandbox.java --- .gitignore | 14 - ...GraphBasedLikelihoodCalculationEngine.java | 4 + .../haplotypecaller/HaplotypeCaller.java | 4 +- .../PairHMMLikelihoodCalculationEngine.java | 65 ++-- .../RandomLikelihoodCalculationEngine.java | 5 + .../utils/pairhmm/DebugJNILoglessPairHMM.java | 141 ++++---- .../utils/pairhmm/JNILoglessPairHMM.java | 0 .../utils/pairhmm/VectorLoglessPairHMM.java | 9 +- public/VectorPairHMM/README.md | 71 ++++ public/VectorPairHMM/pom.xml | 119 +++++++ .../src/main/c++}/.gitignore | 2 + .../src/main/c++}/LoadTimeInitializer.cc | 36 +- .../src/main/c++}/LoadTimeInitializer.h | 33 ++ .../src/main/c++}/Makefile | 37 +- .../src/main/c++}/Sandbox.cc | 29 +- .../src/main/c++}/Sandbox.h | 25 ++ .../src/main/c++}/Sandbox.java | 27 ++ .../Sandbox_JNIHaplotypeDataHolderClass.h | 0 .../c++}/Sandbox_JNIReadDataHolderClass.h | 0 .../main/c++/avx_function_instantiations.cc | 44 +++ .../src/main/c++}/baseline.cc | 56 ++- .../src/main/c++}/define-double.h | 25 ++ .../src/main/c++}/define-float.h | 25 ++ .../src/main/c++}/define-sse-double.h | 25 ++ .../src/main/c++}/define-sse-float.h | 25 ++ .../src/main/c++}/headers.h | 25 ++ .../VectorPairHMM/src/main/c++/jni_common.h | 58 ++++ .../src/main/c++}/jnidebug.h | 25 ++ ...ng_utils_pairhmm_DebugJNILoglessPairHMM.cc | 25 ++ ...ing_utils_pairhmm_DebugJNILoglessPairHMM.h | 25 ++ ...ting_utils_pairhmm_VectorLoglessPairHMM.cc | 83 ++++- ...sting_utils_pairhmm_VectorLoglessPairHMM.h | 25 ++ .../src/main/c++/pairhmm-1-base.cc | 70 ++++ .../src/main/c++}/pairhmm-template-kernel.cc | 28 +- .../src/main/c++}/pairhmm-template-main.cc | 25 ++ .../src/main/c++}/run.sh | 5 +- .../src/main/c++}/shift_template.c | 25 ++ .../main/c++/sse_function_instantiations.cc | 43 +++ public/VectorPairHMM/src/main/c++/template.h | 320 ++++++++++++++++++ .../src/main/c++}/utils.cc | 123 +++++-- .../src/main/c++}/utils.h | 32 ++ .../VectorPairHMM/src/main/c++/vector_defs.h | 55 +++ .../main/c++}/vector_function_prototypes.h | 25 ++ public/c++/VectorPairHMM/JNI_README | 41 --- .../avx_function_instantiations.cc | 19 -- public/c++/VectorPairHMM/jni_common.h | 33 -- public/c++/VectorPairHMM/pairhmm-1-base.cc | 45 --- .../sse_function_instantiations.cc | 18 - public/c++/VectorPairHMM/template.h | 194 ----------- public/c++/VectorPairHMM/vector_defs.h | 30 -- public/sting-root/pom.xml | 6 +- .../utils/pairhmm/libVectorLoglessPairHMM.so | Bin 0 -> 443803 bytes 52 files changed, 1686 insertions(+), 538 deletions(-) rename protected/{java/src => gatk-protected/src/main/java}/org/broadinstitute/sting/utils/pairhmm/DebugJNILoglessPairHMM.java (85%) rename protected/{java/src => gatk-protected/src/main/java}/org/broadinstitute/sting/utils/pairhmm/JNILoglessPairHMM.java (100%) rename protected/{java/src => gatk-protected/src/main/java}/org/broadinstitute/sting/utils/pairhmm/VectorLoglessPairHMM.java (98%) create mode 100644 public/VectorPairHMM/README.md create mode 100644 public/VectorPairHMM/pom.xml rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/.gitignore (88%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/LoadTimeInitializer.cc (79%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/LoadTimeInitializer.h (58%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/Makefile (66%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/Sandbox.cc (68%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/Sandbox.h (60%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/Sandbox.java (89%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/Sandbox_JNIHaplotypeDataHolderClass.h (100%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/Sandbox_JNIReadDataHolderClass.h (100%) create mode 100644 public/VectorPairHMM/src/main/c++/avx_function_instantiations.cc rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/baseline.cc (65%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/define-double.h (82%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/define-float.h (82%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/define-sse-double.h (77%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/define-sse-float.h (77%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/headers.h (56%) create mode 100644 public/VectorPairHMM/src/main/c++/jni_common.h rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/jnidebug.h (84%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc (85%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.h (78%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc (84%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.h (78%) create mode 100644 public/VectorPairHMM/src/main/c++/pairhmm-1-base.cc rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/pairhmm-template-kernel.cc (92%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/pairhmm-template-main.cc (70%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/run.sh (82%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/shift_template.c (71%) create mode 100644 public/VectorPairHMM/src/main/c++/sse_function_instantiations.cc create mode 100644 public/VectorPairHMM/src/main/c++/template.h rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/utils.cc (70%) rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/utils.h (52%) create mode 100644 public/VectorPairHMM/src/main/c++/vector_defs.h rename public/{c++/VectorPairHMM => VectorPairHMM/src/main/c++}/vector_function_prototypes.h (71%) delete mode 100644 public/c++/VectorPairHMM/JNI_README delete mode 100644 public/c++/VectorPairHMM/avx_function_instantiations.cc delete mode 100644 public/c++/VectorPairHMM/jni_common.h delete mode 100644 public/c++/VectorPairHMM/pairhmm-1-base.cc delete mode 100644 public/c++/VectorPairHMM/sse_function_instantiations.cc delete mode 100644 public/c++/VectorPairHMM/template.h delete mode 100644 public/c++/VectorPairHMM/vector_defs.h create mode 100644 public/sting-utils/src/main/resources/org/broadinstitute/sting/utils/pairhmm/libVectorLoglessPairHMM.so diff --git a/.gitignore b/.gitignore index a2ff166e5..ac3b931eb 100644 --- a/.gitignore +++ b/.gitignore @@ -24,19 +24,5 @@ dump/ lib/ out/ /atlassian-ide-plugin.xml -kg_tmp/ -maven-ant-tasks-2.1.3.jar -null-sequenceGraph.25.0.0.raw_readthreading_graph.dot -null-sequenceGraph.25.0.1.cleaned_readthreading_graph.dot -null-sequenceGraph.25.0.1.initial_seqgraph.dot -null-sequenceGraph.3.0.0.raw_readthreading_graph.dot -null-sequenceGraph.3.0.1.cleaned_readthreading_graph.dot -null-sequenceGraph.3.0.1.initial_seqgraph.dot -org/ -package-list -resources/ -velocity.log -perf -verify maven-metadata-local.xml dependency-reduced-pom.xml diff --git a/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/GraphBasedLikelihoodCalculationEngine.java b/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/GraphBasedLikelihoodCalculationEngine.java index 8a35ccb05..8b37e265d 100644 --- a/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/GraphBasedLikelihoodCalculationEngine.java +++ b/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/GraphBasedLikelihoodCalculationEngine.java @@ -155,4 +155,8 @@ public class GraphBasedLikelihoodCalculationEngine implements LikelihoodCalculat } } } + + @Override + public void close() { + } } diff --git a/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/HaplotypeCaller.java b/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/HaplotypeCaller.java index 74164415e..33dfa54ce 100644 --- a/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/HaplotypeCaller.java +++ b/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/HaplotypeCaller.java @@ -433,7 +433,7 @@ public class HaplotypeCaller extends ActiveRegionWalker, In */ @Hidden @Argument(fullName = "pair_hmm_implementation", shortName = "pairHMM", doc = "The PairHMM implementation to use for genotype likelihood calculations", required = false) - public PairHMM.HMM_IMPLEMENTATION pairHMM = PairHMM.HMM_IMPLEMENTATION.LOGLESS_CACHING; + public PairHMM.HMM_IMPLEMENTATION pairHMM = PairHMM.HMM_IMPLEMENTATION.VECTOR_LOGLESS_CACHING; @Hidden @Argument(fullName="keepRG", shortName="keepRG", doc="Only use read from this read group when making calls (but use all reads to build the assembly)", required = false) @@ -1051,7 +1051,7 @@ public class HaplotypeCaller extends ActiveRegionWalker, In referenceConfidenceModel.close(); //TODO remove the need to call close here for debugging, the likelihood output stream should be managed //TODO (open & close) at the walker, not the engine. - likelihoodCalculationEngine.close(); + likelihoodCalculationEngine.close(); logger.info("Ran local assembly on " + result + " active regions"); } diff --git a/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/PairHMMLikelihoodCalculationEngine.java b/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/PairHMMLikelihoodCalculationEngine.java index 9de532dcf..66755542f 100644 --- a/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/PairHMMLikelihoodCalculationEngine.java +++ b/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/PairHMMLikelihoodCalculationEngine.java @@ -80,29 +80,37 @@ public class PairHMMLikelihoodCalculationEngine implements LikelihoodCalculation private final boolean noFpga; private final ThreadLocal pairHMMThreadLocal = new ThreadLocal() { - @Override - protected PairHMM initialValue() { - switch (hmmType) { - case EXACT: return new Log10PairHMM(true); - case ORIGINAL: return new Log10PairHMM(false); - case LOGLESS_CACHING: - if (noFpga || !CnyPairHMM.isAvailable()) - return new LoglessPairHMM(); - else - return new CnyPairHMM(); - case VECTOR_LOGLESS_CACHING: - return new VectorLoglessPairHMM(); - case DEBUG_VECTOR_LOGLESS_CACHING: - return new DebugJNILoglessPairHMM(PairHMM.HMM_IMPLEMENTATION.VECTOR_LOGLESS_CACHING); - case ARRAY_LOGLESS: - if (noFpga || !CnyPairHMM.isAvailable()) - return new ArrayLoglessPairHMM(); - else - return new CnyPairHMM(); - default: - throw new UserException.BadArgumentValue("pairHMM", "Specified pairHMM implementation is unrecognized or incompatible with the HaplotypeCaller. Acceptable options are ORIGINAL, EXACT, CACHING, LOGLESS_CACHING, and ARRAY_LOGLESS."); - } - } + @Override + protected PairHMM initialValue() { + switch (hmmType) { + case EXACT: return new Log10PairHMM(true); + case ORIGINAL: return new Log10PairHMM(false); + case LOGLESS_CACHING: + if (noFpga || !CnyPairHMM.isAvailable()) + return new LoglessPairHMM(); + else + return new CnyPairHMM(); + case VECTOR_LOGLESS_CACHING: + try + { + return new VectorLoglessPairHMM(); + } + catch(UnsatisfiedLinkError ule) + { + logger.warn("Failed to load native library for VectorLoglessPairHMM - using Java implementation of LOGLESS_CACHING"); + return new LoglessPairHMM(); + } + case DEBUG_VECTOR_LOGLESS_CACHING: + return new DebugJNILoglessPairHMM(PairHMM.HMM_IMPLEMENTATION.VECTOR_LOGLESS_CACHING); + case ARRAY_LOGLESS: + if (noFpga || !CnyPairHMM.isAvailable()) + return new ArrayLoglessPairHMM(); + else + return new CnyPairHMM(); + default: + throw new UserException.BadArgumentValue("pairHMM", "Specified pairHMM implementation is unrecognized or incompatible with the HaplotypeCaller. Acceptable options are ORIGINAL, EXACT, CACHING, LOGLESS_CACHING, and ARRAY_LOGLESS."); + } + } }; // Attempted to do as below, to avoid calling pairHMMThreadLocal.get() later on, but it resulted in a NullPointerException // private final PairHMM pairHMM = pairHMMThreadLocal.get(); @@ -169,9 +177,10 @@ public class PairHMMLikelihoodCalculationEngine implements LikelihoodCalculation @Override public void close() { if ( likelihoodsStream != null ) likelihoodsStream.close(); - pairHMMThreadLocal.get().close(); + pairHMMThreadLocal.get().close(); } + private void writeDebugLikelihoods(final GATKSAMRecord processedRead, final Haplotype haplotype, final double log10l){ if ( WRITE_LIKELIHOODS_TO_FILE ) { likelihoodsStream.printf("%s %s %s %s %s %s %f%n", @@ -338,7 +347,7 @@ public class PairHMMLikelihoodCalculationEngine implements LikelihoodCalculation private void finalizePairHMM() { - pairHMMThreadLocal.get().finalizeRegion(); + pairHMMThreadLocal.get().finalizeRegion(); } @@ -358,9 +367,9 @@ public class PairHMMLikelihoodCalculationEngine implements LikelihoodCalculation map.filterPoorlyModelledReads(EXPECTED_ERROR_RATE_PER_BASE); stratifiedReadMap.put(sampleEntry.getKey(), map); } - - //Used mostly by the JNI implementation(s) to free arrays - finalizePairHMM(); + + //Used mostly by the JNI implementation(s) to free arrays + finalizePairHMM(); return stratifiedReadMap; } diff --git a/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/RandomLikelihoodCalculationEngine.java b/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/RandomLikelihoodCalculationEngine.java index b8dba7b86..d5d424ca9 100644 --- a/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/RandomLikelihoodCalculationEngine.java +++ b/protected/gatk-protected/src/main/java/org/broadinstitute/sting/gatk/walkers/haplotypecaller/RandomLikelihoodCalculationEngine.java @@ -79,4 +79,9 @@ public class RandomLikelihoodCalculationEngine implements LikelihoodCalculationE return result; } + + @Override + public void close() { + } + } diff --git a/protected/java/src/org/broadinstitute/sting/utils/pairhmm/DebugJNILoglessPairHMM.java b/protected/gatk-protected/src/main/java/org/broadinstitute/sting/utils/pairhmm/DebugJNILoglessPairHMM.java similarity index 85% rename from protected/java/src/org/broadinstitute/sting/utils/pairhmm/DebugJNILoglessPairHMM.java rename to protected/gatk-protected/src/main/java/org/broadinstitute/sting/utils/pairhmm/DebugJNILoglessPairHMM.java index 7e50f0c30..ea93ebe4a 100644 --- a/protected/java/src/org/broadinstitute/sting/utils/pairhmm/DebugJNILoglessPairHMM.java +++ b/protected/gatk-protected/src/main/java/org/broadinstitute/sting/utils/pairhmm/DebugJNILoglessPairHMM.java @@ -55,7 +55,7 @@ import org.broadinstitute.sting.utils.haplotype.Haplotype; import org.broadinstitute.sting.utils.sam.GATKSAMRecord; import org.broadinstitute.variant.variantcontext.Allele; import org.broadinstitute.sting.utils.exceptions.UserException; - +import static org.broadinstitute.sting.utils.pairhmm.PairHMMModel.*; import java.util.List; import java.util.Map; @@ -75,8 +75,9 @@ import java.io.IOException; */ public class DebugJNILoglessPairHMM extends LoglessPairHMM { + private static final boolean dumpSandboxOnly = false; //simulates ifdef private static final boolean debug = false; //simulates ifdef - private static final boolean verify = debug || true; //simulates ifdef + private static final boolean verify = !dumpSandboxOnly && (debug || true); //simulates ifdef private static final boolean debug0_1 = false; //simulates ifdef private static final boolean debug1 = false; //simulates ifdef private static final boolean debug2 = false; @@ -134,9 +135,11 @@ public class DebugJNILoglessPairHMM extends LoglessPairHMM { public void initialize( final List haplotypes, final Map> perSampleReadList, final int readMaxLength, final int haplotypeMaxLength ) { if(verify) + { super.initialize(haplotypes, perSampleReadList, readMaxLength, haplotypeMaxLength); - jniPairHMM.initialize(haplotypes, perSampleReadList, readMaxLength, haplotypeMaxLength); - haplotypeToHaplotypeListIdxMap = jniPairHMM.getHaplotypeToHaplotypeListIdxMap(); + jniPairHMM.initialize(haplotypes, perSampleReadList, readMaxLength, haplotypeMaxLength); + haplotypeToHaplotypeListIdxMap = jniPairHMM.getHaplotypeToHaplotypeListIdxMap(); + } } /** @@ -145,7 +148,8 @@ public class DebugJNILoglessPairHMM extends LoglessPairHMM { @Override public void finalizeRegion() { - jniPairHMM.finalizeRegion(); + if(!dumpSandboxOnly) + jniPairHMM.finalizeRegion(); } /** @@ -204,46 +208,61 @@ public class DebugJNILoglessPairHMM extends LoglessPairHMM { ++idx; } } - jniPairHMM.computeLikelihoods(reads, alleleHaplotypeMap, GCPArrayMap); - double[] likelihoodArray = jniPairHMM.getLikelihoodArray(); - //to compare values - final PerReadAlleleLikelihoodMap likelihoodMap = super.computeLikelihoods(reads, alleleHaplotypeMap, GCPArrayMap); + double[] likelihoodArray = null; + PerReadAlleleLikelihoodMap likelihoodMap = null; if(verify) { - //re-order values in likelihoodArray - double[] tmpArray = new double[numHaplotypes]; - idx = 0; - int idxInsideHaplotypeList = 0; - int readIdx = 0; - for(GATKSAMRecord read : reads) + jniPairHMM.computeLikelihoods(reads, alleleHaplotypeMap, GCPArrayMap); + likelihoodArray = jniPairHMM.getLikelihoodArray(); + //to compare values + likelihoodMap = super.computeLikelihoods(reads, alleleHaplotypeMap, GCPArrayMap); + } + else + { + likelihoodMap = new PerReadAlleleLikelihoodMap(); + likelihoodArray = new double[numTestcases]; + for(int i=0;i currEntry : alleleHaplotypeMap.entrySet())//order is important - access in same order always + //re-order values in likelihoodArray + double[] tmpArray = new double[numHaplotypes]; + idx = 0; + int idxInsideHaplotypeList = 0; + int readIdx = 0; + for(GATKSAMRecord read : reads) { - idxInsideHaplotypeList = haplotypeToHaplotypeListIdxMap.get(currEntry.getValue()); - likelihoodArray[idx] = tmpArray[idxInsideHaplotypeList]; - ++idx; + for(int j=0;j currEntry : alleleHaplotypeMap.entrySet())//order is important - access in same order always + { + idxInsideHaplotypeList = haplotypeToHaplotypeListIdxMap.get(currEntry.getValue()); + likelihoodArray[idx] = tmpArray[idxInsideHaplotypeList]; + ++idx; + } + readIdx += numHaplotypes; } - readIdx += numHaplotypes; - } - //for floating point values, no exact equality - //check whether numbers are close in terms of abs_error or relative_error - //For very large values, relative_error is relevant - //For very small values, abs_error is relevant - boolean toDump = false; - for(int i=0;i 1e-5 && relative_error > 1e-5) + //for floating point values, no exact equality + //check whether numbers are close in terms of abs_error or relative_error + //For very large values, relative_error is relevant + //For very small values, abs_error is relevant + for(int i=0;i 1e-5 && relative_error > 1e-5) + { + toDump = true; + break; + } } } //if numbers are not close, then dump out the data that produced the inconsistency @@ -251,24 +270,40 @@ public class DebugJNILoglessPairHMM extends LoglessPairHMM { { idx = 0; System.out.println("Dump : Java numReads "+readListSize+" numHaplotypes "+numHaplotypes); + boolean firstLine = true; for(GATKSAMRecord read : reads) { byte [] overallGCP = GCPArrayMap.get(read); + byte[] tmpByteArray = new byte[read.getReadBases().length]; for (Map.Entry currEntry : alleleHaplotypeMap.entrySet()) //order is important - access in same order always { byte[] haplotypeBases = currEntry.getValue().getBases(); debugDump("debug_dump.txt",new String(haplotypeBases)+" ",true); debugDump("debug_dump.txt",new String(read.getReadBases())+" ",true); for(int k=0;k + + 4.0.0 + + + org.broadinstitute.sting + sting-root + 2.8-SNAPSHOT + ../../public/sting-root + + + VectorPairHMM + pom + Vectorized PairHMM native libraries + + Builds a GNU/Linux x86_64 library of VectorPairHMM using icc (Intel C++ compiler). During install, copies it into sting-utils. Neither tested nor expected to work on any other platform. + + + UTF-8 + ${sourceEncoding} + ${sourceEncoding} + ${project.basedir}/../.. + ${sting.basedir}/public/sting-utils + + ${sting-utils.basedir}/src/main/resources/org/broadinstitute/sting/utils/pairhmm + + + + + + + + org.apache.maven.plugins + maven-enforcer-plugin + + + + display-info + + validate + + + + + + + org.codehaus.mojo + exec-maven-plugin + + + + exec + + compile + + make + src/main/c++ + + ${java.home} + ${project.build.directory} + + + + + + + + + org.apache.maven.plugins + maven-install-plugin + + true + + + + + + org.apache.maven.plugins + maven-resources-plugin + + + default-install + + copy-resources + + install + + ${pairhmm.resources.directory} + + + ${project.build.directory} + + **/* + + + + + + + + + + + com.google.code.sortpom + maven-sortpom-plugin + + false + custom_1 + \n + ${sourceEncoding} + true + scope + 4 + false + + + + + diff --git a/public/c++/VectorPairHMM/.gitignore b/public/VectorPairHMM/src/main/c++/.gitignore similarity index 88% rename from public/c++/VectorPairHMM/.gitignore rename to public/VectorPairHMM/src/main/c++/.gitignore index b034f9461..d791ffd80 100644 --- a/public/c++/VectorPairHMM/.gitignore +++ b/public/VectorPairHMM/src/main/c++/.gitignore @@ -12,3 +12,5 @@ reformat subdir_checkout.sh avx/ sse/ +triplicate.sh + diff --git a/public/c++/VectorPairHMM/LoadTimeInitializer.cc b/public/VectorPairHMM/src/main/c++/LoadTimeInitializer.cc similarity index 79% rename from public/c++/VectorPairHMM/LoadTimeInitializer.cc rename to public/VectorPairHMM/src/main/c++/LoadTimeInitializer.cc index f9db149e6..0e3026f65 100644 --- a/public/c++/VectorPairHMM/LoadTimeInitializer.cc +++ b/public/VectorPairHMM/src/main/c++/LoadTimeInitializer.cc @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #include "LoadTimeInitializer.h" #include "utils.h" using namespace std; @@ -14,7 +39,6 @@ char* LoadTimeInitializerStatsNames[] = "dummy" }; - LoadTimeInitializer g_load_time_initializer; LoadTimeInitializer::LoadTimeInitializer() //will be called when library is loaded @@ -47,7 +71,17 @@ LoadTimeInitializer::LoadTimeInitializer() //will be called when library is loa m_filename_to_fptr.clear(); m_written_files_set.clear(); + //Common buffer - 8MB + unsigned size = 1024*1024; + m_buffer = new uint64_t[size]; + m_buffer_size = size*sizeof(uint64_t); + initialize_function_pointers(); + + //Initialize static members of class + Context::initializeStaticMembers(); + Context::initializeStaticMembers(); + cout.flush(); } diff --git a/public/c++/VectorPairHMM/LoadTimeInitializer.h b/public/VectorPairHMM/src/main/c++/LoadTimeInitializer.h similarity index 58% rename from public/c++/VectorPairHMM/LoadTimeInitializer.h rename to public/VectorPairHMM/src/main/c++/LoadTimeInitializer.h index 0eb63849a..e7d9af626 100644 --- a/public/c++/VectorPairHMM/LoadTimeInitializer.h +++ b/public/VectorPairHMM/src/main/c++/LoadTimeInitializer.h @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #ifndef LOAD_TIME_INITIALIZER_H #define LOAD_TIME_INITIALIZER_H #include "headers.h" @@ -22,6 +47,10 @@ class LoadTimeInitializer { public: LoadTimeInitializer(); //will be called when library is loaded + ~LoadTimeInitializer() + { + delete m_buffer; + } void print_profiling(); void debug_dump(std::string filename, std::string s, bool to_append, bool add_newline=true); void debug_close(); @@ -43,6 +72,8 @@ class LoadTimeInitializer uint64_t m_data_transfer_time; //bytes copied uint64_t m_bytes_copied; + unsigned get_buffer_size() { return m_buffer_size; } + char* get_buffer() { return (char*)m_buffer; } private: std::map m_filename_to_fptr; std::set m_written_files_set; @@ -52,6 +83,8 @@ class LoadTimeInitializer double m_sum_square_stats[TOTAL_NUMBER_STATS]; uint64_t m_min_stats[TOTAL_NUMBER_STATS]; uint64_t m_max_stats[TOTAL_NUMBER_STATS]; + unsigned m_buffer_size; + uint64_t* m_buffer; }; extern LoadTimeInitializer g_load_time_initializer; diff --git a/public/c++/VectorPairHMM/Makefile b/public/VectorPairHMM/src/main/c++/Makefile similarity index 66% rename from public/c++/VectorPairHMM/Makefile rename to public/VectorPairHMM/src/main/c++/Makefile index 269cecbd2..354bca0bb 100644 --- a/public/c++/VectorPairHMM/Makefile +++ b/public/VectorPairHMM/src/main/c++/Makefile @@ -1,11 +1,36 @@ +#Copyright (c) 2012 The Broad Institute + +#Permission is hereby granted, free of charge, to any person +#obtaining a copy of this software and associated documentation +#files (the "Software"), to deal in the Software without +#restriction, including without limitation the rights to use, +#copy, modify, merge, publish, distribute, sublicense, and/or sell +#copies of the Software, and to permit persons to whom the +#Software is furnished to do so, subject to the following +#conditions: + +#The above copyright notice and this permission notice shall be +#included in all copies or substantial portions of the Software. + +#THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +#EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +#OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +#NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +#HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +#WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +#FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +#THE USE OR OTHER DEALINGS IN THE SOFTWARE. +# + + #OMPCFLAGS=-fopenmp #OMPLFLAGS=-fopenmp #-openmp-link static #CFLAGS=-O2 -std=c++11 -W -Wall -march=corei7-avx -Wa,-q -pedantic $(OMPCFLAGS) -Wno-unknown-pragmas #CFLAGS=-O2 -W -Wall -march=corei7 -mfpmath=sse -msse4.2 -pedantic $(OMPCFLAGS) -Wno-unknown-pragmas -JAVA_ROOT=/opt/jdk1.7.0_25/ -JNI_COMPILATION_FLAGS=-D_REENTRANT -fPIC -I${JAVA_ROOT}/include -I${JAVA_ROOT}/include/linux +JRE_HOME?=/opt/jdk1.7.0_25/jre +JNI_COMPILATION_FLAGS=-D_REENTRANT -fPIC -I${JRE_HOME}/../include -I${JRE_HOME}/../include/linux COMMON_COMPILATION_FLAGS=$(JNI_COMPILATION_FLAGS) -O3 -W -Wall -pedantic $(OMPCFLAGS) -Wno-unknown-pragmas CC=icc @@ -58,7 +83,7 @@ $(AVX_OBJECTS): CXXFLAGS=$(COMMON_COMPILATION_FLAGS) -xAVX $(SSE_OBJECTS): CXXFLAGS=$(COMMON_COMPILATION_FLAGS) -xSSE4.2 OBJECTS=$(NO_VECTOR_OBJECTS) $(AVX_OBJECTS) $(SSE_OBJECTS) -all: $(BIN) Sandbox.class +all: $(BIN) Sandbox.class copied_lib -include $(addprefix $(DEPDIR)/,$(SOURCES:.cc=.d)) @@ -79,5 +104,11 @@ $(OBJECTS): %.o: %.cc Sandbox.class: Sandbox.java javac Sandbox.java +copied_lib: libVectorLoglessPairHMM.so +ifdef OUTPUT_DIR + mkdir -p $(OUTPUT_DIR) + rsync -a libVectorLoglessPairHMM.so $(OUTPUT_DIR)/ +endif + clean: rm -rf $(BIN) *.o $(DEPDIR) *.class diff --git a/public/c++/VectorPairHMM/Sandbox.cc b/public/VectorPairHMM/src/main/c++/Sandbox.cc similarity index 68% rename from public/c++/VectorPairHMM/Sandbox.cc rename to public/VectorPairHMM/src/main/c++/Sandbox.cc index 7c10e0620..985b19ae9 100644 --- a/public/c++/VectorPairHMM/Sandbox.cc +++ b/public/VectorPairHMM/src/main/c++/Sandbox.cc @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #include "Sandbox.h" #include "org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.h" #include "utils.h" @@ -73,7 +98,9 @@ JNIEXPORT void JNICALL Java_Sandbox_doEverythingNative (JNIEnv* env, jobject thisObject, jstring fileNameString) { const char* fileName = env->GetStringUTFChars(fileNameString, 0); - do_compute((char*)fileName); + char local_array[800]; + strncpy(local_array, fileName, 200); env->ReleaseStringUTFChars(fileNameString, fileName); + do_compute(local_array, true, 10000, false); } diff --git a/public/c++/VectorPairHMM/Sandbox.h b/public/VectorPairHMM/src/main/c++/Sandbox.h similarity index 60% rename from public/c++/VectorPairHMM/Sandbox.h rename to public/VectorPairHMM/src/main/c++/Sandbox.h index 4ac1ea24c..486a1c095 100644 --- a/public/c++/VectorPairHMM/Sandbox.h +++ b/public/VectorPairHMM/src/main/c++/Sandbox.h @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + /* DO NOT EDIT THIS FILE - it is machine generated */ #include /* Header for class Sandbox */ diff --git a/public/c++/VectorPairHMM/Sandbox.java b/public/VectorPairHMM/src/main/c++/Sandbox.java similarity index 89% rename from public/c++/VectorPairHMM/Sandbox.java rename to public/VectorPairHMM/src/main/c++/Sandbox.java index 81a57c0a0..d6b7c2eae 100644 --- a/public/c++/VectorPairHMM/Sandbox.java +++ b/public/VectorPairHMM/src/main/c++/Sandbox.java @@ -1,3 +1,30 @@ +/* +* Copyright (c) 2012 The Broad Institute +* +* Permission is hereby granted, free of charge, to any person +* obtaining a copy of this software and associated documentation +* files (the "Software"), to deal in the Software without +* restriction, including without limitation the rights to use, +* copy, modify, merge, publish, distribute, sublicense, and/or sell +* copies of the Software, and to permit persons to whom the +* Software is furnished to do so, subject to the following +* conditions: +* +* The above copyright notice and this permission notice shall be +* included in all copies or substantial portions of the Software. +* +* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +* THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + +package org.broadinstitute.sting.utils.vectorpairhmm; + import java.util.List; import java.util.LinkedList; import java.util.Map; diff --git a/public/c++/VectorPairHMM/Sandbox_JNIHaplotypeDataHolderClass.h b/public/VectorPairHMM/src/main/c++/Sandbox_JNIHaplotypeDataHolderClass.h similarity index 100% rename from public/c++/VectorPairHMM/Sandbox_JNIHaplotypeDataHolderClass.h rename to public/VectorPairHMM/src/main/c++/Sandbox_JNIHaplotypeDataHolderClass.h diff --git a/public/c++/VectorPairHMM/Sandbox_JNIReadDataHolderClass.h b/public/VectorPairHMM/src/main/c++/Sandbox_JNIReadDataHolderClass.h similarity index 100% rename from public/c++/VectorPairHMM/Sandbox_JNIReadDataHolderClass.h rename to public/VectorPairHMM/src/main/c++/Sandbox_JNIReadDataHolderClass.h diff --git a/public/VectorPairHMM/src/main/c++/avx_function_instantiations.cc b/public/VectorPairHMM/src/main/c++/avx_function_instantiations.cc new file mode 100644 index 000000000..6d90d5070 --- /dev/null +++ b/public/VectorPairHMM/src/main/c++/avx_function_instantiations.cc @@ -0,0 +1,44 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + +#include "template.h" + +#undef SIMD_ENGINE +#undef SIMD_ENGINE_SSE + +#define SIMD_ENGINE avx +#define SIMD_ENGINE_AVX + +#include "define-float.h" +#include "shift_template.c" +#include "pairhmm-template-kernel.cc" + +#include "define-double.h" +#include "shift_template.c" +#include "pairhmm-template-kernel.cc" + +template double compute_full_prob_avxd(testcase* tc, double* nextlog); +template float compute_full_prob_avxs(testcase* tc, float* nextlog); + diff --git a/public/c++/VectorPairHMM/baseline.cc b/public/VectorPairHMM/src/main/c++/baseline.cc similarity index 65% rename from public/c++/VectorPairHMM/baseline.cc rename to public/VectorPairHMM/src/main/c++/baseline.cc index 232df83f2..d6085e661 100644 --- a/public/c++/VectorPairHMM/baseline.cc +++ b/public/VectorPairHMM/src/main/c++/baseline.cc @@ -1,7 +1,34 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #include "headers.h" #include "template.h" #include "utils.h" +#include "LoadTimeInitializer.h" using namespace std; + template NUMBER compute_full_prob(testcase *tc, NUMBER *before_last_log) { @@ -10,12 +37,28 @@ NUMBER compute_full_prob(testcase *tc, NUMBER *before_last_log) int COLS = tc->haplen + 1; Context ctx; - + //#define USE_STACK_ALLOCATION 1 +#ifdef USE_STACK_ALLOCATION + NUMBER M[ROWS][COLS]; + NUMBER X[ROWS][COLS]; + NUMBER Y[ROWS][COLS]; + NUMBER p[ROWS][6]; +#else //allocate on heap in way that simulates a 2D array. Having a 2D array instead of //a straightforward array of pointers ensures that all data lies 'close' in memory, increasing //the chance of being stored together in the cache. Also, prefetchers can learn memory access //patterns for 2D arrays, not possible for array of pointers + //bool locally_allocated = false; + //NUMBER* common_buffer = 0; NUMBER* common_buffer = new NUMBER[3*ROWS*COLS + ROWS*6]; + //unsigned curr_size = sizeof(NUMBER)*(3*ROWS*COLS + ROWS*6); + //if(true) + //{ + //common_buffer = new NUMBER[3*ROWS*COLS + ROWS*6]; + //locally_allocated = true; + //} + //else + //common_buffer = (NUMBER*)(g_load_time_initializer.get_buffer()); //pointers to within the allocated buffer NUMBER** common_pointer_buffer = new NUMBER*[4*ROWS]; NUMBER* ptr = common_buffer; @@ -25,14 +68,11 @@ NUMBER compute_full_prob(testcase *tc, NUMBER *before_last_log) for(;i<4*ROWS;++i, ptr+=6) common_pointer_buffer[i] = ptr; - //NUMBER M[ROWS][COLS]; - //NUMBER X[ROWS][COLS]; - //NUMBER Y[ROWS][COLS]; - //NUMBER p[ROWS][6]; NUMBER** M = common_pointer_buffer; NUMBER** X = M + ROWS; NUMBER** Y = X + ROWS; NUMBER** p = Y + ROWS; +#endif p[0][MM] = ctx._(0.0); @@ -47,7 +87,8 @@ NUMBER compute_full_prob(testcase *tc, NUMBER *before_last_log) int _i = tc->i[r-1] & 127; int _d = tc->d[r-1] & 127; int _c = tc->c[r-1] & 127; - p[r][MM] = ctx._(1.0) - ctx.ph2pr[(_i + _d) & 127]; + //p[r][MM] = ctx._(1.0) - ctx.ph2pr[(_i + _d) & 127]; + SET_MATCH_TO_MATCH_PROB(p[r][MM], _i, _d); p[r][GapM] = ctx._(1.0) - ctx.ph2pr[_c]; p[r][MX] = ctx.ph2pr[_i]; p[r][XX] = ctx.ph2pr[_c]; @@ -111,8 +152,11 @@ NUMBER compute_full_prob(testcase *tc, NUMBER *before_last_log) if (before_last_log != NULL) *before_last_log = result; +#ifndef USE_STACK_ALLOCATION delete common_pointer_buffer; + //if(locally_allocated) delete common_buffer; +#endif return result; //return ctx.LOG10(result) - ctx.LOG10_INITIAL_CONSTANT; diff --git a/public/c++/VectorPairHMM/define-double.h b/public/VectorPairHMM/src/main/c++/define-double.h similarity index 82% rename from public/c++/VectorPairHMM/define-double.h rename to public/VectorPairHMM/src/main/c++/define-double.h index 83589a13d..2067d369c 100644 --- a/public/c++/VectorPairHMM/define-double.h +++ b/public/VectorPairHMM/src/main/c++/define-double.h @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #include #ifdef PRECISION diff --git a/public/c++/VectorPairHMM/define-float.h b/public/VectorPairHMM/src/main/c++/define-float.h similarity index 82% rename from public/c++/VectorPairHMM/define-float.h rename to public/VectorPairHMM/src/main/c++/define-float.h index 87b2b01f3..318f78280 100644 --- a/public/c++/VectorPairHMM/define-float.h +++ b/public/VectorPairHMM/src/main/c++/define-float.h @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #include #ifdef PRECISION diff --git a/public/c++/VectorPairHMM/define-sse-double.h b/public/VectorPairHMM/src/main/c++/define-sse-double.h similarity index 77% rename from public/c++/VectorPairHMM/define-sse-double.h rename to public/VectorPairHMM/src/main/c++/define-sse-double.h index d781d55f3..2d271a854 100644 --- a/public/c++/VectorPairHMM/define-sse-double.h +++ b/public/VectorPairHMM/src/main/c++/define-sse-double.h @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #ifdef PRECISION #undef PRECISION #undef MAIN_TYPE diff --git a/public/c++/VectorPairHMM/define-sse-float.h b/public/VectorPairHMM/src/main/c++/define-sse-float.h similarity index 77% rename from public/c++/VectorPairHMM/define-sse-float.h rename to public/VectorPairHMM/src/main/c++/define-sse-float.h index 7516e6dbf..20af947dd 100644 --- a/public/c++/VectorPairHMM/define-sse-float.h +++ b/public/VectorPairHMM/src/main/c++/define-sse-float.h @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #ifdef PRECISION #undef PRECISION #undef MAIN_TYPE diff --git a/public/c++/VectorPairHMM/headers.h b/public/VectorPairHMM/src/main/c++/headers.h similarity index 56% rename from public/c++/VectorPairHMM/headers.h rename to public/VectorPairHMM/src/main/c++/headers.h index 48bd4d836..4a0d89b57 100644 --- a/public/c++/VectorPairHMM/headers.h +++ b/public/VectorPairHMM/src/main/c++/headers.h @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #ifndef COMMON_HEADERS_H #define COMMON_HEADERS_H diff --git a/public/VectorPairHMM/src/main/c++/jni_common.h b/public/VectorPairHMM/src/main/c++/jni_common.h new file mode 100644 index 000000000..23c323246 --- /dev/null +++ b/public/VectorPairHMM/src/main/c++/jni_common.h @@ -0,0 +1,58 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + +#ifndef JNI_COMMON_H +#define JNI_COMMON_H + +#include +/*#define ENABLE_ASSERTIONS 1*/ +#define DO_PROFILING 1 +/*#define DEBUG 1*/ +/*#define DEBUG0_1 1*/ +/*#define DEBUG3 1*/ +/*#define DUMP_TO_SANDBOX 1*/ + + +#define DIRECT_ACCESS_TO_JAVA_HEAP_MEMORY 1 + +#ifdef DIRECT_ACCESS_TO_JAVA_HEAP_MEMORY +//Gets direct access to Java arrays +#define GET_BYTE_ARRAY_ELEMENTS env->GetPrimitiveArrayCritical +#define RELEASE_BYTE_ARRAY_ELEMENTS env->ReleasePrimitiveArrayCritical +#define JNI_RO_RELEASE_MODE JNI_ABORT +#define GET_DOUBLE_ARRAY_ELEMENTS env->GetPrimitiveArrayCritical +#define RELEASE_DOUBLE_ARRAY_ELEMENTS env->ReleasePrimitiveArrayCritical + +#else +//Likely makes copy of Java arrays to JNI C++ space +#define GET_BYTE_ARRAY_ELEMENTS env->GetByteArrayElements +#define RELEASE_BYTE_ARRAY_ELEMENTS env->ReleaseByteArrayElements +#define JNI_RO_RELEASE_MODE JNI_ABORT +#define GET_DOUBLE_ARRAY_ELEMENTS env->GetDoubleArrayElements +#define RELEASE_DOUBLE_ARRAY_ELEMENTS env->ReleaseDoubleArrayElements + +#endif //ifdef DIRECT_ACCESS_TO_JAVA_HEAP_MEMORY + +#endif //ifndef JNI_COMMON_H diff --git a/public/c++/VectorPairHMM/jnidebug.h b/public/VectorPairHMM/src/main/c++/jnidebug.h similarity index 84% rename from public/c++/VectorPairHMM/jnidebug.h rename to public/VectorPairHMM/src/main/c++/jnidebug.h index 7002b3637..7fcab2a51 100644 --- a/public/c++/VectorPairHMM/jnidebug.h +++ b/public/VectorPairHMM/src/main/c++/jnidebug.h @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #ifndef JNI_DEBUG_H #define JNI_DEBUG_H diff --git a/public/c++/VectorPairHMM/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc b/public/VectorPairHMM/src/main/c++/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc similarity index 85% rename from public/c++/VectorPairHMM/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc rename to public/VectorPairHMM/src/main/c++/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc index 98ebcde92..8a3f8b5bc 100644 --- a/public/c++/VectorPairHMM/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc +++ b/public/VectorPairHMM/src/main/c++/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #include "headers.h" #include "jni_common.h" #include "org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.h" diff --git a/public/c++/VectorPairHMM/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.h b/public/VectorPairHMM/src/main/c++/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.h similarity index 78% rename from public/c++/VectorPairHMM/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.h rename to public/VectorPairHMM/src/main/c++/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.h index 9e099dae0..18c2b0394 100644 --- a/public/c++/VectorPairHMM/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.h +++ b/public/VectorPairHMM/src/main/c++/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.h @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + /* DO NOT EDIT THIS FILE - it is machine generated */ #include /* Header for class org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM */ diff --git a/public/c++/VectorPairHMM/org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc b/public/VectorPairHMM/src/main/c++/org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc similarity index 84% rename from public/c++/VectorPairHMM/org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc rename to public/VectorPairHMM/src/main/c++/org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc index 0f5219dd4..0b54c8a81 100644 --- a/public/c++/VectorPairHMM/org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc +++ b/public/VectorPairHMM/src/main/c++/org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #include "headers.h" #include "jni_common.h" #include "org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.h" @@ -98,6 +123,8 @@ JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_VectorLogless } } +//Create a vector of testcases for computation - copy the references to bytearrays read/readQuals etc into the appropriate +//testcase struct inline JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM_jniInitializeTestcasesVector (JNIEnv* env, jint numReads, jint numHaplotypes, jobjectArray& readDataArray, vector > >& readBasesArrayVector, vector& tc_array) @@ -195,21 +222,26 @@ inline JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_Vector inline void compute_testcases(vector& tc_array, unsigned numTestCases, double* likelihoodDoubleArray, unsigned maxNumThreadsToUse) { -#pragma omp parallel for schedule (dynamic,10000) num_threads(maxNumThreadsToUse) - for(unsigned tc_idx=0;tc_idxGetArrayLength(likelihoodArray) == numTestCases); #endif +#ifdef DO_WARMUP //ignore - only for crazy profiling + vector >& haplotypeBasesArrayVector = g_haplotypeBasesArrayVector; + for(unsigned i=0;iGetArrayLength(haplotypeBasesArrayVector[i].first); + for(unsigned j=0;jGetArrayLength(readBasesArrayVector[i][j].first); + for(unsigned k=0;k /* Header for class org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM */ diff --git a/public/VectorPairHMM/src/main/c++/pairhmm-1-base.cc b/public/VectorPairHMM/src/main/c++/pairhmm-1-base.cc new file mode 100644 index 000000000..7ff219b88 --- /dev/null +++ b/public/VectorPairHMM/src/main/c++/pairhmm-1-base.cc @@ -0,0 +1,70 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + +//#define DEBUG 1 +//#define DEBUG0_1 1 +//#define DEBUG3 1 +#include "headers.h" +#include "utils.h" +#include "LoadTimeInitializer.h" +using namespace std; + +int main(int argc, char** argv) +{ +#define BATCH_SIZE 10000 + if(argc < 2) + { + cerr << "Needs path to input file as argument\n"; + exit(0); + } + bool use_old_read_testcase = false; + if(argc >= 3 && string(argv[2]) == "1") + use_old_read_testcase = true; + unsigned chunk_size = 10000; + bool do_check = true; + uint64_t mask = ~(0ull); + for(int i=3;i @@ -116,7 +141,8 @@ template void CONCAT(CONCAT(initializeVectors,SIMD_ENGINE), PRECIS int _d = tc->d[r-1] & 127; int _c = tc->c[r-1] & 127; - *(ptr_p_MM+r-1) = ctx._(1.0) - ctx.ph2pr[(_i + _d) & 127]; + //*(ptr_p_MM+r-1) = ctx._(1.0) - ctx.ph2pr[(_i + _d) & 127]; + SET_MATCH_TO_MATCH_PROB(*(ptr_p_MM+r-1), _i, _d); *(ptr_p_GAPM+r-1) = ctx._(1.0) - ctx.ph2pr[_c]; *(ptr_p_MX+r-1) = ctx.ph2pr[_i]; *(ptr_p_XX+r-1) = ctx.ph2pr[_c]; diff --git a/public/c++/VectorPairHMM/pairhmm-template-main.cc b/public/VectorPairHMM/src/main/c++/pairhmm-template-main.cc similarity index 70% rename from public/c++/VectorPairHMM/pairhmm-template-main.cc rename to public/VectorPairHMM/src/main/c++/pairhmm-template-main.cc index cbcb8556d..5d9a05830 100644 --- a/public/c++/VectorPairHMM/pairhmm-template-main.cc +++ b/public/VectorPairHMM/src/main/c++/pairhmm-template-main.cc @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #include "headers.h" #include "template.h" #include "vector_defs.h" diff --git a/public/c++/VectorPairHMM/run.sh b/public/VectorPairHMM/src/main/c++/run.sh similarity index 82% rename from public/c++/VectorPairHMM/run.sh rename to public/VectorPairHMM/src/main/c++/run.sh index a6930678b..e821497f1 100755 --- a/public/c++/VectorPairHMM/run.sh +++ b/public/VectorPairHMM/src/main/c++/run.sh @@ -7,8 +7,9 @@ then pair_hmm_implementation=$1; fi -#-Djava.library.path is needed if you are using JNI_LOGLESS_CACHING, else not needed -java -Djava.library.path=${GSA_ROOT_DIR}/public/c++/VectorPairHMM -jar $GSA_ROOT_DIR/dist/GenomeAnalysisTK.jar -T HaplotypeCaller \ +#-Djava.library.path is needed if you wish to override the default 'packed' library +#java -jar $GSA_ROOT_DIR/target/GenomeAnalysisTK.jar -T HaplotypeCaller \ +java -Djava.library.path=${GSA_ROOT_DIR}/public/VectorPairHMM/src/main/c++ -jar $GSA_ROOT_DIR/target/GenomeAnalysisTK.jar -T HaplotypeCaller \ --dbsnp /data/broad/samples/joint_variant_calling/dbSNP/00-All.vcf \ -R /opt/Genomics/ohsu/dnapipeline/humanrefgenome/human_g1k_v37.fasta \ -I /data/simulated/sim1M_pairs_final.bam \ diff --git a/public/c++/VectorPairHMM/shift_template.c b/public/VectorPairHMM/src/main/c++/shift_template.c similarity index 71% rename from public/c++/VectorPairHMM/shift_template.c rename to public/VectorPairHMM/src/main/c++/shift_template.c index 9750e1d8b..c591c6882 100644 --- a/public/c++/VectorPairHMM/shift_template.c +++ b/public/VectorPairHMM/src/main/c++/shift_template.c @@ -1,3 +1,28 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #ifdef PRECISION #ifdef SIMD_ENGINE_AVX diff --git a/public/VectorPairHMM/src/main/c++/sse_function_instantiations.cc b/public/VectorPairHMM/src/main/c++/sse_function_instantiations.cc new file mode 100644 index 000000000..550272494 --- /dev/null +++ b/public/VectorPairHMM/src/main/c++/sse_function_instantiations.cc @@ -0,0 +1,43 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + +#include "template.h" + +#undef SIMD_ENGINE +#undef SIMD_ENGINE_AVX + +#define SIMD_ENGINE sse +#define SIMD_ENGINE_SSE + +#include "define-sse-float.h" +#include "shift_template.c" +#include "pairhmm-template-kernel.cc" + +#include "define-sse-double.h" +#include "shift_template.c" +#include "pairhmm-template-kernel.cc" + +template double compute_full_prob_ssed(testcase* tc, double* nextlog); +template float compute_full_prob_sses(testcase* tc, float* nextlog); diff --git a/public/VectorPairHMM/src/main/c++/template.h b/public/VectorPairHMM/src/main/c++/template.h new file mode 100644 index 000000000..ce4dbfc86 --- /dev/null +++ b/public/VectorPairHMM/src/main/c++/template.h @@ -0,0 +1,320 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + +#ifndef TEMPLATES_H_ +#define TEMPLATES_H_ + +#include "headers.h" + +#define MM 0 +#define GapM 1 +#define MX 2 +#define XX 3 +#define MY 4 +#define YY 5 + +//#define MROWS 500 +//#define MCOLS 1000 + +#define CAT(X,Y) X####Y +#define CONCAT(X,Y) CAT(X,Y) + +#define ALIGNED __attribute__((aligned(32))) + +typedef union __attribute__((aligned(32))) { + ALIGNED __m256 ALIGNED d; + ALIGNED __m128i ALIGNED s[2]; + ALIGNED float ALIGNED f[8]; + ALIGNED __m256i ALIGNED i; +} ALIGNED mix_F ALIGNED; + +typedef union __attribute__((aligned(32))) { + ALIGNED __m128 ALIGNED d; + ALIGNED __m64 ALIGNED s[2]; + ALIGNED float ALIGNED f[4]; + ALIGNED __m128i ALIGNED i; +} ALIGNED mix_F128 ALIGNED; + +typedef union ALIGNED { + __m128i vec ; + __m128 vecf ; + uint32_t masks[4] ; +} MaskVec_F ; + +typedef union ALIGNED { + __m64 vec ; + __m64 vecf ; + uint32_t masks[2] ; +} MaskVec_F128 ; + +typedef union ALIGNED +{ + ALIGNED __m128i ALIGNED i; + ALIGNED __m128 ALIGNED f; +} ALIGNED IF_128f ALIGNED; + +typedef union ALIGNED +{ + ALIGNED int ALIGNED i; + ALIGNED float ALIGNED f; +} ALIGNED IF_32 ALIGNED; + +typedef union __attribute__((aligned(32))) { + ALIGNED __m256d ALIGNED d; + ALIGNED __m128i ALIGNED s[2]; + ALIGNED double ALIGNED f[4]; + ALIGNED __m256i ALIGNED i; +} ALIGNED mix_D ALIGNED; + +typedef union __attribute__((aligned(32))) { + ALIGNED __m128d ALIGNED d; + ALIGNED __m64 ALIGNED s[2]; + ALIGNED double ALIGNED f[2]; + ALIGNED __m128i ALIGNED i; +} ALIGNED mix_D128 ALIGNED; + +typedef union ALIGNED { + __m128i vec ; + __m128d vecf ; + uint64_t masks[2] ; +} MaskVec_D ; + +typedef union ALIGNED { + __m64 vec ; + __m64 vecf ; + uint64_t masks[1] ; +} MaskVec_D128 ; + +typedef union ALIGNED +{ + ALIGNED __m128i ALIGNED i; + ALIGNED __m128d ALIGNED f; +} ALIGNED IF_128d ALIGNED; + +typedef union ALIGNED +{ + ALIGNED int64_t ALIGNED i; + ALIGNED double ALIGNED f; +} ALIGNED IF_64 ALIGNED; + + +#define MAX_QUAL 254 +#define MAX_JACOBIAN_TOLERANCE 8.0 +#define JACOBIAN_LOG_TABLE_STEP 0.0001 +#define JACOBIAN_LOG_TABLE_INV_STEP (1.0 / JACOBIAN_LOG_TABLE_STEP) +#define MAXN 70000 +#define LOG10_CACHE_SIZE (4*MAXN) // we need to be able to go up to 2*(2N) when calculating some of the coefficients +#define JACOBIAN_LOG_TABLE_SIZE ((int) (MAX_JACOBIAN_TOLERANCE / JACOBIAN_LOG_TABLE_STEP) + 1) + +template +struct ContextBase +{ + public: + NUMBER ph2pr[128]; + NUMBER INITIAL_CONSTANT; + NUMBER LOG10_INITIAL_CONSTANT; + NUMBER RESULT_THRESHOLD; + + static bool staticMembersInitializedFlag; + static NUMBER jacobianLogTable[JACOBIAN_LOG_TABLE_SIZE]; + static NUMBER matchToMatchProb[((MAX_QUAL + 1) * (MAX_QUAL + 2)) >> 1]; + + static void initializeStaticMembers() + { + if(!staticMembersInitializedFlag) + { + //Order of calls important - Jacobian first, then MatchToMatch + initializeJacobianLogTable(); + initializeMatchToMatchProb(); + staticMembersInitializedFlag = true; + } + } + + static void deleteStaticMembers() + { + if(staticMembersInitializedFlag) + { + staticMembersInitializedFlag = false; + } + } + + //Called only once during library load - don't bother to optimize with single precision fp + static void initializeJacobianLogTable() + { + for (int k = 0; k < JACOBIAN_LOG_TABLE_SIZE; k++) { + jacobianLogTable[k] = (NUMBER)(log10(1.0 + pow(10.0, -((double) k) * JACOBIAN_LOG_TABLE_STEP))); + } + } + + //Called only once per library load - don't bother optimizing with single fp + static void initializeMatchToMatchProb() + { + double LN10 = log(10); + double INV_LN10 = 1.0/LN10; + for (int i = 0, offset = 0; i <= MAX_QUAL; offset += ++i) + for (int j = 0; j <= i; j++) { + double log10Sum = approximateLog10SumLog10(-0.1*i, -0.1*j); + double matchToMatchLog10 = + log1p(-std::min(1.0,pow(10,log10Sum))) * INV_LN10; + matchToMatchProb[offset + j] = (NUMBER)(pow(10,matchToMatchLog10)); + } + } + //Called during computation - use single precision where possible + static int fastRound(NUMBER d) { + return (d > ((NUMBER)0.0)) ? (int) (d + ((NUMBER)0.5)) : (int) (d - ((NUMBER)0.5)); + } + //Called during computation - use single precision where possible + static NUMBER approximateLog10SumLog10(NUMBER small, NUMBER big) { + // make sure small is really the smaller value + if (small > big) { + NUMBER t = big; + big = small; + small = t; + } + + if (isinf(small) == -1 || isinf(big) == -1) + return big; + + NUMBER diff = big - small; + if (diff >= ((NUMBER)MAX_JACOBIAN_TOLERANCE)) + return big; + + // OK, so |y-x| < tol: we use the following identity then: + // we need to compute log10(10^x + 10^y) + // By Jacobian logarithm identity, this is equal to + // max(x,y) + log10(1+10^-abs(x-y)) + // we compute the second term as a table lookup with integer quantization + // we have pre-stored correction for 0,0.1,0.2,... 10.0 + int ind = fastRound((NUMBER)(diff * ((NUMBER)JACOBIAN_LOG_TABLE_INV_STEP))); // hard rounding + return big + jacobianLogTable[ind]; + } +}; + +template +struct Context : public ContextBase +{}; + +template<> +struct Context : public ContextBase +{ + Context():ContextBase() + { + for (int x = 0; x < 128; x++) + ph2pr[x] = pow(10.0, -((double)x) / 10.0); + + INITIAL_CONSTANT = ldexp(1.0, 1020.0); + LOG10_INITIAL_CONSTANT = log10(INITIAL_CONSTANT); + RESULT_THRESHOLD = 0.0; + } + + double LOG10(double v){ return log10(v); } + inline double POW(double b, double e) { return pow(b,e); } + + static double _(double n){ return n; } + static double _(float n){ return ((double) n); } +}; + +template<> +struct Context : public ContextBase +{ + Context() : ContextBase() + { + for (int x = 0; x < 128; x++) + { + ph2pr[x] = powf(10.f, -((float)x) / 10.f); + } + + INITIAL_CONSTANT = ldexpf(1.f, 120.f); + LOG10_INITIAL_CONSTANT = log10f(INITIAL_CONSTANT); + RESULT_THRESHOLD = ldexpf(1.f, -110.f); + } + + float LOG10(float v){ return log10f(v); } + inline float POW(float b, float e) { return powf(b,e); } + + static float _(double n){ return ((float) n); } + static float _(float n){ return n; } +}; + +#define SET_MATCH_TO_MATCH_PROB(output, insQual, delQual) \ +{ \ + int minQual = delQual; \ + int maxQual = insQual; \ + if (insQual <= delQual) \ + { \ + minQual = insQual; \ + maxQual = delQual; \ + } \ + (output) = (MAX_QUAL < maxQual) ? \ + ((NUMBER)1.0) - ctx.POW(((NUMBER)10), ctx.approximateLog10SumLog10(((NUMBER)-0.1)*minQual, ((NUMBER)-0.1)*maxQual)) \ + : ctx.matchToMatchProb[((maxQual * (maxQual + 1)) >> 1) + minQual]; \ +} + +typedef struct +{ + int rslen, haplen; + /*int *q, *i, *d, *c;*/ + /*int q[MROWS], i[MROWS], d[MROWS], c[MROWS];*/ + char *q, *i, *d, *c; + char *hap, *rs; + int *ihap; + int *irs; +} testcase; + +int normalize(char c); +int read_testcase(testcase *tc, FILE* ifp=0); + + +#define MIN_ACCEPTED 1e-28f +#define NUM_DISTINCT_CHARS 5 +#define AMBIG_CHAR 4 + +class ConvertChar { + + static uint8_t conversionTable[255] ; + +public: + + static void init() { + assert (NUM_DISTINCT_CHARS == 5) ; + assert (AMBIG_CHAR == 4) ; + + conversionTable['A'] = 0 ; + conversionTable['C'] = 1 ; + conversionTable['T'] = 2 ; + conversionTable['G'] = 3 ; + conversionTable['N'] = 4 ; + } + + static inline uint8_t get(uint8_t input) { + return conversionTable[input] ; + } + +}; + + +#endif + + diff --git a/public/c++/VectorPairHMM/utils.cc b/public/VectorPairHMM/src/main/c++/utils.cc similarity index 70% rename from public/c++/VectorPairHMM/utils.cc rename to public/VectorPairHMM/src/main/c++/utils.cc index 76cee7327..9f83cffa2 100644 --- a/public/c++/VectorPairHMM/utils.cc +++ b/public/VectorPairHMM/src/main/c++/utils.cc @@ -1,14 +1,48 @@ +/*Copyright (c) 2012 The Broad Institute + +*Permission is hereby granted, free of charge, to any person +*obtaining a copy of this software and associated documentation +*files (the "Software"), to deal in the Software without +*restriction, including without limitation the rights to use, +*copy, modify, merge, publish, distribute, sublicense, and/or sell +*copies of the Software, and to permit persons to whom the +*Software is furnished to do so, subject to the following +*conditions: + +*The above copyright notice and this permission notice shall be +*included in all copies or substantial portions of the Software. + +*THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +*EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +*OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +*NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +*HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +*WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +*FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR +*THE USE OR OTHER DEALINGS IN THE SOFTWARE. +*/ + + #include "headers.h" #include "template.h" #include "utils.h" #include "vector_defs.h" #include "LoadTimeInitializer.h" +using namespace std; +//static members from ConvertChar uint8_t ConvertChar::conversionTable[255]; +//Global function pointers in utils.h float (*g_compute_full_prob_float)(testcase *tc, float* before_last_log) = 0; double (*g_compute_full_prob_double)(testcase *tc, double* before_last_log) = 0; +//Static members in ContextBase +bool ContextBase::staticMembersInitializedFlag = false; +double ContextBase::jacobianLogTable[JACOBIAN_LOG_TABLE_SIZE]; +double ContextBase::matchToMatchProb[((MAX_QUAL + 1) * (MAX_QUAL + 2)) >> 1]; +bool ContextBase::staticMembersInitializedFlag = false; +float ContextBase::jacobianLogTable[JACOBIAN_LOG_TABLE_SIZE]; +float ContextBase::matchToMatchProb[((MAX_QUAL + 1) * (MAX_QUAL + 2)) >> 1]; -using namespace std; bool is_avx_supported() { @@ -61,6 +95,7 @@ uint64_t get_machine_capabilities() void initialize_function_pointers(uint64_t mask) { //mask = 0ull; + //mask = (1 << SSE41_CUSTOM_IDX); if(is_avx_supported() && (mask & (1<< AVX_CUSTOM_IDX))) { cout << "Using AVX accelerated implementation of PairHMM\n"; @@ -286,6 +321,13 @@ double getCurrClk() { return (double)tv.tv_sec + (double)tv.tv_usec / 1000000.0; } +inline unsigned long long rdtsc(void) +{ + unsigned hi, lo; + __asm__ __volatile__ ("rdtsc" : "=a"(lo), "=d"(hi)); + return ( (unsigned long long)lo)|( ((unsigned long long)hi)<<32 ); +} + void get_time(struct timespec* store_struct) { clock_gettime(CLOCK_REALTIME, store_struct); @@ -298,9 +340,12 @@ uint64_t diff_time(struct timespec& prev_time) return (uint64_t)((curr_time.tv_sec-prev_time.tv_sec)*1000000000+(curr_time.tv_nsec-prev_time.tv_nsec)); } -//#define DUMP_COMPUTE_VALUES 1 -#define BATCH_SIZE 10000 -#define RUN_HYBRID + +#ifdef USE_PAPI +#include "papi.h" +#define NUM_PAPI_COUNTERS 4 +#endif + void do_compute(char* filename, bool use_old_read_testcase, unsigned chunk_size, bool do_check) { FILE* fptr = 0; @@ -321,7 +366,22 @@ void do_compute(char* filename, bool use_old_read_testcase, unsigned chunk_size, uint64_t vector_compute_time = 0; uint64_t baseline_compute_time = 0; uint64_t num_double_calls = 0; + unsigned num_testcases = 0; bool all_ok = do_check ? true : false; +#ifdef USE_PAPI + uint32_t all_mask = (0); + uint32_t no_usr_mask = (1 << 16); //bit 16 user mode, bit 17 kernel mode + uint32_t no_kernel_mask = (1 << 17); //bit 16 user mode, bit 17 kernel mode + PAPI_num_counters(); + int events[NUM_PAPI_COUNTERS] = { 0, 0, 0, 0 }; + char* eventnames[NUM_PAPI_COUNTERS]= { "cycles", "itlb_walk_cycles", "dtlb_load_walk_cycles", "dtlb_store_walk_cycles" }; + assert(PAPI_event_name_to_code("UNHALTED_REFERENCE_CYCLES:u=1:k=1",&(events[0])) == PAPI_OK); + assert(PAPI_event_name_to_code("ITLB_MISSES:WALK_DURATION", &(events[1])) == PAPI_OK); + assert(PAPI_event_name_to_code("DTLB_LOAD_MISSES:WALK_DURATION", &(events[2])) == PAPI_OK); + assert(PAPI_event_name_to_code("DTLB_STORE_MISSES:WALK_DURATION", &(events[3])) == PAPI_OK); + long long values[NUM_PAPI_COUNTERS] = { 0, 0, 0, 0 }; + long long accum_values[NUM_PAPI_COUNTERS] = { 0, 0, 0, 0 }; +#endif while(1) { int break_value = use_old_read_testcase ? read_testcase(&tc, fptr) : read_mod_testcase(ifptr,&tc,true); @@ -336,26 +396,42 @@ void do_compute(char* filename, bool use_old_read_testcase, unsigned chunk_size, results_vec.resize(tc_vector.size()); baseline_results_vec.resize(tc_vector.size()); struct timespec start_time; +#ifdef USE_PAPI + assert(PAPI_start_counters(events, NUM_PAPI_COUNTERS) == PAPI_OK); +#endif get_time(&start_time); #pragma omp parallel for schedule(dynamic,chunk_size) num_threads(12) - for(unsigned i=0;i/bin/compilervars.sh intel64 -See run.sh in this directory on how to invoke HaplotypeCaller with the native library. The -argument -Djava.library.path is needed if the native implementation is selected, else -unnecessary. diff --git a/public/c++/VectorPairHMM/avx_function_instantiations.cc b/public/c++/VectorPairHMM/avx_function_instantiations.cc deleted file mode 100644 index 8f0de827d..000000000 --- a/public/c++/VectorPairHMM/avx_function_instantiations.cc +++ /dev/null @@ -1,19 +0,0 @@ -#include "template.h" - -#undef SIMD_ENGINE -#undef SIMD_ENGINE_SSE - -#define SIMD_ENGINE avx -#define SIMD_ENGINE_AVX - -#include "define-float.h" -#include "shift_template.c" -#include "pairhmm-template-kernel.cc" - -#include "define-double.h" -#include "shift_template.c" -#include "pairhmm-template-kernel.cc" - -template double compute_full_prob_avxd(testcase* tc, double* nextlog); -template float compute_full_prob_avxs(testcase* tc, float* nextlog); - diff --git a/public/c++/VectorPairHMM/jni_common.h b/public/c++/VectorPairHMM/jni_common.h deleted file mode 100644 index 1cffea1cc..000000000 --- a/public/c++/VectorPairHMM/jni_common.h +++ /dev/null @@ -1,33 +0,0 @@ -#ifndef JNI_COMMON_H -#define JNI_COMMON_H - -#include -/*#define ENABLE_ASSERTIONS 1*/ -#define DO_PROFILING 1 -/*#define DEBUG 1*/ -/*#define DEBUG0_1 1*/ -/*#define DEBUG3 1*/ -/*#define DUMP_TO_SANDBOX 1*/ - - -#define DIRECT_ACCESS_TO_JAVA_HEAP_MEMORY 1 - -#ifdef DIRECT_ACCESS_TO_JAVA_HEAP_MEMORY -//Gets direct access to Java arrays -#define GET_BYTE_ARRAY_ELEMENTS env->GetPrimitiveArrayCritical -#define RELEASE_BYTE_ARRAY_ELEMENTS env->ReleasePrimitiveArrayCritical -#define JNI_RO_RELEASE_MODE JNI_ABORT -#define GET_DOUBLE_ARRAY_ELEMENTS env->GetPrimitiveArrayCritical -#define RELEASE_DOUBLE_ARRAY_ELEMENTS env->ReleasePrimitiveArrayCritical - -#else -//Likely makes copy of Java arrays to JNI C++ space -#define GET_BYTE_ARRAY_ELEMENTS env->GetByteArrayElements -#define RELEASE_BYTE_ARRAY_ELEMENTS env->ReleaseByteArrayElements -#define JNI_RO_RELEASE_MODE JNI_ABORT -#define GET_DOUBLE_ARRAY_ELEMENTS env->GetDoubleArrayElements -#define RELEASE_DOUBLE_ARRAY_ELEMENTS env->ReleaseDoubleArrayElements - -#endif //ifdef DIRECT_ACCESS_TO_JAVA_HEAP_MEMORY - -#endif //ifndef JNI_COMMON_H diff --git a/public/c++/VectorPairHMM/pairhmm-1-base.cc b/public/c++/VectorPairHMM/pairhmm-1-base.cc deleted file mode 100644 index 78dbe296d..000000000 --- a/public/c++/VectorPairHMM/pairhmm-1-base.cc +++ /dev/null @@ -1,45 +0,0 @@ -//#define DEBUG 1 -//#define DEBUG0_1 1 -//#define DEBUG3 1 -#include "headers.h" -#include "utils.h" -#include "LoadTimeInitializer.h" -using namespace std; - -int main(int argc, char** argv) -{ -#define BATCH_SIZE 10000 - if(argc < 2) - { - cerr << "Needs path to input file as argument\n"; - exit(0); - } - bool use_old_read_testcase = false; - if(argc >= 3 && string(argv[2]) == "1") - use_old_read_testcase = true; - unsigned chunk_size = 10000; - bool do_check = true; - uint64_t mask = ~(0ull); - for(int i=3;i(testcase* tc, double* nextlog); -template float compute_full_prob_sses(testcase* tc, float* nextlog); diff --git a/public/c++/VectorPairHMM/template.h b/public/c++/VectorPairHMM/template.h deleted file mode 100644 index f91c2300e..000000000 --- a/public/c++/VectorPairHMM/template.h +++ /dev/null @@ -1,194 +0,0 @@ -#ifndef TEMPLATES_H_ -#define TEMPLATES_H_ - -#include "headers.h" - -#define MM 0 -#define GapM 1 -#define MX 2 -#define XX 3 -#define MY 4 -#define YY 5 - -//#define MROWS 500 -//#define MCOLS 1000 - -#define CAT(X,Y) X####Y -#define CONCAT(X,Y) CAT(X,Y) - -#define ALIGNED __attribute__((aligned(32))) - -typedef union __attribute__((aligned(32))) { - ALIGNED __m256 ALIGNED d; - ALIGNED __m128i ALIGNED s[2]; - ALIGNED float ALIGNED f[8]; - ALIGNED __m256i ALIGNED i; -} ALIGNED mix_F ALIGNED; - -typedef union __attribute__((aligned(32))) { - ALIGNED __m128 ALIGNED d; - ALIGNED __m64 ALIGNED s[2]; - ALIGNED float ALIGNED f[4]; - ALIGNED __m128i ALIGNED i; -} ALIGNED mix_F128 ALIGNED; - -typedef union ALIGNED { - __m128i vec ; - __m128 vecf ; - uint32_t masks[4] ; -} MaskVec_F ; - -typedef union ALIGNED { - __m64 vec ; - __m64 vecf ; - uint32_t masks[2] ; -} MaskVec_F128 ; - -typedef union ALIGNED -{ - ALIGNED __m128i ALIGNED i; - ALIGNED __m128 ALIGNED f; -} ALIGNED IF_128f ALIGNED; - -typedef union ALIGNED -{ - ALIGNED int ALIGNED i; - ALIGNED float ALIGNED f; -} ALIGNED IF_32 ALIGNED; - -typedef union __attribute__((aligned(32))) { - ALIGNED __m256d ALIGNED d; - ALIGNED __m128i ALIGNED s[2]; - ALIGNED double ALIGNED f[4]; - ALIGNED __m256i ALIGNED i; -} ALIGNED mix_D ALIGNED; - -typedef union __attribute__((aligned(32))) { - ALIGNED __m128d ALIGNED d; - ALIGNED __m64 ALIGNED s[2]; - ALIGNED double ALIGNED f[2]; - ALIGNED __m128i ALIGNED i; -} ALIGNED mix_D128 ALIGNED; - -typedef union ALIGNED { - __m128i vec ; - __m128d vecf ; - uint64_t masks[2] ; -} MaskVec_D ; - -typedef union ALIGNED { - __m64 vec ; - __m64 vecf ; - uint64_t masks[1] ; -} MaskVec_D128 ; - -typedef union ALIGNED -{ - ALIGNED __m128i ALIGNED i; - ALIGNED __m128d ALIGNED f; -} ALIGNED IF_128d ALIGNED; - -typedef union ALIGNED -{ - ALIGNED int64_t ALIGNED i; - ALIGNED double ALIGNED f; -} ALIGNED IF_64 ALIGNED; - -template -struct Context{}; - -template<> -struct Context -{ - Context() - { - for (int x = 0; x < 128; x++) - ph2pr[x] = pow(10.0, -((double)x) / 10.0); - - INITIAL_CONSTANT = ldexp(1.0, 1020.0); - LOG10_INITIAL_CONSTANT = log10(INITIAL_CONSTANT); - RESULT_THRESHOLD = 0.0; - } - - double LOG10(double v){ return log10(v); } - - static double _(double n){ return n; } - static double _(float n){ return ((double) n); } - double ph2pr[128]; - double INITIAL_CONSTANT; - double LOG10_INITIAL_CONSTANT; - double RESULT_THRESHOLD; -}; - -template<> -struct Context -{ - Context() - { - for (int x = 0; x < 128; x++) - { - ph2pr[x] = powf(10.f, -((float)x) / 10.f); - } - - INITIAL_CONSTANT = ldexpf(1.f, 120.f); - LOG10_INITIAL_CONSTANT = log10f(INITIAL_CONSTANT); - RESULT_THRESHOLD = ldexpf(1.f, -110.f); - } - - float LOG10(float v){ return log10f(v); } - - static float _(double n){ return ((float) n); } - static float _(float n){ return n; } - float ph2pr[128]; - float INITIAL_CONSTANT; - float LOG10_INITIAL_CONSTANT; - float RESULT_THRESHOLD; -}; - - -typedef struct -{ - int rslen, haplen; - /*int *q, *i, *d, *c;*/ - /*int q[MROWS], i[MROWS], d[MROWS], c[MROWS];*/ - char *q, *i, *d, *c; - char *hap, *rs; - int *ihap; - int *irs; -} testcase; - -int normalize(char c); -int read_testcase(testcase *tc, FILE* ifp=0); - - -#define MIN_ACCEPTED 1e-28f -#define NUM_DISTINCT_CHARS 5 -#define AMBIG_CHAR 4 - -class ConvertChar { - - static uint8_t conversionTable[255] ; - -public: - - static void init() { - assert (NUM_DISTINCT_CHARS == 5) ; - assert (AMBIG_CHAR == 4) ; - - conversionTable['A'] = 0 ; - conversionTable['C'] = 1 ; - conversionTable['T'] = 2 ; - conversionTable['G'] = 3 ; - conversionTable['N'] = 4 ; - } - - static inline uint8_t get(uint8_t input) { - return conversionTable[input] ; - } - -}; - - -#endif - - diff --git a/public/c++/VectorPairHMM/vector_defs.h b/public/c++/VectorPairHMM/vector_defs.h deleted file mode 100644 index 80b48ae99..000000000 --- a/public/c++/VectorPairHMM/vector_defs.h +++ /dev/null @@ -1,30 +0,0 @@ -#undef SIMD_ENGINE -#undef SIMD_ENGINE_AVX -#undef SIMD_ENGINE_SSE - -#define SIMD_ENGINE avx -#define SIMD_ENGINE_AVX - -#include "define-float.h" -#include "vector_function_prototypes.h" - -#include "define-double.h" -#include "vector_function_prototypes.h" - -#undef SIMD_ENGINE -#undef SIMD_ENGINE_AVX - -#define SIMD_ENGINE sse -#define SIMD_ENGINE_SSE - - -#include "define-sse-float.h" -#include "vector_function_prototypes.h" - -#include "define-sse-double.h" -#include "vector_function_prototypes.h" - -#undef SIMD_ENGINE -#undef SIMD_ENGINE_AVX -#undef SIMD_ENGINE_SSE - diff --git a/public/sting-root/pom.xml b/public/sting-root/pom.xml index 171eb7620..549a99ae6 100644 --- a/public/sting-root/pom.xml +++ b/public/sting-root/pom.xml @@ -335,7 +335,11 @@ maven-assembly-plugin 2.4 - + + org.apache.maven.plugins + maven-enforcer-plugin + 1.3.1 +