Mercurial > repos > rliterman > csp2

comparison CSP2/CSP2_env/env-d9b9114564458d9d-741b3de822f2aaca6c6caa4325c4afce/opt/bbmap-39.01-1/current/icecream/IceCreamAlignerJNI.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 import shared.KillSwitch;
5 import shared.Shared;
6
7 public final class IceCreamAlignerJNI extends IceCreamAligner {
8
9 static {
10 Shared.loadJNI();
11 }
12
13 IceCreamAlignerJNI(){}
14
15 /**
16 * @param query
17 * @param ref
18 * @param qstart
19 * @param rstart
20 * @param rstop
21 * @param minScore Quit early if score drops below this
22 * @param minRatio Don't return results if max score is less than this fraction of max possible score
23 * @return
24 */
25 @Override
26 public AlignmentResult alignForward(final byte[] query, final byte[] ref, final int rstart, final int rstop, final int minScore,
27 final float minRatio) {
28
29 final int qlen=query.length;
30 final int rlen=rstop-rstart+1;
31 final int[] retVec=KillSwitch.allocInt1D(4);
32
33 // alignForwardJNI(qInt, rInt, retVec, qlen, rlen, minScore, minRatio);
34 // alignForwardPseudo(qInt, rInt, retVec, qlen, rlen, minScore, minRatio);
35
36 if((qlen+rlen+32)*2>Short.MAX_VALUE) {
37 final int[] qInt=new int[query.length];
38 final int[] rInt=new int[rlen];
39 for(int i=0; i<query.length; i++){qInt[i]=query[i];}
40 for(int i=0; i<rlen; i++){rInt[i]=ref[i+rstart];}
41
42 // alignForwardPseudo(qInt, rInt, retVec, qlen, rlen, minScore, minRatio);
43 alignForwardJNI(qInt, rInt, retVec, qlen, rlen, minScore, minRatio);
44 // alignForwardShortJNI(qInt, rInt, retVec, qlen, rlen);
45 }else{
46 final short[] qInt=new short[query.length];
47 final short[] rInt=new short[rlen];
48 for(int i=0; i<query.length; i++){qInt[i]=query[i];}
49 for(int i=0; i<rlen; i++){rInt[i]=ref[i+rstart];}
50
51 alignForward16JNI(qInt, rInt, retVec, (short)qlen, (short)rlen, (short)minScore, minRatio);
52 // alignForwardShort16JNI(qInt, rInt, retVec, (short)qlen, (short)rlen);
53 }
54
55 int maxScore=retVec[0];
56 final int maxQpos=retVec[1];
57 final int maxRpos=retVec[2]+rstart;
58 final int jniIters=retVec[3];
59
60 iters+=jniIters;
61
62 maxScore=(maxScore-pointsSub*query.length)/(pointsMatch-pointsSub);//Rescale 0 to length
63 final float ratio=maxScore/(float)query.length;
64
65 if(ratio<minRatio){return null;}
66
67 return new AlignmentResult(maxScore, maxQpos, maxRpos, query.length, ref.length, rstart, rstop, ratio);
68 }
69
70 /**
71 * @param query
72 * @param ref
73 * @param qstart
74 * @param rstart
75 * @param rstop
76 * @param minScore Quit early if score drops below this
77 * @param minRatio Don't return results if max score is less than this fraction of max possible score
78 * @return
79 */
80 @Override
81 public AlignmentResult alignForwardShort(final byte[] query, final byte[] ref, final int rstart, final int rstop, final int minScore,
82 final float minRatio) {
83
84 final int qlen=query.length;
85 final int rlen=rstop-rstart+1;
86 final int[] retVec=KillSwitch.allocInt1D(4);
87
88 if((qlen+rlen+32)*2>Short.MAX_VALUE) {
89 final int[] qInt=new int[query.length];
90 final int[] rInt=new int[rlen];
91 for(int i=0; i<query.length; i++){qInt[i]=query[i];}
92 for(int i=0; i<rlen; i++){rInt[i]=ref[i+rstart];}
93
94 // alignForwardShortPseudo(qInt, rInt, retVec, qlen, rlen);
95 alignForwardShortJNI(qInt, rInt, retVec, qlen, rlen);
96 }else{
97 final short[] qInt=new short[query.length];
98 final short[] rInt=new short[rlen];
99 for(int i=0; i<query.length; i++){qInt[i]=query[i];}
100 for(int i=0; i<rlen; i++){rInt[i]=ref[i+rstart];}
101
102 alignForwardShort16JNI(qInt, rInt, retVec, (short)qlen, (short)rlen);
103 }
104 int maxScore=retVec[0];
105 final int maxQpos=retVec[1];
106 final int maxRpos=retVec[2]+rstart;
107 final int jniIters=retVec[3];
108
109 itersShort+=jniIters;
110
111 maxScore=(maxScore-pointsSub*query.length)/(pointsMatch-pointsSub);//Rescale 0 to length
112 final float ratio=maxScore/(float)query.length;
113
114 if(ratio<minRatio){return null;}
115
116 return new AlignmentResult(maxScore, maxQpos, maxRpos, query.length, ref.length, rstart, rstop, ratio);
117 }
118
119 /*--------------------------------------------------------------*/
120 /*---------------- JNI ----------------*/
121 /*--------------------------------------------------------------*/
122
123 private static void alignForwardPseudo(final int[] query, final int[] ref, final int[] retArray,
124 final int qlen, final int rlen, final int minScore, final float minRatio) {
125
126 final int arrayLength=rlen;
127 final int arrayLength2=rlen+1;
128
129 //Stack allocated; faster.
130 int[] array1=new int[arrayLength2];
131 int[] array2=new int[arrayLength2];
132
133 int[] prev=array1;
134 int[] next=array2;
135 int maxScore=-32000;
136 int maxQpos=-1;
137 int maxRpos=-1;
138 int iters=0;
139
140 final int minPassingScore=(int)(qlen*minRatio*pointsMatch);
141 final int minPassingScore3=minPassingScore-qlen*pointsMatch;
142
143 for(int i=0; i<=arrayLength-qlen; i++){prev[i]=0;}
144 for(int i=arrayLength-qlen, score=0; i<=arrayLength; i++, score+=pointsDel) {
145 prev[i]=score;
146 next[i]=0;
147 }
148
149 for(int qpos=0; qpos<qlen; qpos++){
150 prev[0]=pointsIns*qpos;
151
152 final int q=query[qpos];
153 int maxScoreThisPass=-32000;
154 final int remainingBases=(qlen-qpos);
155 final int remainingPoints=remainingBases*pointsMatch;
156 final int minViableScore=minPassingScore3-remainingPoints;
157
158 // for(int rpos=0, apos=1; rpos<rlen; rpos++, apos++){
159 // final int r=ref[rpos];
160 // final boolean match=(q==r);
161 // final int vScore=prev[apos]+pointsIns;
162 // final int hScore=next[apos-1]+pointsDel;
163 // final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
164 //
165 // //Should be branchless conditional moves
166 // int score=(dScore>=vScore ? dScore : vScore);
167 // score=(hScore>score ? hScore : score);
168 // next[apos]=score;
169 // maxScoreThisPass=(score>maxScoreThisPass ? score : maxScoreThisPass);
170 // }
171
172 for(int rpos=0, apos=1; rpos<rlen; rpos++, apos++){
173 final int r=ref[rpos];
174 final boolean match=(q==r);
175 final int vScore=prev[apos]+pointsIns;
176 final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
177
178 //Should be branchless conditional moves
179 int score=(dScore>=vScore ? dScore : vScore);
180 next[apos]=score;
181 }
182
183 for(int apos=1; apos<arrayLength2; apos++){
184 final int hScore=next[apos-1]+pointsDel;
185
186 //Should be branchless conditional moves
187 int score=next[apos];
188 score=(hScore>score ? hScore : score);
189 next[apos]=score;
190 maxScoreThisPass=(score>maxScoreThisPass ? score : maxScoreThisPass);
191 }
192 iters+=arrayLength;
193
194 //Aggressive early exit
195 if(maxScoreThisPass<minScore){return;}
196
197 //Safe early exit
198 if(maxScoreThisPass<minViableScore){return;}
199
200 int[] temp=prev;
201 prev=next;
202 next=temp;
203 }
204
205 maxScore=-32000;
206 for(int apos=1; apos<arrayLength2; apos++){//Grab high score from last iteration
207 int score=prev[apos];
208 if(score>=maxScore){
209 maxScore=score;
210 maxQpos=qlen;
211 maxRpos=apos-1;
212 }
213 }
214 retArray[0]=maxScore;
215 retArray[1]=maxQpos;
216 retArray[2]=maxRpos;
217 retArray[3]=iters;
218 }
219
220 private static void alignForwardShortPseudo(int[] query, int[] ref, int[] retArray, int qlen, int rlen) {
221
222 final int arrayLength=qlen;
223 final int arrayLength2=qlen+1;
224
225 //Stack allocated; faster.
226 int[] array1=new int[arrayLength2];
227 int[] array2=new int[arrayLength2];
228
229 int[] prev=array1;
230 int[] next=array2;
231 int maxScore=-999999;
232 int maxQpos=-1;
233 int maxRpos=-1;
234 int itersShort=0;
235
236 for(int i=0; i<arrayLength2; i++) {
237 prev[i]=pointsIns*i;
238 next[i]=0;//For C version initialization
239 }
240
241 for(int rpos=0; rpos<rlen; rpos++){
242 if(rlen-rpos<arrayLength){prev[0]=next[0]+pointsDel;}
243
244 final int r=ref[rpos];
245 int score=0;
246
247 // //Inner loop
248 // for(int qpos=0, apos=1; qpos<arrayLength; qpos++, apos++){
249 // final int q=query[qpos];
250 // final boolean match=(q==r);
251 // final int vScore=prev[apos]+pointsIns;
252 // final int hScore=next[apos-1]+pointsDel;
253 // final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
254 //
255 // score=(dScore>=vScore ? dScore : vScore);
256 // score=(hScore>score ? hScore : score);
257 //
258 // next[apos]=score;
259 // }
260
261 //Inner DV loop
262 for(int qpos=0, apos=1; qpos<arrayLength; qpos++, apos++){
263 final int q=query[qpos];
264 final boolean match=(q==r);
265 final int vScore=prev[apos]+pointsIns;
266 final int dScore=(match ? pointsMatch : pointsSub)+prev[apos-1];
267
268 score=(dScore>=vScore ? dScore : vScore);
269
270 next[apos]=score;
271 }
272
273 //Inner I loop
274 for(int qpos=0, apos=1; qpos<arrayLength; qpos++, apos++){
275 final int hScore=next[apos-1]+pointsDel;
276
277 score=next[apos];
278 score=(hScore>score ? hScore : score);
279
280 next[apos]=score;
281 }
282
283 itersShort+=arrayLength;
284
285 if(score>=maxScore){
286 maxScore=score;
287 maxQpos=arrayLength-1;
288 maxRpos=rpos;
289 }
290
291 int[] temp=prev;
292 prev=next;
293 next=temp;
294 }
295 retArray[0]=maxScore;
296 retArray[1]=maxQpos;
297 retArray[2]=maxRpos;
298 retArray[3]=itersShort;
299 }
300
301 private static native void alignForwardJNI(int[] query, int[] ref, int[] retArray, int qlen, int rlen, int minScore, float minRatio);
302 private static native void alignForward16JNI(short[] query, short[] ref, int[] retArray, short qlen, short rlen, short minScore, float minRatio);
303 private static native void alignForwardShortJNI(int[] query, int[] ref, int[] retArray, int qlen, int rlen);
304 private static native void alignForwardShort16JNI(short[] query, short[] ref, int[] retArray, short qlen, short rlen);
305
306 /*--------------------------------------------------------------*/
307 /*---------------- Getters ----------------*/
308 /*--------------------------------------------------------------*/
309
310 @Override
311 long iters(){return iters;}
312
313 @Override
314 long itersShort(){return itersShort;}
315
316 /*--------------------------------------------------------------*/
317 /*---------------- Fields ----------------*/
318 /*--------------------------------------------------------------*/
319
320 long iters = 0;
321 long itersShort = 0;
322
323 /*--------------------------------------------------------------*/
324 /*---------------- Constants ----------------*/
325 /*--------------------------------------------------------------*/
326
327 public static final int pointsMatch = 1;
328 public static final int pointsSub = -1;
329 public static final int pointsDel = -2;
330 public static final int pointsIns = -2;
331
332 }