29 using namespace shogun;
51 -1, 0, 0, 2, -5, 3, 7,
52 0, -3, -1, -2, -4, -3, -3, 8,
53 -2, 0, 1, -2, -4, 1, 0, -3, 11,
54 -2, -5, -5, -5, -2, -4, -5, -6, -5, 6,
55 -2, -3, -5, -5, -2, -3, -4, -5, -4, 2, 6,
56 -1, 3, 0, -1, -5, 2, 1, -2, -1, -4, -4, 7,
57 -1, -2, -3, -5, -2, -1, -3, -4, -2, 2, 3, -2, 8,
58 -3, -4, -5, -5, -4, -5, -5, -5, -2, 0, 1, -5, 0, 9,
59 -1, -3, -3, -2, -4, -2, -2, -3, -3, -4, -4, -2, -4, -5, 11,
60 2, -1, 1, 0, -1, 0, 0, 0, -1, -4, -4, 0, -2, -4, -1, 6,
61 0, -2, 0, -2, -1, -1, -1, -2, -3, -1, -2, -1, -1, -3, -2, 2, 7,
62 -4, -4, -6, -6, -3, -3, -4, -4, -4, -4, -2, -4, -2, 1, -6, -4, -4, 16,
63 -3, -3, -3, -5, -4, -2, -3, -5, 3, -2, -2, -3, -1, 4, -4, -3, -2, 3, 10,
64 0, -4, -4, -5, -1, -3, -4, -5, -5, 4, 1, -3, 1, -1, -4, -2, 0, -4, -2, 6};
67 static int32_t
BINDEX(int32_t i, int32_t j)
69 return (((i)>(j))?(j)+(((i)*(i+1))/2):(i)+(((j)*(j+1))/2));
84 #define LOGP(x,y) LogSum(x,y)
102 init_static_variables();
114 init_static_variables();
146 void CLocalAlignmentStringKernel::init_logsum(){
153 int32_t CLocalAlignmentStringKernel::LogSum(int32_t p1, int32_t p2)
156 static int32_t firsttime=1;
165 if (diff>=LOGSUM_TBL)
167 else if (diff<=-LOGSUM_TBL)
172 return p2+logsum_lookup[-diff];
179 return (p1-p2>50.) ? p1 : p1+log(1.+exp(p2-p1));
181 return (p2-p1>50.) ? p2 : p2+log(1.+exp(p1-p2));
185 void CLocalAlignmentStringKernel::init_static_variables()
191 if ((
aaIndex=(int32_t *)calloc(
NLET,
sizeof(int32_t)))==NULL)
197 if ((
isAA=(int32_t *)calloc(256,
sizeof(int32_t)))==NULL)
203 for (i=0 ; i<NAA*(NAA+1)/2; i++)
217 float64_t CLocalAlignmentStringKernel::LAkernelcompute(
218 int32_t* aaX, int32_t* aaY,
219 int32_t nX, int32_t nY )
259 for (j=0; j<cl; j++) {
277 for (i=1; i<=nX; i++) {
288 for (j=1; j<=nY; j++) {
302 logX2[curpos]=
LOGP(logM[frompos], logX2[frompos]);
314 aux=
LOGP(logM[frompos], logY2[frompos]);
315 logY2[curpos]=
LOGP(aux, logX2[frompos]);
322 aux=
LOGP(logX[frompos], logY[frompos]);
323 aux2=
LOGP(0, logM[frompos]);
344 aux=
LOGP(logX2[curpos], logY2[curpos]);
345 aux2=
LOGP(0, logM[curpos]);
377 if ( (lx<1) || (ly<1) )
382 if ((aax=(int32_t *)calloc(lx,
sizeof(int32_t)))==NULL)
384 if ((aay=(int32_t *)calloc(ly,
sizeof(int32_t)))==NULL)
391 if (
isAA[toupper(x[i])])
392 aax[j++]=
aaIndex[toupper(x[i])-
'A'];
396 if (
isAA[toupper(y[i])])
397 aay[j++]=
aaIndex[toupper(y[i])-
'A'];
402 float64_t result=LAkernelcompute(aax, aay, lx, ly);
414 void CLocalAlignmentStringKernel::init()