jpayne@68: package icecream;
jpayne@68: 
jpayne@68: import aligner.AlignmentResult;
jpayne@68: import shared.KillSwitch;
jpayne@68: import shared.Shared;
jpayne@68: 
jpayne@68: public final class IceCreamAlignerJNI extends IceCreamAligner {
jpayne@68: 
jpayne@68: 	static {
jpayne@68: 		Shared.loadJNI();
jpayne@68: 	}
jpayne@68: 
jpayne@68: 	IceCreamAlignerJNI(){}
jpayne@68: 	
jpayne@68: 	/**
jpayne@68: 	 * @param query
jpayne@68: 	 * @param ref
jpayne@68: 	 * @param qstart
jpayne@68: 	 * @param rstart
jpayne@68: 	 * @param rstop
jpayne@68: 	 * @param minScore Quit early if score drops below this
jpayne@68: 	 * @param minRatio Don't return results if max score is less than this fraction of max possible score
jpayne@68: 	 * @return
jpayne@68: 	 */
jpayne@68: 	@Override
jpayne@68: 	public AlignmentResult alignForward(final byte[] query, final byte[] ref, final int rstart, final int rstop, final int minScore,
jpayne@68: 			final float minRatio) {
jpayne@68: 		
jpayne@68: 		final int qlen=query.length;
jpayne@68: 		final int rlen=rstop-rstart+1;
jpayne@68: 		final int[] retVec=KillSwitch.allocInt1D(4);
jpayne@68: 		
jpayne@68: //		alignForwardJNI(qInt, rInt, retVec, qlen, rlen, minScore, minRatio);
jpayne@68: //		alignForwardPseudo(qInt, rInt, retVec, qlen, rlen, minScore, minRatio);
jpayne@68: 		
jpayne@68: 		if((qlen+rlen+32)*2>Short.MAX_VALUE) {
jpayne@68: 			final int[] qInt=new int[query.length];
jpayne@68: 			final int[] rInt=new int[rlen];
jpayne@68: 			for(int i=0; i<query.length; i++){qInt[i]=query[i];}
jpayne@68: 			for(int i=0; i<rlen; i++){rInt[i]=ref[i+rstart];}
jpayne@68: 
jpayne@68: //			alignForwardPseudo(qInt, rInt, retVec, qlen, rlen, minScore, minRatio);
jpayne@68: 			alignForwardJNI(qInt, rInt, retVec, qlen, rlen, minScore, minRatio);
jpayne@68: //			alignForwardShortJNI(qInt, rInt, retVec, qlen, rlen);
jpayne@68: 		}else{
jpayne@68: 			final short[] qInt=new short[query.length];
jpayne@68: 			final short[] rInt=new short[rlen];
jpayne@68: 			for(int i=0; i<query.length; i++){qInt[i]=query[i];}
jpayne@68: 			for(int i=0; i<rlen; i++){rInt[i]=ref[i+rstart];}
jpayne@68: 			
jpayne@68: 			alignForward16JNI(qInt, rInt, retVec, (short)qlen, (short)rlen, (short)minScore, minRatio);
jpayne@68: //			alignForwardShort16JNI(qInt, rInt, retVec, (short)qlen, (short)rlen);
jpayne@68: 		}
jpayne@68: 		
jpayne@68: 		int maxScore=retVec[0];
jpayne@68: 		final int maxQpos=retVec[1];
jpayne@68: 		final int maxRpos=retVec[2]+rstart;
jpayne@68: 		final int jniIters=retVec[3];
jpayne@68: 		
jpayne@68: 		iters+=jniIters;
jpayne@68: 		
jpayne@68: 		maxScore=(maxScore-pointsSub*query.length)/(pointsMatch-pointsSub);//Rescale 0 to length
jpayne@68: 		final float ratio=maxScore/(float)query.length;
jpayne@68: 		
jpayne@68: 		if(ratio<minRatio){return null;}
jpayne@68: 		
jpayne@68: 		return new AlignmentResult(maxScore, maxQpos, maxRpos, query.length, ref.length, rstart, rstop, ratio);
jpayne@68: 	}
jpayne@68: 
jpayne@68: 	/**
jpayne@68: 	 * @param query
jpayne@68: 	 * @param ref
jpayne@68: 	 * @param qstart
jpayne@68: 	 * @param rstart
jpayne@68: 	 * @param rstop
jpayne@68: 	 * @param minScore Quit early if score drops below this
jpayne@68: 	 * @param minRatio Don't return results if max score is less than this fraction of max possible score
jpayne@68: 	 * @return
jpayne@68: 	 */
jpayne@68: 	@Override
jpayne@68: 	public AlignmentResult alignForwardShort(final byte[] query, final byte[] ref, final int rstart, final int rstop, final int minScore,
jpayne@68: 			final float minRatio) {
jpayne@68: 		
jpayne@68: 		final int qlen=query.length;
jpayne@68: 		final int rlen=rstop-rstart+1;
jpayne@68: 		final int[] retVec=KillSwitch.allocInt1D(4);
jpayne@68: 		
jpayne@68: 		if((qlen+rlen+32)*2>Short.MAX_VALUE) {
jpayne@68: 			final int[] qInt=new int[query.length];
jpayne@68: 			final int[] rInt=new int[rlen];
jpayne@68: 			for(int i=0; i<query.length; i++){qInt[i]=query[i];}
jpayne@68: 			for(int i=0; i<rlen; i++){rInt[i]=ref[i+rstart];}
jpayne@68: 
jpayne@68: //			alignForwardShortPseudo(qInt, rInt, retVec, qlen, rlen);
jpayne@68: 			alignForwardShortJNI(qInt, rInt, retVec, qlen, rlen);
jpayne@68: 		}else{
jpayne@68: 			final short[] qInt=new short[query.length];
jpayne@68: 			final short[] rInt=new short[rlen];
jpayne@68: 			for(int i=0; i<query.length; i++){qInt[i]=query[i];}
jpayne@68: 			for(int i=0; i<rlen; i++){rInt[i]=ref[i+rstart];}
jpayne@68: 			
jpayne@68: 			alignForwardShort16JNI(qInt, rInt, retVec, (short)qlen, (short)rlen);
jpayne@68: 		}
jpayne@68: 		int maxScore=retVec[0];
jpayne@68: 		final int maxQpos=retVec[1];
jpayne@68: 		final int maxRpos=retVec[2]+rstart;
jpayne@68: 		final int jniIters=retVec[3];
jpayne@68: 		
jpayne@68: 		itersShort+=jniIters;
jpayne@68: 		
jpayne@68: 		maxScore=(maxScore-pointsSub*query.length)/(pointsMatch-pointsSub);//Rescale 0 to length
jpayne@68: 		final float ratio=maxScore/(float)query.length;
jpayne@68: 		
jpayne@68: 		if(ratio<minRatio){return null;}
jpayne@68: 		
jpayne@68: 		return new AlignmentResult(maxScore, maxQpos, maxRpos, query.length, ref.length, rstart, rstop, ratio);
jpayne@68: 	}
jpayne@68: 
jpayne@68: 	/*--------------------------------------------------------------*/
jpayne@68: 	/*----------------             JNI              ----------------*/
jpayne@68: 	/*--------------------------------------------------------------*/
jpayne@68: 	
jpayne@68: 	private static void alignForwardPseudo(final int[] query, final int[] ref, final int[] retArray, 
jpayne@68: 			final int qlen, final int rlen, final int minScore, final float minRatio) {
jpayne@68: 		
jpayne@68: 		final int arrayLength=rlen;
jpayne@68: 		final int arrayLength2=rlen+1;
jpayne@68: 
jpayne@68: 		//Stack allocated; faster.
jpayne@68: 		int[] array1=new int[arrayLength2];
jpayne@68: 		int[] array2=new int[arrayLength2];
jpayne@68: 
jpayne@68: 		int[] prev=array1;
jpayne@68: 		int[] next=array2;
jpayne@68: 		int maxScore=-32000;
jpayne@68: 		int maxQpos=-1;
jpayne@68: 		int maxRpos=-1;
jpayne@68: 		int iters=0;
jpayne@68: 
jpayne@68: 		final int minPassingScore=(int)(qlen*minRatio*pointsMatch);
jpayne@68: 		final int minPassingScore3=minPassingScore-qlen*pointsMatch;
jpayne@68: 
jpayne@68: 		for(int i=0; i<=arrayLength-qlen; i++){prev[i]=0;}
jpayne@68: 		for(int i=arrayLength-qlen, score=0; i<=arrayLength; i++, score+=pointsDel) {
jpayne@68: 			prev[i]=score;
jpayne@68: 			next[i]=0;
jpayne@68: 		}
jpayne@68: 		
jpayne@68: 		for(int qpos=0; qpos<qlen; qpos++){
jpayne@68: 			prev[0]=pointsIns*qpos;
jpayne@68: 
jpayne@68: 			final int q=query[qpos];
jpayne@68: 			int maxScoreThisPass=-32000;
jpayne@68: 			final int remainingBases=(qlen-qpos);
jpayne@68: 			final int remainingPoints=remainingBases*pointsMatch;
jpayne@68: 			final int minViableScore=minPassingScore3-remainingPoints;
jpayne@68: 			
jpayne@68: //			for(int rpos=0, apos=1; rpos<rlen; rpos++, apos++){
jpayne@68: //				final int r=ref[rpos];
jpayne@68: //				final boolean match=(q==r);
jpayne@68: //				final int vScore=prev[apos]+pointsIns;
jpayne@68: //				final int hScore=next[apos-1]+pointsDel;
jpayne@68: //				final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
jpayne@68: //
jpayne@68: //				//Should be branchless conditional moves
jpayne@68: //				int score=(dScore>=vScore ? dScore : vScore);
jpayne@68: //				score=(hScore>score ? hScore : score);
jpayne@68: //				next[apos]=score;
jpayne@68: //				maxScoreThisPass=(score>maxScoreThisPass ? score : maxScoreThisPass);
jpayne@68: //			}
jpayne@68: 			
jpayne@68: 			for(int rpos=0, apos=1; rpos<rlen; rpos++, apos++){
jpayne@68: 				final int r=ref[rpos];
jpayne@68: 				final boolean match=(q==r);
jpayne@68: 				final int vScore=prev[apos]+pointsIns;
jpayne@68: 				final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
jpayne@68: 
jpayne@68: 				//Should be branchless conditional moves
jpayne@68: 				int score=(dScore>=vScore ? dScore : vScore);
jpayne@68: 				next[apos]=score;
jpayne@68: 			}
jpayne@68: 			
jpayne@68: 			for(int apos=1; apos<arrayLength2; apos++){
jpayne@68: 				final int hScore=next[apos-1]+pointsDel;
jpayne@68: 
jpayne@68: 				//Should be branchless conditional moves
jpayne@68: 				int score=next[apos];
jpayne@68: 				score=(hScore>score ? hScore : score);
jpayne@68: 				next[apos]=score;
jpayne@68: 				maxScoreThisPass=(score>maxScoreThisPass ? score : maxScoreThisPass);
jpayne@68: 			}
jpayne@68: 			iters+=arrayLength;
jpayne@68: 			
jpayne@68: 			//Aggressive early exit
jpayne@68: 			if(maxScoreThisPass<minScore){return;}
jpayne@68: 
jpayne@68: 			//Safe early exit
jpayne@68: 			if(maxScoreThisPass<minViableScore){return;}
jpayne@68: 
jpayne@68: 			int[] temp=prev;
jpayne@68: 			prev=next;
jpayne@68: 			next=temp;
jpayne@68: 		}
jpayne@68: 
jpayne@68: 		maxScore=-32000;
jpayne@68: 		for(int apos=1; apos<arrayLength2; apos++){//Grab high score from last iteration
jpayne@68: 			int score=prev[apos];
jpayne@68: 			if(score>=maxScore){
jpayne@68: 				maxScore=score;
jpayne@68: 				maxQpos=qlen;
jpayne@68: 				maxRpos=apos-1;
jpayne@68: 			}
jpayne@68: 		}
jpayne@68: 		retArray[0]=maxScore;
jpayne@68: 		retArray[1]=maxQpos;
jpayne@68: 		retArray[2]=maxRpos;
jpayne@68: 		retArray[3]=iters;
jpayne@68: 	}
jpayne@68: 	
jpayne@68: 	private static void alignForwardShortPseudo(int[] query, int[] ref, int[] retArray, int qlen, int rlen) {
jpayne@68: 		
jpayne@68: 		final int arrayLength=qlen;
jpayne@68: 		final int arrayLength2=qlen+1;
jpayne@68: 
jpayne@68: 		//Stack allocated; faster.
jpayne@68: 		int[] array1=new int[arrayLength2];
jpayne@68: 		int[] array2=new int[arrayLength2];
jpayne@68: 
jpayne@68: 		int[] prev=array1;
jpayne@68: 		int[] next=array2;
jpayne@68: 		int maxScore=-999999;
jpayne@68: 		int maxQpos=-1;
jpayne@68: 		int maxRpos=-1;
jpayne@68: 		int itersShort=0;
jpayne@68: 
jpayne@68: 		for(int i=0; i<arrayLength2; i++) {
jpayne@68: 			prev[i]=pointsIns*i;
jpayne@68: 			next[i]=0;//For C version initialization
jpayne@68: 		}
jpayne@68: 
jpayne@68: 		for(int rpos=0; rpos<rlen; rpos++){
jpayne@68: 			if(rlen-rpos<arrayLength){prev[0]=next[0]+pointsDel;}
jpayne@68: 			
jpayne@68: 			final int r=ref[rpos];
jpayne@68: 			int score=0;
jpayne@68: 			
jpayne@68: //			//Inner loop
jpayne@68: //			for(int qpos=0, apos=1; qpos<arrayLength; qpos++, apos++){
jpayne@68: //				final int q=query[qpos];
jpayne@68: //				final boolean match=(q==r);
jpayne@68: //				final int vScore=prev[apos]+pointsIns;
jpayne@68: //				final int hScore=next[apos-1]+pointsDel;
jpayne@68: //				final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
jpayne@68: //
jpayne@68: //				score=(dScore>=vScore ? dScore : vScore);
jpayne@68: //				score=(hScore>score ? hScore : score);
jpayne@68: //
jpayne@68: //				next[apos]=score;
jpayne@68: //			}
jpayne@68: 			
jpayne@68: 			//Inner DV loop
jpayne@68: 			for(int qpos=0, apos=1; qpos<arrayLength; qpos++, apos++){
jpayne@68: 				final int q=query[qpos];
jpayne@68: 				final boolean match=(q==r);
jpayne@68: 				final int vScore=prev[apos]+pointsIns;
jpayne@68: 				final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
jpayne@68: 
jpayne@68: 				score=(dScore>=vScore ? dScore : vScore);
jpayne@68: 
jpayne@68: 				next[apos]=score;
jpayne@68: 			}
jpayne@68: 			
jpayne@68: 			//Inner I loop
jpayne@68: 			for(int qpos=0, apos=1; qpos<arrayLength; qpos++, apos++){
jpayne@68: 				final int hScore=next[apos-1]+pointsDel;
jpayne@68: 
jpayne@68: 				score=next[apos];
jpayne@68: 				score=(hScore>score ? hScore : score);
jpayne@68: 
jpayne@68: 				next[apos]=score;
jpayne@68: 			}
jpayne@68: 
jpayne@68: 			itersShort+=arrayLength;
jpayne@68: 
jpayne@68: 			if(score>=maxScore){
jpayne@68: 				maxScore=score;
jpayne@68: 				maxQpos=arrayLength-1;
jpayne@68: 				maxRpos=rpos;
jpayne@68: 			}
jpayne@68: 
jpayne@68: 			int[] temp=prev;
jpayne@68: 			prev=next;
jpayne@68: 			next=temp;
jpayne@68: 		}
jpayne@68: 		retArray[0]=maxScore;
jpayne@68: 		retArray[1]=maxQpos;
jpayne@68: 		retArray[2]=maxRpos;
jpayne@68: 		retArray[3]=itersShort;
jpayne@68: 	}
jpayne@68: 
jpayne@68: 	private static native void alignForwardJNI(int[] query, int[] ref, int[] retArray, int qlen, int rlen, int minScore, float minRatio);
jpayne@68: 	private static native void alignForward16JNI(short[] query, short[] ref, int[] retArray, short qlen, short rlen, short minScore, float minRatio);
jpayne@68: 	private static native void alignForwardShortJNI(int[] query, int[] ref, int[] retArray, int qlen, int rlen);
jpayne@68: 	private static native void alignForwardShort16JNI(short[] query, short[] ref, int[] retArray, short qlen, short rlen);
jpayne@68: 	
jpayne@68: 	/*--------------------------------------------------------------*/
jpayne@68: 	/*----------------           Getters            ----------------*/
jpayne@68: 	/*--------------------------------------------------------------*/
jpayne@68: 
jpayne@68: 	@Override
jpayne@68: 	long iters(){return iters;}
jpayne@68: 
jpayne@68: 	@Override
jpayne@68: 	long itersShort(){return itersShort;}
jpayne@68: 
jpayne@68: 	/*--------------------------------------------------------------*/
jpayne@68: 	/*----------------            Fields            ----------------*/
jpayne@68: 	/*--------------------------------------------------------------*/
jpayne@68: 	
jpayne@68: 	long iters = 0;
jpayne@68: 	long itersShort = 0;
jpayne@68: 	
jpayne@68: 	/*--------------------------------------------------------------*/
jpayne@68: 	/*----------------           Constants          ----------------*/
jpayne@68: 	/*--------------------------------------------------------------*/
jpayne@68: 	
jpayne@68: 	public static final int pointsMatch = 1;
jpayne@68: 	public static final int pointsSub = -1;
jpayne@68: 	public static final int pointsDel = -2;
jpayne@68: 	public static final int pointsIns = -2;
jpayne@68: 	
jpayne@68: }