Mercurial > repos > rliterman > csp2

comparison CSP2/CSP2_env/env-d9b9114564458d9d-741b3de822f2aaca6c6caa4325c4afce/opt/bbmap-39.01-1/current/icecream/IceCreamAlignerJava.java @ 68:5028fdace37b

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload commit 2e9511a184a1ca667c7be0c6321a36dc4e3d116d
author jpayne
date Tue, 18 Mar 2025 16:23:26 -0400
parents
children
comparison
equal deleted inserted replaced
67:0e9998148a16 68:5028fdace37b
1 package icecream;
2
3 import aligner.AlignmentResult;
4
5 public final class IceCreamAlignerJava extends IceCreamAligner {
6
7 IceCreamAlignerJava(){}
8
9 /**
10 * @param query
11 * @param ref
12 * @param rstart
13 * @param rstop
14 * @param minScore Quit early if score drops below this
15 * @param minRatio Don't return results if max score is less than this fraction of max possible score
16 * @return
17 */
18 @Override
19 public AlignmentResult alignForward(final byte[] query, final byte[] ref, final int rstart, final int rstop, final int minScore,
20 final float minRatio) {
21 final int arrayLength=rstop-rstart+1;
22
23 //Stack allocated; faster.
24 final int[] array1=new int[arrayLength+1], array2=new int[arrayLength+1];
25
26 final int minPassingScore=(int)(query.length*minRatio*pointsMatch);
27 final int minPassingScore3=minPassingScore-query.length*pointsMatch;
28
29 int[] prev=array1, next=array2;
30 int maxScore=-999999;
31 int maxQpos=-1;
32 int maxRpos=-1;
33
34 for(int i=0; i<=arrayLength-query.length; i++){prev[i]=0;}
35 for(int i=arrayLength-query.length, score=0; i<=arrayLength; i++, score+=pointsDel) {
36 prev[i]=score;
37 }
38
39 int currentRstart=rstart;
40 for(int qpos=0; qpos<query.length; qpos++){
41 prev[0]=pointsIns*qpos;
42
43 final byte q=query[qpos];
44 int maxScoreThisPass=-9999;
45 final int remainingBases=(query.length-qpos);
46 final int remainingPoints=remainingBases*pointsMatch;
47 final int minViableScore=minPassingScore3-remainingPoints;
48 // int minViableRstart=rstop+1;
49 for(int rpos=currentRstart, apos=1+currentRstart-rstart; rpos<=rstop; rpos++, apos++){
50 final byte r=ref[rpos];
51 final boolean match=(q==r);
52 final int vScore=prev[apos]+pointsIns;
53 final int hScore=next[apos-1]+pointsDel;
54 final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
55
56 //Slow branchy code
57 // final int score=Tools.max(vScore, hScore, dScore);
58 // next[apos]=score;
59 // if(score>=maxScoreThisPass){
60 // maxScoreThisPass=score;
61 // if(score>=maxScore){
62 // maxScore=score;
63 // maxQpos=qpos;
64 // maxRpos=rpos;
65 // }
66 // }
67
68 //Should be branchless conditional moves
69 int score=(dScore>=vScore ? dScore : vScore);
70 score=(hScore>score ? hScore : score);
71 next[apos]=score;
72 maxScoreThisPass=(score>maxScoreThisPass ? score : maxScoreThisPass);
73
74 // minViableRstart=((score<minViableScore || rpos>minViableRstart) ? minViableRstart : rpos);
75 }
76 iters+=arrayLength;
77 if(maxScoreThisPass<minScore){//Aggressive early exit
78 return null;
79 }
80
81 {//Safe early exit
82 if(maxScoreThisPass<minViableScore){return null;}
83 }
84
85 int[] temp=prev;
86 prev=next;
87 next=temp;
88 }
89
90 maxScore=-999999;
91 for(int rpos=rstart, apos=1; rpos<=rstop; rpos++, apos++){//Grab high score from last iteration
92 int score=prev[apos];
93 if(score>=maxScore){
94 maxScore=score;
95 maxQpos=query.length;
96 maxRpos=rpos;
97 }
98 }
99
100 maxScore=(maxScore-pointsSub*query.length)/(pointsMatch-pointsSub);//Rescale 0 to length
101 final float ratio=maxScore/(float)query.length;
102
103 if(ratio<minRatio){return null;}
104 return new AlignmentResult(maxScore, maxQpos, maxRpos, query.length, ref.length, rstart, rstop, ratio);
105 }
106
107 /**
108 * @param query
109 * @param ref
110 * @param qstart
111 * @param rstart
112 * @param rstop
113 * @param minScore Quit early if score drops below this
114 * @param minRatio Don't return results if max score is less than this fraction of max possible score
115 * @return
116 */
117 @Override
118 public AlignmentResult alignForwardShort(final byte[] query, final byte[] ref, final int rstart, final int rstop, final int minScore,
119 final float minRatio) {
120
121 final int arrayLength=query.length;
122
123 //Stack allocated; faster.
124 final int[] array1=new int[arrayLength+1], array2=new int[arrayLength+1];
125
126 final int minPassingScore=(int)(query.length*minRatio);
127 final int minPassingScore3=minPassingScore-query.length;
128
129 int[] prev=array1, next=array2;
130 int maxScore=-999999;
131 int maxQpos=-1;
132 int maxRpos=-1;
133
134 for(int i=0; i<prev.length; i++) {
135 prev[i]=pointsIns*i;
136 }
137
138 for(int rpos=rstart; rpos<=rstop; rpos++){
139 if(rstop-rpos<arrayLength){prev[0]=next[0]+pointsDel;}
140
141 final byte r=ref[rpos];
142 int score=0;
143
144 //In Java, the unsplit loop is faster
145 //-XX:+UseSuperWord seems to have no impact on speed.
146 //Might be worth packing the query and ref into int arrays and trying again.
147 for(int qpos=0, apos=1; qpos<arrayLength; qpos++, apos++){
148 final byte q=query[qpos];
149 final boolean match=(q==r);
150 final int vScore=prev[apos]+pointsIns;
151 final int hScore=next[apos-1]+pointsDel;
152 final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
153
154 score=(dScore>=vScore ? dScore : vScore);
155 score=(hScore>score ? hScore : score);
156
157 next[apos]=score;
158 }
159
160 // //Inner DV loop
161 // for(int qpos=0, apos=1; qpos<arrayLength; qpos++, apos++){
162 // final int q=query[qpos];
163 // final boolean match=(q==r);
164 // final int vScore=prev[apos]+pointsIns;
165 // final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
166 //
167 // score=(dScore>=vScore ? dScore : vScore);
168 //
169 // next[apos]=score;
170 // }
171 //
172 // //Inner I loop
173 // for(int qpos=0, apos=1; qpos<arrayLength; qpos++, apos++){
174 // final int hScore=next[apos-1]+pointsDel;
175 //
176 // score=next[apos];
177 // score=(hScore>score ? hScore : score);
178 //
179 // next[apos]=score;
180 // }
181
182 itersShort+=arrayLength;
183
184 if(score>=maxScore){
185 maxScore=score;
186 maxQpos=arrayLength-1;
187 maxRpos=rpos;
188 }
189
190 int[] temp=prev;
191 prev=next;
192 next=temp;
193 }
194
195 maxScore=(maxScore-pointsSub*query.length)/(pointsMatch-pointsSub);//Rescale 0 to length
196 final float ratio=maxScore/(float)query.length;
197
198 if(ratio<minRatio){return null;}
199
200 return new AlignmentResult(maxScore, maxQpos, maxRpos, query.length, ref.length, rstart, rstop, ratio);
201 }
202
203 /*--------------------------------------------------------------*/
204 /*---------------- Getters ----------------*/
205 /*--------------------------------------------------------------*/
206
207 @Override
208 long iters(){return iters;}
209
210 @Override
211 long itersShort(){return itersShort;}
212
213 /*--------------------------------------------------------------*/
214 /*---------------- Fields ----------------*/
215 /*--------------------------------------------------------------*/
216
217 long iters = 0;
218 long itersShort = 0;
219
220 /*--------------------------------------------------------------*/
221 /*---------------- Constants ----------------*/
222 /*--------------------------------------------------------------*/
223
224 public static final int pointsMatch = 1;
225 public static final int pointsSub = -1;
226 public static final int pointsDel = -2;
227 public static final int pointsIns = -2;
228
229 }