From 018e9e2c5f3ea7e1bb531c12be2371da02e6a440 Mon Sep 17 00:00:00 2001 From: Karthik Gururaj Date: Sun, 26 Jan 2014 19:18:12 -0800 Subject: [PATCH] 1. Cleaned up code 2. Split into DebugJNILoglessPairHMM and VectorLoglessPairHMM with base class JNILoglessPairHMM. DebugJNILoglessPairHMM can, in principle, invoke any other child class of JNILoglessPairHMM. 3. Added more profiling code for Java parts of LoglessPairHMM --- PairHMM_JNI/Makefile | 8 +- PairHMM_JNI/jni_common.h | 32 ++++ PairHMM_JNI/jnidebug.h | 166 ++++++++++++++++++ ...ng_utils_pairhmm_DebugJNILoglessPairHMM.cc | 151 ++++++++++++++++ ...te_sting_utils_pairhmm_JNILoglessPairHMM.h | 111 ------------ ...ing_utils_pairhmm_VectorLoglessPairHMM.cc} | 139 +++------------ PairHMM_JNI/run.sh | 10 +- .../PairHMMLikelihoodCalculationEngine.java | 4 +- .../utils/pairhmm/DebugJNILoglessPairHMM.java | 2 +- .../utils/pairhmm/JNILoglessPairHMM.java | 1 + .../utils/pairhmm/VectorLoglessPairHMM.java | 14 +- .../sting/utils/pairhmm/PairHMM.java | 39 ++-- 12 files changed, 423 insertions(+), 254 deletions(-) create mode 100644 PairHMM_JNI/jni_common.h create mode 100644 PairHMM_JNI/jnidebug.h create mode 100644 PairHMM_JNI/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc delete mode 100644 PairHMM_JNI/org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM.h rename PairHMM_JNI/{org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM.cc => org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc} (72%) diff --git a/PairHMM_JNI/Makefile b/PairHMM_JNI/Makefile index 9b461a833..e207dd4cc 100644 --- a/PairHMM_JNI/Makefile +++ b/PairHMM_JNI/Makefile @@ -13,7 +13,7 @@ CXX=icc LDFLAGS=-lm -lrt $(OMPLDFLAGS) -BIN=libJNILoglessPairHMM.so pairhmm-template-main checker +BIN=libVectorLoglessPairHMM.so pairhmm-template-main checker #BIN=checker DEPDIR=.deps @@ -22,14 +22,14 @@ DF=$(DEPDIR)/$(*).d #Common across libJNI and sandbox COMMON_SOURCES=utils.cc avx_function_instantiations.cc baseline.cc sse_function_instantiations.cc LoadTimeInitializer.cc #Part of libJNI -LIBSOURCES=org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM.cc $(COMMON_SOURCES) +LIBSOURCES=org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc $(COMMON_SOURCES) SOURCES=$(LIBSOURCES) pairhmm-template-main.cc pairhmm-1-base.cc LIBOBJECTS=$(LIBSOURCES:.cc=.o) COMMON_OBJECTS=$(COMMON_SOURCES:.cc=.o) #No vectorization for these files -NO_VECTOR_SOURCES=org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM.cc pairhmm-template-main.cc pairhmm-1-base.cc utils.cc baseline.cc LoadTimeInitializer.cc +NO_VECTOR_SOURCES=org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc pairhmm-template-main.cc pairhmm-1-base.cc utils.cc baseline.cc LoadTimeInitializer.cc #Use -xAVX for these files AVX_SOURCES=avx_function_instantiations.cc #Use -xSSE4.2 for these files @@ -53,7 +53,7 @@ checker: pairhmm-1-base.o $(COMMON_OBJECTS) pairhmm-template-main: pairhmm-template-main.o $(COMMON_OBJECTS) $(CXX) $(OMPLFLAGS) -o $@ $^ $(LDFLAGS) -libJNILoglessPairHMM.so: $(LIBOBJECTS) +libVectorLoglessPairHMM.so: $(LIBOBJECTS) $(CXX) $(OMPLFLAGS) -shared -o $@ $(LIBOBJECTS) ${LDFLAGS} diff --git a/PairHMM_JNI/jni_common.h b/PairHMM_JNI/jni_common.h new file mode 100644 index 000000000..d1ce86e9c --- /dev/null +++ b/PairHMM_JNI/jni_common.h @@ -0,0 +1,32 @@ +#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 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/PairHMM_JNI/jnidebug.h b/PairHMM_JNI/jnidebug.h new file mode 100644 index 000000000..7002b3637 --- /dev/null +++ b/PairHMM_JNI/jnidebug.h @@ -0,0 +1,166 @@ +#ifndef JNI_DEBUG_H +#define JNI_DEBUG_H + +template +class DataHolder +{ +#define INIT_MATRIX(X) \ + X = new NUMBER*[m_paddedMaxReadLength]; \ + for(int i=0;i ctx; + for (int r = 1; r <= length;r++) //in original code, r < ROWS (where ROWS = paddedReadLength) + { + int _i = insertionGOP[r-1]; //insertionGOP + int _d = deletionGOP[r-1]; //deletionGOP + int _c = overallGCP[r-1]; //overallGCP + m_transition[r][MM] = ctx._(1.0) - ctx.ph2pr[(_i + _d) & 127]; //lines 161-162 + m_transition[r][GapM] = ctx._(1.0) - ctx.ph2pr[_c]; //line 163 + m_transition[r][MX] = ctx.ph2pr[_i]; //164 + m_transition[r][XX] = ctx.ph2pr[_c]; //165 + m_transition[r][MY] = ctx.ph2pr[_d];//last row seems different, compared to line 166 + m_transition[r][YY] = ctx.ph2pr[_c];//same as above for line 167 + //m_transition[r][MY] = (r == length) ? ctx._(1.0) : ctx.ph2pr[_d];//last row seems different, compared to line 166 + //m_transition[r][YY] = (r == length) ? ctx._(1.0) : ctx.ph2pr[_c];//same as above for line 167 +#ifdef DEBUG3 + for(int j=0;j<6;++j) + debug_dump("transitions_jni.txt", to_string(m_transition[r][j]),true); +#endif + } + ++g_num_prob_init; + } + bool m_is_initialized; + int m_readMaxLength; + int m_haplotypeMaxLength; + int m_paddedMaxReadLength; + int m_paddedMaxHaplotypeLength; + NUMBER** m_matchMatrix; + NUMBER** m_insertionMatrix; + NUMBER** m_deletionMatrix; + NUMBER** m_prior; + NUMBER (*m_transition)[6]; +}; +extern DataHolder g_double_dataholder; + +template +NUMBER compute_full_prob(testcase *tc, NUMBER** M, NUMBER** X, NUMBER** Y, NUMBER (*p)[6], + bool do_initialization, jint hapStartIndex, NUMBER *before_last_log = NULL) +{ + int r, c; + int ROWS = tc->rslen + 1; //ROWS = paddedReadLength + int COLS = tc->haplen + 1; //COLS = paddedHaplotypeLength + + Context ctx; + //////NOTES + ////ctx.ph2pr[quality]; //This quantity is QualityUtils.qualToErrorProb(quality) + ////1-ctx.ph2pr[quality]; //This corresponds to QualityUtils.qualToProb(quality); + + //Initialization + if(do_initialization) + { + for (c = 0; c < COLS; c++) + { + M[0][c] = ctx._(0.0); + X[0][c] = ctx._(0.0); + Y[0][c] = ctx.INITIAL_CONSTANT / (tc->haplen); //code from 87-90 in LoglessPairHMM + } + + for (r = 1; r < ROWS; r++) + { + M[r][0] = ctx._(0.0); + //deletionMatrix row 0 in above nest is initialized in the Java code + //However, insertionMatrix column 0 is not initialized in Java code, could it be that + //values are re-used from a previous iteration? + //Why even do this, X[0][0] = 0 from above loop nest, X[idx][0] = 0 from this computation + X[r][0] = X[r-1][0] * p[r][XX]; + Y[r][0] = ctx._(0.0); + } + } + + for (r = 1; r < ROWS; r++) + for (c = hapStartIndex+1; c < COLS; c++) + { + //The following lines correspond to initializePriors() + char _rs = tc->rs[r-1]; //line 137 + char _hap = tc->hap[c-1]; //line 140 + //int _q = tc->q[r-1] & 127; //line 138 - q is the quality (qual), should be byte hence int ANDed with 0xFF + int _q = tc->q[r-1]; //line 138 - q is the quality (qual), should be byte hence int ANDed with 0xFF + NUMBER distm = ctx.ph2pr[_q]; //This quantity is QualityUtils.qualToErrorProb(_q) + //The assumption here is that doNotUseTristateCorrection is true + //TOASK + if (_rs == _hap || _rs == 'N' || _hap == 'N') + distm = ctx._(1.0) - distm; //This is the quantity QualityUtils.qualToProb(qual) + else + distm = distm/3; +#ifdef DEBUG3 + debug_dump("priors_jni.txt",to_string(distm),true); +#endif + + //Computation inside updateCell + M[r][c] = distm * (M[r-1][c-1] * p[r][MM] + X[r-1][c-1] * p[r][GapM] + Y[r-1][c-1] * p[r][GapM]); + X[r][c] = M[r-1][c] * p[r][MX] + X[r-1][c] * p[r][XX]; + Y[r][c] = M[r][c-1] * p[r][MY] + Y[r][c-1] * p[r][YY]; +#ifdef DEBUG3 + debug_dump("matrices_jni.txt",to_string(M[r][c]),true); + debug_dump("matrices_jni.txt",to_string(X[r][c]),true); + debug_dump("matrices_jni.txt",to_string(Y[r][c]),true); +#endif + } + + NUMBER result = ctx._(0.0); + for (c = 0; c < COLS; c++) + result += M[ROWS-1][c] + X[ROWS-1][c]; + + if (before_last_log != NULL) + *before_last_log = result; + +#ifdef DEBUG + debug_dump("return_values_jni.txt",to_string(ctx.LOG10(result) - ctx.LOG10_INITIAL_CONSTANT),true); +#endif + return ctx.LOG10(result) - ctx.LOG10_INITIAL_CONSTANT; +} + +#endif diff --git a/PairHMM_JNI/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc b/PairHMM_JNI/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc new file mode 100644 index 000000000..30f118c31 --- /dev/null +++ b/PairHMM_JNI/org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.cc @@ -0,0 +1,151 @@ +#include "headers.h" +#include "jni_common.h" +#include "org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM.h" +#include "template.h" +#include "utils.h" +#include "LoadTimeInitializer.h" +#include "jnidebug.h" +DataHolder g_double_dataholder; + +using namespace std; + +JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM_jniInitialize +(JNIEnv* env, jobject thisObject, + jint readMaxLength, jint haplotypeMaxLength) +{ + static int g_num_init_calls = 0; +#ifdef DEBUG3 + cout << "Entered alloc initialized .. readMaxLength "<GetArrayLength(insertionGOP); +#ifdef DEBUG3 + cout << "Entered initializeProbabilities .. length "<GetByteArrayElements(insertionGOP, &is_copy); + jbyte* deletionGOPArray = (env)->GetByteArrayElements(deletionGOP, &is_copy); + jbyte* overallGCPArray = (env)->GetByteArrayElements(overallGCP, &is_copy); +#ifdef DEBUG + if(insertionGOPArray == 0) + cerr << "insertionGOP array not initialized in JNI\n"; + ////assert(insertionGOPArray && "insertionGOP array not initialized in JNI"); + if(deletionGOPArray == 0) + cerr << "deletionGOP array not initialized in JNI\n"; + ////assert(deletionGOPArray && "deletionGOP array not initialized in JNI"); + assert(overallGCPArray && "OverallGCP array not initialized in JNI"); +#endif + + g_double_dataholder.initializeProbabilities(length, insertionGOPArray, deletionGOPArray, overallGCPArray); + + env->ReleaseByteArrayElements(overallGCP, overallGCPArray, JNI_ABORT); + env->ReleaseByteArrayElements(deletionGOP, deletionGOPArray, JNI_ABORT); + env->ReleaseByteArrayElements(insertionGOP, insertionGOPArray, JNI_ABORT); +} + +JNIEXPORT jdouble JNICALL +Java_org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM_jniInitializePriorsAndUpdateCells( + JNIEnv* env, jobject thisObject, + jboolean doInitialization, jint paddedReadLength, jint paddedHaplotypeLength, + jbyteArray readBases, jbyteArray haplotypeBases, jbyteArray readQuals, + jint hapStartIndex + ) +{ +#ifdef DEBUG3 + cout << "Entered mainCompute .. doInitialization "<<(doInitialization == JNI_TRUE)<<" hapStartIndex "<GetByteArrayElements(readBases, &is_copy); + jbyte* haplotypeBasesArray = (env)->GetByteArrayElements(haplotypeBases, &is_copy); + jbyte* readQualsArray = (env)->GetByteArrayElements(readQuals, &is_copy); +#ifdef DEBUG + assert(readBasesArray && "readBasesArray not initialized in JNI"); + assert(haplotypeBasesArray && "haplotypeBasesArray not initialized in JNI"); + assert(readQualsArray && "readQualsArray not initialized in JNI"); +#endif + testcase tc; + + tc.rslen = paddedReadLength-1; + tc.haplen = paddedHaplotypeLength-1; + + tc.rs = (char*)readBasesArray; + tc.hap = (char*)haplotypeBasesArray; + tc.q = (char*)readQualsArray; //TOASK - q is now char* + + compute_full_prob(&tc, g_double_dataholder.m_matchMatrix, g_double_dataholder.m_insertionMatrix, + g_double_dataholder.m_deletionMatrix, g_double_dataholder.m_transition, + doInitialization == JNI_TRUE, hapStartIndex, NULL); + + env->ReleaseByteArrayElements(readBases, readBasesArray, JNI_ABORT); + env->ReleaseByteArrayElements(haplotypeBases, haplotypeBasesArray, JNI_ABORT); + env->ReleaseByteArrayElements(readQuals, readQualsArray, JNI_ABORT); + return 0.0; +} + +JNIEXPORT jdouble JNICALL +Java_org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM_jniSubComputeReadLikelihoodGivenHaplotypeLog10( + JNIEnv* env, jobject thisObject, + jint readLength, jint haplotypeLength, + jbyteArray readBases, jbyteArray haplotypeBases, jbyteArray readQuals, + jbyteArray insertionGOP, jbyteArray deletionGOP, jbyteArray overallGCP, + jint hapStartIndex + ) +{ + jboolean is_copy = JNI_FALSE; + jbyte* readBasesArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(readBases, &is_copy); + jbyte* haplotypeBasesArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(haplotypeBases, &is_copy); + jbyte* readQualsArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(readQuals, &is_copy); + jbyte* insertionGOPArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(insertionGOP, &is_copy); + jbyte* deletionGOPArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(deletionGOP, &is_copy); + jbyte* overallGCPArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(overallGCP, &is_copy); +#ifdef DEBUG + assert(readBasesArray && "readBasesArray not initialized in JNI"); + assert(haplotypeBasesArray && "haplotypeBasesArray not initialized in JNI"); + assert(readQualsArray && "readQualsArray not initialized in JNI"); + assert(insertionGOPArray && "insertionGOP array not initialized in JNI"); + assert(deletionGOPArray && "deletionGOP array not initialized in JNI"); + assert(overallGCPArray && "OverallGCP array not initialized in JNI"); + assert(readLength < MROWS); +#endif + testcase tc; + tc.rslen = readLength; + tc.haplen = haplotypeLength; + tc.rs = (char*)readBasesArray; + tc.hap = (char*)haplotypeBasesArray; + for(unsigned i=0;i -/* Header for class org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM */ - -#ifndef _Included_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM -#define _Included_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM -#ifdef __cplusplus -extern "C" { -#endif -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_TRISTATE_CORRECTION -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_TRISTATE_CORRECTION 3.0 -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_matchToMatch -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_matchToMatch 0L -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_indelToMatch -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_indelToMatch 1L -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_matchToInsertion -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_matchToInsertion 2L -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_insertionToInsertion -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_insertionToInsertion 3L -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_matchToDeletion -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_matchToDeletion 4L -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_deletionToDeletion -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_deletionToDeletion 5L -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_debug -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_debug 0L -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_verify -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_verify 0L -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_debug0_1 -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_debug0_1 0L -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_debug1 -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_debug1 0L -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_debug2 -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_debug2 0L -#undef org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_debug3 -#define org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_debug3 0L -/* - * Class: org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM - * Method: jniGlobalInit - * Signature: (Ljava/lang/Class;Ljava/lang/Class;)V - */ -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniGlobalInit - (JNIEnv *, jobject, jclass, jclass); - -/* - * Class: org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM - * Method: jniClose - * Signature: ()V - */ -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniClose - (JNIEnv *, jobject); - -/* - * Class: org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM - * Method: jniInitialize - * Signature: (II)V - */ -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniInitialize - (JNIEnv *, jobject, jint, jint); - -/* - * Class: org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM - * Method: jniInitializeProbabilities - * Signature: ([[D[B[B[B)V - */ -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniInitializeProbabilities - (JNIEnv *, jclass, jobjectArray, jbyteArray, jbyteArray, jbyteArray); - -/* - * Class: org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM - * Method: jniInitializePriorsAndUpdateCells - * Signature: (ZII[B[B[BI)D - */ -JNIEXPORT jdouble JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniInitializePriorsAndUpdateCells - (JNIEnv *, jobject, jboolean, jint, jint, jbyteArray, jbyteArray, jbyteArray, jint); - -/* - * Class: org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM - * Method: jniSubComputeReadLikelihoodGivenHaplotypeLog10 - * Signature: (II[B[B[B[B[B[BI)D - */ -JNIEXPORT jdouble JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniSubComputeReadLikelihoodGivenHaplotypeLog10 - (JNIEnv *, jobject, jint, jint, jbyteArray, jbyteArray, jbyteArray, jbyteArray, jbyteArray, jbyteArray, jint); - -/* - * Class: org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM - * Method: jniInitializeHaplotypes - * Signature: (I[Lorg/broadinstitute/sting/utils/pairhmm/JNILoglessPairHMM/JNIHaplotypeDataHolderClass;)V - */ -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniInitializeHaplotypes - (JNIEnv *, jobject, jint, jobjectArray); - -/* - * Class: org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM - * Method: jniFinalizeRegion - * Signature: ()V - */ -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniFinalizeRegion - (JNIEnv *, jobject); - -/* - * Class: org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM - * Method: jniComputeLikelihoods - * Signature: (II[Lorg/broadinstitute/sting/utils/pairhmm/JNILoglessPairHMM/JNIReadDataHolderClass;[Lorg/broadinstitute/sting/utils/pairhmm/JNILoglessPairHMM/JNIHaplotypeDataHolderClass;[DI)V - */ -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniComputeLikelihoods - (JNIEnv *, jobject, jint, jint, jobjectArray, jobjectArray, jdoubleArray, jint); - -#ifdef __cplusplus -} -#endif -#endif diff --git a/PairHMM_JNI/org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM.cc b/PairHMM_JNI/org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc similarity index 72% rename from PairHMM_JNI/org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM.cc rename to PairHMM_JNI/org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc index a4cd63bc4..54df4265f 100644 --- a/PairHMM_JNI/org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM.cc +++ b/PairHMM_JNI/org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.cc @@ -1,120 +1,16 @@ #include "headers.h" -#include -#include "org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM.h" - - -#define ENABLE_ASSERTIONS 1 -#define DO_PROFILING 1 -//#define DEBUG 1 -//#define DEBUG0_1 1 -//#define DEBUG3 1 - +#include "jni_common.h" +#include "org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM.h" #include "template.h" #include "utils.h" #include "LoadTimeInitializer.h" using namespace std; -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniInitializeProbabilities -(JNIEnv* env, jclass thisObject, - jobjectArray transition, jbyteArray insertionGOP, jbyteArray deletionGOP, jbyteArray overallGCP - ) -{} - - JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniInitialize -(JNIEnv* env, jobject thisObject, - jint readMaxLength, jint haplotypeMaxLength) -{} - -JNIEXPORT jdouble JNICALL -Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniInitializePriorsAndUpdateCells( - JNIEnv* env, jobject thisObject, - jboolean doInitialization, jint paddedReadLength, jint paddedHaplotypeLength, - jbyteArray readBases, jbyteArray haplotypeBases, jbyteArray readQuals, - jint hapStartIndex - ) - -{ return 0.0; } - -#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 - -JNIEXPORT jdouble JNICALL -Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniSubComputeReadLikelihoodGivenHaplotypeLog10( - JNIEnv* env, jobject thisObject, - jint readLength, jint haplotypeLength, - jbyteArray readBases, jbyteArray haplotypeBases, jbyteArray readQuals, - jbyteArray insertionGOP, jbyteArray deletionGOP, jbyteArray overallGCP, - jint hapStartIndex - ) -{ - jboolean is_copy = JNI_FALSE; - jbyte* readBasesArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(readBases, &is_copy); - jbyte* haplotypeBasesArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(haplotypeBases, &is_copy); - jbyte* readQualsArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(readQuals, &is_copy); - jbyte* insertionGOPArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(insertionGOP, &is_copy); - jbyte* deletionGOPArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(deletionGOP, &is_copy); - jbyte* overallGCPArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(overallGCP, &is_copy); -#ifdef DEBUG - assert(readBasesArray && "readBasesArray not initialized in JNI"); - assert(haplotypeBasesArray && "haplotypeBasesArray not initialized in JNI"); - assert(readQualsArray && "readQualsArray not initialized in JNI"); - assert(insertionGOPArray && "insertionGOP array not initialized in JNI"); - assert(deletionGOPArray && "deletionGOP array not initialized in JNI"); - assert(overallGCPArray && "OverallGCP array not initialized in JNI"); - assert(readLength < MROWS); -#endif - testcase tc; - tc.rslen = readLength; - tc.haplen = haplotypeLength; - tc.rs = (char*)readBasesArray; - tc.hap = (char*)haplotypeBasesArray; - for(unsigned i=0;i > g_haplotypeBasesArrayVector; -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniInitializeHaplotypes +vector g_haplotypeBasesLengths; +JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM_jniInitializeHaplotypes (JNIEnv * env, jobject thisObject, jint numHaplotypes, jobjectArray haplotypeDataArray) { jboolean is_copy = JNI_FALSE; //To ensure, GET_BYTE_ARRAY_ELEMENTS is invoked only once for each haplotype, store bytearrays in a vector vector >& haplotypeBasesArrayVector = g_haplotypeBasesArrayVector; haplotypeBasesArrayVector.clear(); + g_haplotypeBasesLengths.clear(); haplotypeBasesArrayVector.resize(numHaplotypes); + g_haplotypeBasesLengths.resize(numHaplotypes); + jsize haplotypeBasesLength = 0; for(unsigned j=0;jGetObjectArrayElement(haplotypeDataArray, j); @@ -165,20 +65,22 @@ JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPai #endif env->DeleteLocalRef(haplotypeBases); //free the local reference jbyte* haplotypeBasesArray = (jbyte*)GET_BYTE_ARRAY_ELEMENTS(haplotypeBasesGlobalRef, &is_copy); + haplotypeBasesLength = env->GetArrayLength(haplotypeBasesGlobalRef); #ifdef ENABLE_ASSERTIONS assert(haplotypeBasesArray && "haplotypeBasesArray not initialized in JNI"); - assert(env->GetArrayLength(haplotypeBasesGlobalRef) < MCOLS); + assert(haplotypeBasesLength < MCOLS); #endif #ifdef DEBUG0_1 - cout << "JNI haplotype length "<GetArrayLength(haplotypeBasesGlobalRef)<<"\n"; + cout << "JNI haplotype length "<GetArrayLength(haplotypeBases);++k) + for(unsigned k=0;kGetArrayLength(haplotypeBasesGlobalRef); + g_load_time_initializer.m_sumHaplotypeLengths += haplotypeBasesLength; #endif } } @@ -188,7 +90,7 @@ JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPai //haplotypeDataArray - array of JNIHaplotypeDataHolderClass objects which contain the haplotypeBases //likelihoodArray - array of doubles to return results back to Java. Memory allocated by Java prior to JNI call //maxNumThreadsToUse - Max number of threads that OpenMP can use for the HMM computation -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniComputeLikelihoods +JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM_jniComputeLikelihoods (JNIEnv* env, jobject thisObject, jint numReads, jint numHaplotypes, jobjectArray readDataArray, jobjectArray haplotypeDataArray, jdoubleArray likelihoodArray, jint maxNumThreadsToUse) { @@ -265,7 +167,7 @@ JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPai for(unsigned j=0;jGetArrayLength(haplotypeBasesArrayVector[j].first); + jsize haplotypeLength = (jsize)g_haplotypeBasesLengths[j]; jbyte* haplotypeBasesArray = haplotypeBasesArrayVector[j].second; tc_array[tc_idx].rslen = (int)readLength; tc_array[tc_idx].haplen = (int)haplotypeLength; @@ -363,7 +265,7 @@ JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPai } //Release haplotypes at the end of a region -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniFinalizeRegion +JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM_jniFinalizeRegion (JNIEnv * env, jobject thisObject) { vector >& haplotypeBasesArrayVector = g_haplotypeBasesArrayVector; @@ -373,11 +275,12 @@ JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPai RELEASE_BYTE_ARRAY_ELEMENTS(haplotypeBasesArrayVector[j].first, haplotypeBasesArrayVector[j].second, JNI_RO_RELEASE_MODE); env->DeleteGlobalRef(haplotypeBasesArrayVector[j].first); //free the global reference } - haplotypeBasesArrayVector.clear(); + haplotypeBasesArrayVector.clear(); + g_haplotypeBasesLengths.clear(); } -JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniClose +JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM_jniClose (JNIEnv* env, jobject thisObject) { #ifdef DO_PROFILING diff --git a/PairHMM_JNI/run.sh b/PairHMM_JNI/run.sh index 65a53289a..a1eacdefa 100755 --- a/PairHMM_JNI/run.sh +++ b/PairHMM_JNI/run.sh @@ -1,14 +1,20 @@ #!/bin/bash rm -f *.txt *.log GSA_ROOT_DIR=/home/karthikg/broad/gsa-unstable +pair_hmm_implementation="VECTOR_LOGLESS_CACHING"; +if [ "$#" -ge 1 ]; +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}/PairHMM_JNI/ -jar $GSA_ROOT_DIR/dist/GenomeAnalysisTK.jar -T HaplotypeCaller \ +--dbsnp /data/broad/samples/joint_variant_calling/dbSNP/00-All.vcf \ -R /data/broad/samples/joint_variant_calling/broad_reference/human_g1k_v37_decoy.fasta \ -I /data/broad/samples/joint_variant_calling/NA12878_low_coverage_alignment/NA12878.chrom11.ILLUMINA.bwa.CEU.low_coverage.20121211.bam \ ---dbsnp /data/broad/samples/joint_variant_calling/dbSNP/00-All.vcf \ -stand_call_conf 50.0 \ -stand_emit_conf 10.0 \ ---pair_hmm_implementation JNI_LOGLESS_CACHING \ +--pair_hmm_implementation $pair_hmm_implementation \ -o output.raw.snps.indels.vcf #--pair_hmm_implementation JNI_LOGLESS_CACHING \ diff --git a/protected/java/src/org/broadinstitute/sting/gatk/walkers/haplotypecaller/PairHMMLikelihoodCalculationEngine.java b/protected/java/src/org/broadinstitute/sting/gatk/walkers/haplotypecaller/PairHMMLikelihoodCalculationEngine.java index 66d8ee33c..c534764a5 100644 --- a/protected/java/src/org/broadinstitute/sting/gatk/walkers/haplotypecaller/PairHMMLikelihoodCalculationEngine.java +++ b/protected/java/src/org/broadinstitute/sting/gatk/walkers/haplotypecaller/PairHMMLikelihoodCalculationEngine.java @@ -93,8 +93,8 @@ public class PairHMMLikelihoodCalculationEngine implements LikelihoodCalculation return new CnyPairHMM(); case VECTOR_LOGLESS_CACHING: return new VectorLoglessPairHMM(); - case DEBUG_JNI_LOGLESS_CACHING: - return new DebugJNILoglessPairHMM(hmmType); + case DEBUG_VECTOR_LOGLESS_CACHING: + return new DebugJNILoglessPairHMM(PairHMM.HMM_IMPLEMENTATION.VECTOR_LOGLESS_CACHING); case ARRAY_LOGLESS: if (noFpga || !CnyPairHMM.isAvailable()) return new ArrayLoglessPairHMM(); diff --git a/protected/java/src/org/broadinstitute/sting/utils/pairhmm/DebugJNILoglessPairHMM.java b/protected/java/src/org/broadinstitute/sting/utils/pairhmm/DebugJNILoglessPairHMM.java index 294c5c451..55491bcd8 100644 --- a/protected/java/src/org/broadinstitute/sting/utils/pairhmm/DebugJNILoglessPairHMM.java +++ b/protected/java/src/org/broadinstitute/sting/utils/pairhmm/DebugJNILoglessPairHMM.java @@ -76,7 +76,7 @@ import java.io.IOException; public class DebugJNILoglessPairHMM extends LoglessPairHMM { private static final boolean debug = false; //simulates ifdef - private static final boolean verify = debug || false; //simulates ifdef + private static final boolean verify = 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; diff --git a/protected/java/src/org/broadinstitute/sting/utils/pairhmm/JNILoglessPairHMM.java b/protected/java/src/org/broadinstitute/sting/utils/pairhmm/JNILoglessPairHMM.java index 1c7adf4c3..9b0772941 100644 --- a/protected/java/src/org/broadinstitute/sting/utils/pairhmm/JNILoglessPairHMM.java +++ b/protected/java/src/org/broadinstitute/sting/utils/pairhmm/JNILoglessPairHMM.java @@ -58,4 +58,5 @@ import java.util.HashMap; */ public abstract class JNILoglessPairHMM extends LoglessPairHMM { public abstract HashMap getHaplotypeToHaplotypeListIdxMap(); + protected long setupTime = 0; } diff --git a/protected/java/src/org/broadinstitute/sting/utils/pairhmm/VectorLoglessPairHMM.java b/protected/java/src/org/broadinstitute/sting/utils/pairhmm/VectorLoglessPairHMM.java index 9d6812c0c..6207c9679 100644 --- a/protected/java/src/org/broadinstitute/sting/utils/pairhmm/VectorLoglessPairHMM.java +++ b/protected/java/src/org/broadinstitute/sting/utils/pairhmm/VectorLoglessPairHMM.java @@ -70,6 +70,7 @@ public class VectorLoglessPairHMM extends JNILoglessPairHMM { //For machine capabilities public static final long sse42Mask = 1; public static final long avxMask = 2; + public static final long enableAll = 0xFFFFFFFFFFFFFFFFl; //Used to copy references to byteArrays to JNI from reads protected class JNIReadDataHolderClass { @@ -95,13 +96,15 @@ public class VectorLoglessPairHMM extends JNILoglessPairHMM { if(!isVectorLoglessPairHMMLibraryLoaded) { System.loadLibrary("VectorLoglessPairHMM"); isVectorLoglessPairHMMLibraryLoaded = true; - jniGlobalInit(JNIReadDataHolderClass.class, JNIHaplotypeDataHolderClass.class, 0xFFFFFFFFFFFFFFFFl); //need to do this only once + jniGlobalInit(JNIReadDataHolderClass.class, JNIHaplotypeDataHolderClass.class, enableAll); //need to do this only once } } @Override public void close() { + System.out.println("Time spent in setup for JNI call : "+setupTime); + super.close(); jniClose(); } @@ -154,6 +157,8 @@ public class VectorLoglessPairHMM extends JNILoglessPairHMM { */ @Override public PerReadAlleleLikelihoodMap computeLikelihoods( final List reads, final Map alleleHaplotypeMap, final Map GCPArrayMap ) { + if(doProfiling) + startTime = System.nanoTime(); int readListSize = reads.size(); int numHaplotypes = alleleHaplotypeMap.size(); int numTestcases = readListSize*numHaplotypes; @@ -170,12 +175,13 @@ public class VectorLoglessPairHMM extends JNILoglessPairHMM { ++idx; } - JNIHaplotypeDataHolderClass[] haplotypeDataArray = new JNIHaplotypeDataHolderClass[numHaplotypes]; mLikelihoodArray = new double[readListSize*numHaplotypes]; //to store results + if(doProfiling) + setupTime += (System.nanoTime() - startTime); //for(reads) // for(haplotypes) // compute_full_prob() - jniComputeLikelihoods(readListSize, numHaplotypes, readDataArray, haplotypeDataArray, mLikelihoodArray, 12); + jniComputeLikelihoods(readListSize, numHaplotypes, readDataArray, null, mLikelihoodArray, 12); final PerReadAlleleLikelihoodMap likelihoodMap = new PerReadAlleleLikelihoodMap(); idx = 0; @@ -193,6 +199,8 @@ public class VectorLoglessPairHMM extends JNILoglessPairHMM { } readIdx += numHaplotypes; } + if(doProfiling) + computeTime += (System.nanoTime() - startTime); return likelihoodMap; } } diff --git a/public/java/src/org/broadinstitute/sting/utils/pairhmm/PairHMM.java b/public/java/src/org/broadinstitute/sting/utils/pairhmm/PairHMM.java index 8b2123751..4ab17d5f3 100644 --- a/public/java/src/org/broadinstitute/sting/utils/pairhmm/PairHMM.java +++ b/public/java/src/org/broadinstitute/sting/utils/pairhmm/PairHMM.java @@ -59,8 +59,8 @@ public abstract class PairHMM { LOGLESS_CACHING, /* Optimized AVX implementation of LOGLESS_CACHING called through JNI */ VECTOR_LOGLESS_CACHING, - /* Debugging for any JNI implementation of LOGLESS_CACHING */ - DEBUG_JNI_LOGLESS_CACHING, + /* Debugging for vector implementation of LOGLESS_CACHING */ + DEBUG_VECTOR_LOGLESS_CACHING, /* Logless caching PairHMM that stores computations in 1D arrays instead of matrices, and which proceeds diagonally over the (read x haplotype) intersection matrix */ ARRAY_LOGLESS } @@ -77,6 +77,11 @@ public abstract class PairHMM { //debug array protected double[] mLikelihoodArray; + //profiling information + protected static final boolean doProfiling = true; + protected long computeTime = 0; + protected long startTime = 0; + /** * Initialize this PairHMM, making it suitable to run against a read and haplotype with given lengths * @@ -107,7 +112,7 @@ public abstract class PairHMM { */ public void finalizeRegion() { - ; + ; } /** @@ -119,7 +124,7 @@ public abstract class PairHMM { * @param readMaxLength the max length of reads we want to use with this PairHMM */ public void initialize( final List haplotypes, final Map> perSampleReadList, final int readMaxLength, final int haplotypeMaxLength ) { - initialize(readMaxLength, haplotypeMaxLength); + initialize(readMaxLength, haplotypeMaxLength); } protected int findMaxReadLength(final List reads) { @@ -152,6 +157,8 @@ public abstract class PairHMM { * said read coming from the said haplotype under the provided error model */ public PerReadAlleleLikelihoodMap computeLikelihoods(final List reads, final Map alleleHaplotypeMap, final Map GCPArrayMap) { + if(doProfiling) + startTime = System.nanoTime(); // (re)initialize the pairHMM only if necessary final int readMaxLength = findMaxReadLength(reads); @@ -159,8 +166,8 @@ public abstract class PairHMM { if (!initialized || readMaxLength > maxReadLength || haplotypeMaxLength > maxHaplotypeLength) { initialize(readMaxLength, haplotypeMaxLength); } final PerReadAlleleLikelihoodMap likelihoodMap = new PerReadAlleleLikelihoodMap(); - mLikelihoodArray = new double[reads.size()*alleleHaplotypeMap.size()]; - int idx = 0; + mLikelihoodArray = new double[reads.size()*alleleHaplotypeMap.size()]; + int idx = 0; for(GATKSAMRecord read : reads){ final byte[] readBases = read.getReadBases(); final byte[] readQuals = read.getBaseQualities(); @@ -173,16 +180,16 @@ public abstract class PairHMM { boolean isFirstHaplotype = true; Allele currentAllele = null; double log10l; - //for (final Allele allele : alleleHaplotypeMap.keySet()){ + //for (final Allele allele : alleleHaplotypeMap.keySet()){ for (Map.Entry currEntry : alleleHaplotypeMap.entrySet()){ - //final Haplotype haplotype = alleleHaplotypeMap.get(allele); - final Allele allele = currEntry.getKey(); - final Haplotype haplotype = currEntry.getValue(); + //final Haplotype haplotype = alleleHaplotypeMap.get(allele); + final Allele allele = currEntry.getKey(); + final Haplotype haplotype = currEntry.getValue(); final byte[] nextHaplotypeBases = haplotype.getBases(); if (currentHaplotypeBases != null) { log10l = computeReadLikelihoodGivenHaplotypeLog10(currentHaplotypeBases, readBases, readQuals, readInsQuals, readDelQuals, overallGCP, isFirstHaplotype, nextHaplotypeBases); - mLikelihoodArray[idx++] = log10l; + mLikelihoodArray[idx++] = log10l; likelihoodMap.add(read, currentAllele, log10l); } // update the current haplotype @@ -196,9 +203,11 @@ public abstract class PairHMM { log10l = computeReadLikelihoodGivenHaplotypeLog10(currentHaplotypeBases, readBases, readQuals, readInsQuals, readDelQuals, overallGCP, isFirstHaplotype, null); likelihoodMap.add(read, currentAllele, log10l); - mLikelihoodArray[idx++] = log10l; + mLikelihoodArray[idx++] = log10l; } } + if(doProfiling) + computeTime += (System.nanoTime() - startTime); return likelihoodMap; } @@ -305,5 +314,9 @@ public abstract class PairHMM { public double[] getLikelihoodArray() { return mLikelihoodArray; } //Called at the end of all HC calls - public void close() { ; } + public void close() + { + if(doProfiling) + System.out.println("Total compute time in PairHMM computeLikelihoods() : "+computeTime); + } }