diff CSP2/CSP2_env/env-d9b9114564458d9d-741b3de822f2aaca6c6caa4325c4afce/opt/bbmap-39.01-1/current/aligner/FlatAligner2.java @ 68:5028fdace37b

planemo upload commit 2e9511a184a1ca667c7be0c6321a36dc4e3d116d
author jpayne
date Tue, 18 Mar 2025 16:23:26 -0400
parents
children
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/CSP2/CSP2_env/env-d9b9114564458d9d-741b3de822f2aaca6c6caa4325c4afce/opt/bbmap-39.01-1/current/aligner/FlatAligner2.java	Tue Mar 18 16:23:26 2025 -0400
@@ -0,0 +1,248 @@
+package aligner;
+
+
+/** FlatAligner with flatter weights */
+public final class FlatAligner2 {
+
+	public FlatAligner2(){}
+	
+	/**
+	 * @param query
+	 * @param ref
+	 * @param qstart
+	 * @param rstart
+	 * @param rstop
+	 * @param minScore Quit early if score drops below this
+	 * @param minRatio Don't return results if max score is less than this fraction of max possible score
+	 * @return
+	 */
+	public AlignmentResult alignForward(final byte[] query, final byte[] ref, final int rstart, final int rstop, final int minScore,
+			final float minRatio) {
+//		if(ref.length<16300){return ica16.alignForward(query, ref, qstart, rstart, rstop, minScore, minRatio);}
+		final int arrayLength=rstop-rstart+1;
+		//		if(array1==null || arrayLength+1>array1.length){
+		//			array1=new int[(int)((arrayLength+2)*1.2)];
+		//			array2=new int[(int)((arrayLength+2)*1.2)];
+		//		}
+
+		//Stack allocated; faster.
+		final int[] array1=new int[arrayLength+1], array2=new int[arrayLength+1];
+
+		final int minPassingScore=(int)(query.length*minRatio*pointsMatch);
+		final int minPassingScore3=minPassingScore-query.length*pointsMatch;
+
+		int[] prev=array1, next=array2;
+		int maxScore=-999999;
+		int maxQpos=-1;
+		int maxRpos=-1;
+
+		for(int i=0; i<=arrayLength-query.length; i++){prev[i]=0;}
+		for(int i=arrayLength-query.length, score=0; i<=arrayLength; i++, score+=pointsDel) {
+			prev[i]=score;
+		}
+
+		int currentRstart=rstart;
+		for(int qpos=0; qpos<query.length; qpos++){
+			prev[0]=pointsIns*qpos;
+
+			final byte q=query[qpos];
+			int maxScoreThisPass=-9999;
+			final int remainingBases=(query.length-qpos);
+			final int remainingPoints=remainingBases*pointsMatch;
+			final int minViableScore=minPassingScore3-remainingPoints;
+			//			int minViableRstart=rstop+1;
+			for(int rpos=currentRstart, apos=1+currentRstart-rstart; rpos<=rstop; rpos++, apos++){
+				final byte r=ref[rpos];
+				final boolean match=(q==r);
+				final int vScore=prev[apos]+pointsIns;
+				final int hScore=next[apos-1]+pointsDel;
+				final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
+
+				//Slow branchy code
+				//				final int score=Tools.max(vScore, hScore, dScore);
+				//				next[apos]=score;
+				//				if(score>=maxScoreThisPass){
+				//					maxScoreThisPass=score;
+				//					if(score>=maxScore){
+				//						maxScore=score;
+				//						maxQpos=qpos;
+				//						maxRpos=rpos;
+				//					}
+				//				}
+
+				//Should be branchless conditional moves
+				int score=(dScore>=vScore ? dScore : vScore);
+				score=(hScore>score ? hScore : score);
+				next[apos]=score;
+				maxScoreThisPass=(score>maxScoreThisPass ? score : maxScoreThisPass);
+
+				//				minViableRstart=((score<minViableScore || rpos>minViableRstart) ? minViableRstart : rpos);
+			}
+			iters+=arrayLength;
+			//			currentRstart=minViableRstart;
+			//			System.err.println("qPos="+qpos+", maxScoreThisPass="+maxScoreThisPass+", maxScore="+maxScore+", minScore="+minScore);
+			//			System.err.println(Arrays.toString(prev));
+			//			System.err.println(Arrays.toString(next));
+			if(maxScoreThisPass<minScore){//Aggressive early exit
+				//				System.err.print('.');
+				//				System.err.println("qPos="+qpos+", maxScoreThisPass="+maxScoreThisPass+", maxScore="+maxScore+", minScore="+minScore);
+				return null;
+			}
+
+			{//Safe early exit
+				//				if((maxScoreThisPass+remaining+query.length)/2<minPassingScore){return null;}
+				//				if((maxScoreThisPass+remaining+query.length)<minPassingScore2){return null;}
+				if(maxScoreThisPass<minViableScore){return null;}
+			}
+
+			int[] temp=prev;
+			prev=next;
+			next=temp;
+		}
+		//		System.err.print('*');
+
+		maxScore=-999999;
+		for(int rpos=rstart, apos=1; rpos<=rstop; rpos++, apos++){//Grab high score from last iteration
+			int score=prev[apos];
+			if(score>=maxScore){
+				maxScore=score;
+				maxQpos=query.length;
+				maxRpos=rpos;
+			}
+		}
+
+
+		//		maxScore=(maxScore+query.length)/2;//Rescale 0 to length
+		maxScore=(maxScore-pointsSub*query.length)/(pointsMatch-pointsSub);//Rescale 0 to length
+		final float ratio=maxScore/(float)query.length;
+		//		System.err.println("maxqPos="+maxQpos+", maxScore="+maxScore+", minScore="+(minRatio*query.length));
+		//					if(minRatio*query.length>maxScore){return null;}
+		if(ratio<minRatio){return null;}
+		//		System.err.print('^');
+		return new AlignmentResult(maxScore, maxQpos, maxRpos, query.length, ref.length, rstart, rstop, ratio);
+	}
+
+	/**
+	 * @param query
+	 * @param ref
+	 * @param qstart
+	 * @param rstart
+	 * @param rstop
+	 * @param minScore Quit early if score drops below this
+	 * @param minRatio Don't return results if max score is less than this fraction of max possible score
+	 * @return
+	 */
+	public AlignmentResult alignForwardShort(final byte[] query, final byte[] ref, final int rstart, final int rstop,
+			final float minRatio) {
+//		if(ref.length<16300){return ica16.alignForwardShort(query, ref, qstart, rstart, rstop, minScore, minRatio);}
+
+		final int arrayLength=query.length;
+		//		if(array1==null || arrayLength+1>array1.length){
+		//			array1=new int[(int)((arrayLength+2)*1.2)];
+		//			array2=new int[(int)((arrayLength+2)*1.2)];
+		//		}
+
+		//Stack allocated; faster.
+		final int[] array1=new int[arrayLength+1], array2=new int[arrayLength+1];
+
+		final int minPassingScore=(int)(pointsMatch*query.length*minRatio);
+
+		int[] prev=array1, next=array2;
+		int maxScore=-1;
+		int maxQpos=-1;
+		int maxRpos=-1;
+
+		for(int i=0; i<prev.length; i++) {
+			prev[i]=pointsIns*i;
+		}
+		
+		for(int rpos=rstart; rpos<=rstop && maxScore<minPassingScore; rpos++){
+			if(rstop-rpos<arrayLength){prev[0]=next[0]+pointsDel;}
+
+			final byte r=ref[rpos];
+			int score=0;
+			for(int qpos=0, apos=1; qpos<arrayLength; qpos++, apos++){
+				final byte q=query[qpos];
+				final boolean match=(q==r);
+				final int vScore=prev[apos]+pointsIns;
+				final int hScore=next[apos-1]+pointsDel;
+				final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
+
+				//				score=Tools.max(vScore, hScore, dScore);
+				score=(dScore>=vScore ? dScore : vScore);
+				score=(hScore>score ? hScore : score);
+
+				next[apos]=score;
+			}
+			itersShort+=arrayLength;
+
+			if(score>=maxScore){
+				maxScore=score;
+				maxQpos=arrayLength-1;
+				maxRpos=rpos;
+			}
+
+			//			System.err.println("qPos="+qpos+", maxScoreThisPass="+maxScoreThisPass+", maxScore="+maxScore+", minScore="+minScore);
+			//			System.err.println(Arrays.toString(prev));
+			//			System.err.println(Arrays.toString(next));
+			//			if(maxScoreThisPass<minScore){//Aggressive early exit
+			////				System.err.print('.');
+			////				System.err.println("qPos="+qpos+", maxScoreThisPass="+maxScoreThisPass+", maxScore="+maxScore+", minScore="+minScore);
+			//				return null;
+			//			}
+
+			//			{//Safe early exit
+			//				int remaining=query.length-qpos;
+			////				if((maxScoreThisPass+remaining+query.length)/2<minPassingScore){return null;}
+			////				if((maxScoreThisPass+remaining+query.length)<minPassingScore2){return null;}
+			//				if((maxScoreThisPass+remaining)<minPassingScore3){return null;}
+			//			}
+
+			int[] temp=prev;
+			prev=next;
+			next=temp;
+		}
+		//		System.err.print('*');
+
+		//		maxScore=-999999;
+		//		for(int rpos=rstart, apos=1; rpos<=rstop; rpos++, apos++){//Grab high score from last iteration
+		//			int score=prev[apos];
+		//			if(score>=maxScore){
+		//				maxScore=score;
+		//				maxQpos=query.length;
+		//				maxRpos=rpos;
+		//			}
+		//		}
+		final float ratio=maxScore/(float)(pointsMatch*query.length);
+//		System.err.println(ratio+"\t"+maxScore+"\t"+query.length);
+		//		System.err.println("maxqPos="+maxQpos+", maxScore="+maxScore+", minScore="+(minRatio*query.length));
+		//					if(minRatio*query.length>maxScore){return null;}
+		if(ratio<minRatio){return null;}
+		//		System.err.print('^');
+		return new AlignmentResult(maxScore, maxQpos, maxRpos, query.length, ref.length, rstart, rstop, ratio);
+	}
+
+	/*--------------------------------------------------------------*/
+	/*----------------           Getters            ----------------*/
+	/*--------------------------------------------------------------*/
+
+	long iters(){return iters;}
+	long itersShort(){return itersShort;}
+
+	/*--------------------------------------------------------------*/
+	/*----------------            Fields            ----------------*/
+	/*--------------------------------------------------------------*/
+	
+	long iters = 0;
+	long itersShort = 0;
+	
+	/*--------------------------------------------------------------*/
+	/*----------------           Constants          ----------------*/
+	/*--------------------------------------------------------------*/
+	
+	public static final int pointsMatch = 10;
+	public static final int pointsSub = -9;
+	public static final int pointsDel = -11;
+	public static final int pointsIns = -11;
+	
+}