1 /******************************************************************
\r
2 Copyright 2006 by Michael Farrar. All rights reserved.
\r
3 This program may not be sold or incorporated into a commercial product,
\r
4 in whole or in part, without written consent of Michael Farrar. For
\r
5 further information regarding permission for use or reproduction, please
\r
6 contact: Michael Farrar at farrar.michael@gmail.com.
\r
7 *******************************************************************/
\r
10 Written by Michael Farrar, 2006.
\r
11 Please send bug reports and/or suggestions to farrar.michael@gmail.com.
\r
18 #include "swscalar.h"
\r
28 swScalar (unsigned char *querySeq,
\r
30 unsigned char *dbSeq,
\r
32 unsigned short gapOpen,
\r
33 unsigned short gapExtend,
\r
39 swScalarInit(unsigned char *querySeq,
\r
41 signed char *matrix)
\r
48 SwScalarData *pSwData;
\r
50 /* remove unreferenced warnings */
\r
53 pSwData = (SwScalarData *) malloc (sizeof (SwScalarData));
\r
55 fprintf (stderr, "Unable to allocate memory for SW data\n");
\r
59 nCount = ALPHA_SIZE * ALPHA_SIZE + /* matrix */
\r
60 (queryLength * 3); /* vH1, vH2 and vE */
\r
62 pSwData->pData = (int *) calloc (nCount, sizeof (int));
\r
63 if (!pSwData->pData) {
\r
64 fprintf (stderr, "Unable to allocate memory for SW data buffers\n");
\r
68 pSwData->pMatrix = pSwData->pData;
\r
70 pSwData->pH = pSwData->pMatrix + ALPHA_SIZE * ALPHA_SIZE;
\r
71 pSwData->pE = pSwData->pH + queryLength;
\r
73 /* Conver the scoring matrix to ints */
\r
74 pi = pSwData->pMatrix;
\r
75 for (i = 0; i < ALPHA_SIZE; ++i) {
\r
76 for (j = 0; j < ALPHA_SIZE; ++j) {
\r
77 *pi++ = (int) *matrix++;
\r
84 void swScalarScan (unsigned char *querySeq,
\r
88 SEARCH_OPTIONS *options,
\r
93 int threshold = options->threshold;
\r
95 unsigned char *dbSeq;
\r
98 int gapInit = -(options->gapInit + options->gapExt);
\r
99 int gapExt = -options->gapExt;
\r
101 SwScalarData *scalarData = (SwScalarData *) swData;
\r
103 dbSeq = nextSeq (dbLib, &dbLen);
\r
104 while (dbLen > 0) {
\r
106 score = swScalar (querySeq, queryLength,
\r
109 scalarData->pMatrix,
\r
113 if (score >= threshold) {
\r
114 int minScore = insertList (scores, score, seqName (dbLib));
\r
115 if (minScore >= threshold) {
\r
116 threshold = minScore;
\r
120 dbSeq = nextSeq (dbLib, &dbLen);
\r
125 swScalarComplete(void *pSwData)
\r
127 SwScalarData *pScalarData = (SwScalarData *) pSwData;
\r
129 free (pScalarData->pData);
\r
130 free (pScalarData);
\r
134 swScalar(unsigned char *querySeq,
\r
136 unsigned char *dbSeq,
\r
138 unsigned short gapOpen,
\r
139 unsigned short gapExtend,
\r
153 /* Zero out the storage vector */
\r
154 for (i = 0; i < queryLength; i++)
\r
160 for (i = 0; i < dbLength; ++i)
\r
162 /* fetch first data asap. */
\r
163 pScore = pMatrix + dbSeq[i] * ALPHA_SIZE;
\r
168 /* load the next h value */
\r
171 for (j = 0; j < queryLength; j++)
\r
173 /* load E values */
\r
176 /* add score to H */
\r
177 H = prevH + *(pScore + querySeq[j]);
\r
182 /* Update highest score encountered this far */
\r
183 if (H > maxScore) {
\r
187 /* get max from H, E and F */
\r
195 /* save H values */
\r
201 /* update E value */
\r
210 /* update vF value */
\r
219 /* save vE values */
\r
224 /* return largest score */
\r