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
This commit is contained in:
Karthik Gururaj
parent
81bdfbd00d
commit
018e9e2c5f
|
|
@ -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}
|
||||
|
||||
|
||||
|
|
|
|||
|
|
@ -0,0 +1,32 @@
|
|||
#ifndef JNI_COMMON_H
|
||||
#define JNI_COMMON_H
|
||||
|
||||
#include <jni.h>
|
||||
#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
|
||||
|
|
@ -0,0 +1,166 @@
|
|||
#ifndef JNI_DEBUG_H
|
||||
#define JNI_DEBUG_H
|
||||
|
||||
template<class NUMBER>
|
||||
class DataHolder
|
||||
{
|
||||
#define INIT_MATRIX(X) \
|
||||
X = new NUMBER*[m_paddedMaxReadLength]; \
|
||||
for(int i=0;i<m_paddedMaxReadLength;++i) \
|
||||
{ \
|
||||
X[i] = new NUMBER[m_paddedMaxHaplotypeLength]; \
|
||||
for(int j=0;j<m_paddedMaxHaplotypeLength;++j) \
|
||||
X[i][j] = (NUMBER)0; \
|
||||
}
|
||||
|
||||
#define FREE_MATRIX(X) \
|
||||
for(int i=0;i<m_paddedMaxReadLength;++i) \
|
||||
delete X[i]; \
|
||||
delete X;
|
||||
|
||||
public:
|
||||
DataHolder() { m_is_initialized = false; }
|
||||
void initialize(int readMaxLength, int haplotypeMaxLength)
|
||||
{
|
||||
if(m_is_initialized)
|
||||
{
|
||||
FREE_MATRIX(m_matchMatrix);
|
||||
FREE_MATRIX(m_insertionMatrix);
|
||||
FREE_MATRIX(m_deletionMatrix);
|
||||
FREE_MATRIX(m_prior);
|
||||
delete m_transition;
|
||||
}
|
||||
|
||||
m_readMaxLength = readMaxLength;
|
||||
m_haplotypeMaxLength = haplotypeMaxLength;
|
||||
m_paddedMaxReadLength = readMaxLength + 1;
|
||||
m_paddedMaxHaplotypeLength = haplotypeMaxLength + 1;
|
||||
|
||||
INIT_MATRIX(m_matchMatrix);
|
||||
INIT_MATRIX(m_insertionMatrix);
|
||||
INIT_MATRIX(m_deletionMatrix);
|
||||
INIT_MATRIX(m_prior);
|
||||
m_transition = new NUMBER[m_paddedMaxReadLength][6];
|
||||
for(int i=0;i<m_paddedMaxReadLength;++i)
|
||||
for(int j=0;j<6;++j)
|
||||
m_transition[i][j] = (NUMBER)0;
|
||||
m_is_initialized = true;
|
||||
}
|
||||
|
||||
//Corresponds to initializeProbabilities
|
||||
void initializeProbabilities(jint length, jbyte* insertionGOP, jbyte* deletionGOP, jbyte* overallGCP)
|
||||
{
|
||||
static unsigned g_num_prob_init = 0;
|
||||
Context<NUMBER> 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<double> g_double_dataholder;
|
||||
|
||||
template<class NUMBER>
|
||||
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<NUMBER> 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
|
||||
|
|
@ -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<double> 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 "<<readMaxLength<<" haplotypeMaxLength "<<haplotypeMaxLength<<"\n";
|
||||
#endif
|
||||
g_double_dataholder.initialize(readMaxLength, haplotypeMaxLength);
|
||||
#ifdef DEBUG3
|
||||
debug_dump("lengths_jni.txt", to_string(readMaxLength)+" "+to_string(haplotypeMaxLength),true);
|
||||
#endif
|
||||
++g_num_init_calls;
|
||||
}
|
||||
|
||||
JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_DebugJNILoglessPairHMM_jniInitializeProbabilities
|
||||
(JNIEnv* env, jclass thisObject,
|
||||
jobjectArray transition, jbyteArray insertionGOP, jbyteArray deletionGOP, jbyteArray overallGCP
|
||||
)
|
||||
{
|
||||
jboolean is_copy = JNI_FALSE;
|
||||
jsize length = (env)->GetArrayLength(insertionGOP);
|
||||
#ifdef DEBUG3
|
||||
cout << "Entered initializeProbabilities .. length "<<length<<"\n";
|
||||
#endif
|
||||
jbyte* insertionGOPArray = (env)->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 "<<hapStartIndex<<"\n";
|
||||
cout << "mainCompute padded lengths "<< paddedReadLength << " " << paddedHaplotypeLength <<"\n";
|
||||
#endif
|
||||
jboolean is_copy = JNI_FALSE;
|
||||
jbyte* readBasesArray = (env)->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<double>(&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<readLength;++i)
|
||||
{
|
||||
tc.q[i] = (int)readQualsArray[i];
|
||||
tc.i[i] = (int)insertionGOPArray[i];
|
||||
tc.d[i] = (int)deletionGOPArray[i];
|
||||
tc.c[i] = (int)overallGCPArray[i];
|
||||
}
|
||||
|
||||
double result_avxd = g_compute_full_prob_double(&tc, 0);
|
||||
double result = log10(result_avxd) - log10(ldexp(1.0, 1020));
|
||||
#ifdef DEBUG
|
||||
g_load_time_initializer.debug_dump("return_values_jni.txt",to_string(result),true);
|
||||
#endif
|
||||
|
||||
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(overallGCP, overallGCPArray, JNI_RO_RELEASE_MODE);
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(deletionGOP, deletionGOPArray, JNI_RO_RELEASE_MODE);
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(insertionGOP, insertionGOPArray, JNI_RO_RELEASE_MODE);
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(readQuals, readQualsArray, JNI_RO_RELEASE_MODE);
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(haplotypeBases, haplotypeBasesArray, JNI_RO_RELEASE_MODE);
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(readBases, readBasesArray, JNI_RO_RELEASE_MODE);
|
||||
|
||||
return 0.0;
|
||||
}
|
||||
|
||||
|
|
@ -1,111 +0,0 @@
|
|||
/* DO NOT EDIT THIS FILE - it is machine generated */
|
||||
#include <jni.h>
|
||||
/* 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
|
||||
|
|
@ -1,120 +1,16 @@
|
|||
#include "headers.h"
|
||||
#include <jni.h>
|
||||
#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<readLength;++i)
|
||||
{
|
||||
tc.q[i] = (int)readQualsArray[i];
|
||||
tc.i[i] = (int)insertionGOPArray[i];
|
||||
tc.d[i] = (int)deletionGOPArray[i];
|
||||
tc.c[i] = (int)overallGCPArray[i];
|
||||
}
|
||||
|
||||
double result_avxd = g_compute_full_prob_double(&tc, 0);
|
||||
double result = log10(result_avxd) - log10(ldexp(1.0, 1020));
|
||||
#ifdef DEBUG
|
||||
g_load_time_initializer.debug_dump("return_values_jni.txt",to_string(result),true);
|
||||
#endif
|
||||
|
||||
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(overallGCP, overallGCPArray, JNI_RO_RELEASE_MODE);
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(deletionGOP, deletionGOPArray, JNI_RO_RELEASE_MODE);
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(insertionGOP, insertionGOPArray, JNI_RO_RELEASE_MODE);
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(readQuals, readQualsArray, JNI_RO_RELEASE_MODE);
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(haplotypeBases, haplotypeBasesArray, JNI_RO_RELEASE_MODE);
|
||||
RELEASE_BYTE_ARRAY_ELEMENTS(readBases, readBasesArray, JNI_RO_RELEASE_MODE);
|
||||
|
||||
return 0.0;
|
||||
}
|
||||
|
||||
//Should be called only once for the whole Java process - initializes field ids for the classes JNIReadDataHolderClass
|
||||
//and JNIHaplotypeDataHolderClass
|
||||
JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniGlobalInit
|
||||
(JNIEnv* env, jobject thisObject, jclass readDataHolderClass, jclass haplotypeDataHolderClass)
|
||||
JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_VectorLoglessPairHMM_jniGlobalInit
|
||||
(JNIEnv* env, jobject thisObject, jclass readDataHolderClass, jclass haplotypeDataHolderClass, jlong mask)
|
||||
{
|
||||
assert(readDataHolderClass);
|
||||
assert(haplotypeDataHolderClass);
|
||||
|
|
@ -143,14 +39,18 @@ JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPai
|
|||
//Since the list of haplotypes against which the reads are evaluated in PairHMM is the same for a region,
|
||||
//transfer the list only once
|
||||
vector<pair<jbyteArray, jbyte*> > g_haplotypeBasesArrayVector;
|
||||
JNIEXPORT void JNICALL Java_org_broadinstitute_sting_utils_pairhmm_JNILoglessPairHMM_jniInitializeHaplotypes
|
||||
vector<unsigned> 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<pair<jbyteArray, jbyte*> >& haplotypeBasesArrayVector = g_haplotypeBasesArrayVector;
|
||||
haplotypeBasesArrayVector.clear();
|
||||
g_haplotypeBasesLengths.clear();
|
||||
haplotypeBasesArrayVector.resize(numHaplotypes);
|
||||
g_haplotypeBasesLengths.resize(numHaplotypes);
|
||||
jsize haplotypeBasesLength = 0;
|
||||
for(unsigned j=0;j<numHaplotypes;++j)
|
||||
{
|
||||
jobject haplotypeObject = env->GetObjectArrayElement(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 "<<env->GetArrayLength(haplotypeBasesGlobalRef)<<"\n";
|
||||
cout << "JNI haplotype length "<<haplotypeBasesLength<<"\n";
|
||||
#endif
|
||||
haplotypeBasesArrayVector[j] = make_pair(haplotypeBasesGlobalRef, haplotypeBasesArray);
|
||||
g_haplotypeBasesLengths[j] = haplotypeBasesLength;
|
||||
#ifdef DEBUG3
|
||||
for(unsigned k=0;k<env->GetArrayLength(haplotypeBases);++k)
|
||||
for(unsigned k=0;k<haplotypeBasesLength;++k)
|
||||
g_load_time_initializer.debug_dump("haplotype_bases_jni.txt",to_string((int)haplotypeBasesArray[k]),true);
|
||||
#endif
|
||||
#ifdef DO_PROFILING
|
||||
g_load_time_initializer.m_sumHaplotypeLengths += env->GetArrayLength(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;j<numHaplotypes;++j)
|
||||
{
|
||||
jsize haplotypeLength = env->GetArrayLength(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<pair<jbyteArray, jbyte*> >& 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
|
||||
|
|
@ -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 \
|
||||
|
|
|
|||
|
|
@ -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();
|
||||
|
|
|
|||
|
|
@ -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;
|
||||
|
|
|
|||
|
|
@ -58,4 +58,5 @@ import java.util.HashMap;
|
|||
*/
|
||||
public abstract class JNILoglessPairHMM extends LoglessPairHMM {
|
||||
public abstract HashMap<Haplotype, Integer> getHaplotypeToHaplotypeListIdxMap();
|
||||
protected long setupTime = 0;
|
||||
}
|
||||
|
|
|
|||
|
|
@ -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<GATKSAMRecord> reads, final Map<Allele, Haplotype> alleleHaplotypeMap, final Map<GATKSAMRecord, byte[]> 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;
|
||||
}
|
||||
}
|
||||
|
|
|
|||
|
|
@ -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<Haplotype> haplotypes, final Map<String, List<GATKSAMRecord>> perSampleReadList, final int readMaxLength, final int haplotypeMaxLength ) {
|
||||
initialize(readMaxLength, haplotypeMaxLength);
|
||||
initialize(readMaxLength, haplotypeMaxLength);
|
||||
}
|
||||
|
||||
protected int findMaxReadLength(final List<GATKSAMRecord> 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<GATKSAMRecord> reads, final Map<Allele, Haplotype> alleleHaplotypeMap, final Map<GATKSAMRecord, byte[]> 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<Allele,Haplotype> 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);
|
||||
}
|
||||
}
|
||||
|
|
|
|||
Loading…
Reference in New Issue