1 #include <U2Core/disable-warnings.h>
2 U2_DISABLE_WARNINGS
3
4 #include "muscle.h"
5 #include "muscle_context.h"
6
7 extern SCOREMATRIX VTML_LA;
8 extern SCOREMATRIX PAM200;
9 extern SCOREMATRIX PAM200NoCenter;
10 extern SCOREMATRIX VTML_SP;
11 extern SCOREMATRIX VTML_SPNoCenter;
12 extern SCOREMATRIX NUC_SP;
13
InitializePPTerm()14 bool MuscleContext::glbalndimer_stuct::InitializePPTerm()
15 {
16 PPTerm.m_bAllGaps = false;
17 PPTerm.m_LL = 1;
18 PPTerm.m_LG = 0;
19 PPTerm.m_GL = 0;
20 PPTerm.m_GG = 0;
21 PPTerm.m_fOcc = 1;
22 return true;
23 }
24
MuscleContext(int _nThreads)25 MuscleContext::MuscleContext(int _nThreads)
26 {
27 memset(&intmath, 0, sizeof(intmath));
28 memset(¶ms, 0, sizeof(params));
29 memset(&alpha, 0, sizeof(alpha));
30 memset(&muscle, 0, sizeof(muscle));
31 memset(&nwdasimple, 0, sizeof(nwdasimple));
32 memset(&ppscore, 0, sizeof(ppscore));
33 memset(&scoredist, 0, sizeof(scoredist));
34 memset(&objscore2, 0, sizeof(objscore2));
35 memset(&nwsmalls, 0, sizeof(nwsmalls));
36 memset(&fastdistnuc, 0, sizeof(fastdistnuc));
37 memset(&fastdistmafft, 0, sizeof(fastdistmafft));
38 memset(&glbalignsp, 0, sizeof(glbalignsp));
39 memset(&glbalignspn, 0, sizeof(glbalignspn));
40 memset(&glbalignss, 0, sizeof(glbalignss));
41 memset(&glbalignle, 0, sizeof(glbalignle));
42 memset(&setnewhandler, 0, sizeof(setnewhandler));
43 memset(&options, 0, sizeof(options));
44 memset(&globals, 0, sizeof(globals));
45 memset(&progress, 0, sizeof(progress));
46 memset(&scoregaps, 0, sizeof(scoregaps));
47 memset(&glbaligndiag, 0, sizeof(glbaligndiag));
48 memset(&upgma2, 0, sizeof(upgma2));
49 memset(&msa2, 0, sizeof(msa2));
50 memset(&validateids, 0, sizeof(validateids));
51 memset(&refinehoriz, 0, sizeof(refinehoriz));
52 memset(&spfast, 0, sizeof(spfast));
53 memset(&readmx, 0, sizeof(readmx));
54 memset(&mhack, 0, sizeof(mhack));
55 memset(&glbalndimer, 0, sizeof(glbalndimer));
56 memset(&subfams, 0, sizeof(subfams));
57 memset(&savebest, 0, sizeof(savebest));
58 memset(&enumtostr, 0, sizeof(enumtostr));
59 memset(&finddiags, 0, sizeof(finddiags));
60 memset(&finddiagsn, 0, sizeof(finddiagsn));
61 memset(&globalswin32, 0, sizeof(globalswin32));
62
63 nThreads = _nThreads;
64 progressStub = 0;
65 progressPercent = &progressStub;
66 cancelStub = 0;
67 cancelFlag = &cancelStub;
68
69 m_uIdCount = 0;
70
71 //intmath.cpp
72 intmath.bInit = false;
73
74 //params
75 params.g_scoreCenter = 0;
76 params.g_scoreGapExtend = 0;
77 params.g_scoreGapOpen2 = MINUS_INFINITY;
78 params.g_scoreGapExtend2 = MINUS_INFINITY;
79 params.g_scoreGapAmbig = 0;
80 params.g_scoreAmbigFactor = 0;
81
82 params.g_pstrInFileName = "-";
83 params.g_pstrOutFileName = "-";
84 params.g_pstrFASTAOutFileName = 0;
85 params.g_pstrMSFOutFileName = 0;
86 params.g_pstrClwOutFileName = 0;
87 params.g_pstrClwStrictOutFileName = 0;
88 params.g_pstrHTMLOutFileName = 0;
89 params.g_pstrPHYIOutFileName = 0;
90 params.g_pstrPHYSOutFileName = 0;
91 params.g_pstrDistMxFileName1 = 0;
92 params.g_pstrDistMxFileName2 = 0;
93
94 params.g_pstrFileName1 = 0;
95 params.g_pstrFileName2 = 0;
96
97 params.g_pstrSPFileName = 0;
98 params.g_pstrMatrixFileName = 0;
99
100 params.g_pstrUseTreeFileName = 0;
101
102 ////
103 params.g_pstrInFileName = "-";
104 params.g_pstrOutFileName = "-";
105 params.g_pstrFASTAOutFileName = 0;
106 params.g_pstrMSFOutFileName = 0;
107 params.g_pstrClwOutFileName = 0;
108 params.g_pstrClwStrictOutFileName = 0;
109 params.g_pstrHTMLOutFileName = 0;
110 params.g_pstrPHYIOutFileName = 0;
111 params.g_pstrPHYSOutFileName = 0;
112 params.g_pstrDistMxFileName1 = 0;
113 params.g_pstrDistMxFileName2 = 0;
114
115 params.g_pstrComputeWeightsFileName = 0;
116 params.g_pstrFileName1 = 0;
117 params.g_pstrFileName2 = 0;
118
119 params.g_pstrSPFileName = 0;
120 params.g_pstrMatrixFileName = 0;
121
122 params.g_pstrUseTreeFileName = 0;
123
124 params.g_pstrComputeWeightsFileName = NULL;
125
126 params.g_pstrProf1FileName = 0;
127 params.g_pstrProf2FileName = 0;
128
129 params.g_bUseTreeNoWarn = false;
130
131 params.g_uSmoothWindowLength = 7;
132 params.g_uAnchorSpacing = 32;
133 params.g_uMaxTreeRefineIters = 1;
134
135 params.g_uRefineWindow = 200;
136 params.g_uWindowFrom = 0;
137 params.g_uWindowTo = 0;
138 params.g_uSaveWindow = uInsane;
139 params.g_uWindowOffset = 0;
140
141 params.g_uMaxSubFamCount = 5;
142
143 params.g_uHydrophobicRunLength = 5;
144 params.g_dHydroFactor = (float) 1.2;
145
146 params.g_uMinDiagLength = 24; // TODO alpha -- should depend on alphabet?
147 params.g_uMaxDiagBreak = 1;
148 params.g_uDiagMargin = 5;
149
150 params.g_dSUEFF = (float) 0.1;
151
152 params.g_bPrecompiledCenter = true;
153 params.g_bNormalizeCounts = false;
154 params.g_bDiags1 = false;
155 params.g_bDiags2 = false;
156 params.g_bAnchors = true;
157 params.g_bQuiet = false;
158 params.g_bVerbose = false;
159 params.g_bRefine = false;
160 params.g_bRefineW = false;
161 params.g_bProfDB = false;
162 params.g_bLow = false;
163 params.g_bSW = false;
164 params.g_bClusterOnly = false;
165 params.g_bProfile = false;
166 params.g_bPPScore = false;
167 params.g_bBrenner = false;
168 params.g_bDimer = false;
169 params.g_bVersion = false;
170 params.g_bStable = false;
171 params.g_bFASTA = false;
172 params.g_bPAS = false;
173 params.g_bTomHydro = false;
174 params.g_bMakeTree = false;
175
176 #if DEBUG
177 params.g_bCatchExceptions = false;
178 #else
179 params.g_bCatchExceptions = true;
180 #endif
181
182 params.g_bMSF = false;
183 params.g_bAln = false;
184 params.g_bClwStrict = false;
185 params.g_bHTML = false;
186 params.g_bPHYI = false;
187 params.g_bPHYS = false;
188
189 params.g_uMaxIters = 8;
190 params.g_ulMaxSecs = 0;
191 params.g_uMaxMB = 500;
192
193 params.g_PPScore = PPSCORE_LE;
194 params.g_ObjScore = OBJSCORE_SPM;
195
196 params.g_SeqWeight1 = SEQWEIGHT_ClustalW;
197 params.g_SeqWeight2 = SEQWEIGHT_ClustalW;
198
199 params.g_Distance1 = DISTANCE_Kmer6_6;
200 params.g_Distance2 = DISTANCE_PctIdKimura;
201
202 params.g_Cluster1 = CLUSTER_UPGMB;
203 params.g_Cluster2 = CLUSTER_UPGMB;
204
205 params.g_Root1 = ROOT_Pseudo;
206 params.g_Root2 = ROOT_Pseudo;
207
208 // params.g_bDiags;
209
210 params.g_SeqType = SEQTYPE_Auto;
211
212 params.g_TermGaps = TERMGAPS_Half;
213
214 //------------------------------------------------------
215 // These parameters depending on the chosen prof-prof
216 // (params.g_PP), initialized to "Undefined".
217 params.g_dSmoothScoreCeil = fInsane;
218 params.g_dMinBestColScore = fInsane;
219 params.g_dMinSmoothScore = fInsane;
220 params.g_scoreGapOpen = fInsane;
221
222 //alpha.cpp
223 alpha.g_Alpha = ALPHA_Undefined;
224 alpha.g_AlphaSize = 0;
225 alpha.InvalidLetterCount = 0;
226
227 //ppscore.cpp
228 ppscore.g_bTracePPScore = false;
229 ppscore.g_ptrPPScoreMSA1 = 0;
230 ppscore.g_ptrPPScoreMSA2 = 0;
231
232 //setnewhandler.cpp
233 setnewhandler.EmergencyReserve = 0;
234
235 //options.cpp
236 static options_struct::VALUE_OPT ValueOptsSample[] =
237 {
238 {"in", 0},
239 {"in1", 0},
240 {"in2", 0},
241 {"out", 0},
242 {"MaxIters", 0},
243 {"MaxHours", 0},
244 {"GapOpen", 0},
245 {"GapOpen2", 0},
246 {"GapExtend", 0},
247 {"GapExtend2", 0},
248 {"GapAmbig", 0},
249 {"Center", 0},
250 {"SmoothScoreCeil", 0},
251 {"MinBestColScore", 0},
252 {"MinSmoothScore", 0},
253 {"ObjScore", 0},
254 {"SmoothWindow", 0},
255 {"RefineWindow", 0},
256 {"FromWindow", 0},
257 {"ToWindow", 0},
258 {"SaveWindow", 0},
259 {"WindowOffset", 0},
260 {"FirstWindow", 0},
261 {"AnchorSpacing", 0},
262 {"Log", 0},
263 {"LogA", 0},
264 {"MaxTrees", 0},
265 {"SUEFF", 0},
266 {"Distance", 0},
267 {"Distance1", 0},
268 {"Distance2", 0},
269 {"Weight", 0},
270 {"Weight1", 0},
271 {"Weight2", 0},
272 {"Cluster", 0},
273 {"Cluster1", 0},
274 {"Cluster2", 0},
275 {"Root1", 0},
276 {"Root2", 0},
277 {"Tree1", 0},
278 {"Tree2", 0},
279 {"UseTree", 0},
280 {"UseTree_NoWarn", 0},
281 {"DiagLength", 0},
282 {"DiagMargin", 0},
283 {"DiagBreak", 0},
284 {"Hydro", 0},
285 {"HydroFactor", 0},
286 {"SPScore", 0},
287 {"SeqType", 0},
288 {"MaxMB", 0},
289 {"ComputeWeights", 0},
290 {"MaxSubFam", 0},
291 {"ScoreFile", 0},
292 {"TermGaps", 0},
293 {"FASTAOut", 0},
294 {"CLWOut", 0},
295 {"CLWStrictOut", 0},
296 {"HTMLOut", 0},
297 {"MSFOut", 0},
298 {"PHYIOut", 0},
299 {"PHYSOut", 0},
300 {"Matrix", 0},
301 {"DistMx1", 0},
302 {"DistMx2", 0},
303 {"Weight", 0},
304 };
305 options.ValueOptCount = sizeof(ValueOptsSample)/sizeof(ValueOptsSample[0]);
306 options.ValueOpts = new options_struct::VALUE_OPT[options.ValueOptCount];
307 memcpy(options.ValueOpts,ValueOptsSample,sizeof(ValueOptsSample));
308
309 static options_struct::FLAG_OPT FlagOptsSample[] =
310 {
311 {"LE", false},
312 {"SP", false},
313 {"SV", false},
314 {"SPN", false},
315 {"Core", false},
316 {"NoCore", false},
317 {"Diags1", false},
318 {"Diags2", false},
319 {"Diags", false},
320 {"Quiet", false},
321 {"MSF", false},
322 {"Verbose", false},
323 {"Anchors", false},
324 {"NoAnchors", false},
325 {"Refine", false},
326 {"RefineW", false},
327 {"SW", false},
328 {"Profile", false},
329 {"PPScore", false},
330 {"ClusterOnly", false},
331 {"Brenner", false},
332 {"Dimer", false},
333 {"clw", false},
334 {"clwstrict", false},
335 {"HTML", false},
336 {"Version", false},
337 {"Stable", false},
338 {"Group", false},
339 {"FASTA", false},
340 {"ProfDB", false},
341 {"PAS", false},
342 {"PHYI", false},
343 {"PHYS", false},
344 {"TomHydro", false},
345 {"MakeTree", false},
346 };
347 options.FlagOptCount = sizeof(FlagOptsSample)/sizeof(FlagOptsSample[0]);
348 options.FlagOpts = new options_struct::FLAG_OPT[options.FlagOptCount];
349 memcpy(options.FlagOpts,FlagOptsSample,sizeof(FlagOptsSample));
350
351 //globals.cpp
352 globals.g_bListFileAppend = false;
353 globals.g_SeqWeight = SEQWEIGHT_Undefined;
354
355 //progress.cpp
356 progress.g_uIter = 0; // Main MUSCLE iteration 1, 2..
357 progress.g_uLocalMaxIters = 0; // Max iters
358 progress.g_fProgress = stderr; // Default to standard error
359 progress.pr_printf = fprintf;
360 progress.g_bWipeDesc = false;
361
362 //glbaligndiag.cpp
363 glbaligndiag.g_dDPAreaWithoutDiags = 0.0;
364 glbaligndiag.g_dDPAreaWithDiags = 0.0;
365
366 //glbalndimer.cpp
367 glbalndimer.PPTermInitialized = glbalndimer.InitializePPTerm();
368
369 //validateids.cpp
370 #if _DEBUG
371 validateids.g_ptrMuscleSeqVect = 0;
372 #endif
373 validateids.g_uTreeSplitNode1 = NULL_NEIGHBOR;
374 validateids.g_uTreeSplitNode2 = NULL_NEIGHBOR;
375
376 //msa2.cpp
377 msa2.g_ptrMuscleTree = 0;
378 msa2.g_uTreeSplitNode1 = NULL_NEIGHBOR;
379 msa2.g_uTreeSplitNode2 = NULL_NEIGHBOR;
380 msa2.g_MuscleWeights = NULL;
381
382 //mhack.cpp
383 mhack.M = NULL;
384
385 //nwsmall->cpp
386 nwsmalls = new nwsmall_struct[nThreads];
387 for(int i=0;i<nThreads;i++) {
388 nwsmall_struct* nwsmall = &nwsmalls[i];
389 nwsmall->uCachePrefixCountA = 0;
390 nwsmall->uCachePrefixCountB = 0;
391 nwsmall->CacheTB = NULL;
392 nwsmall->CacheDRow = NULL;
393 nwsmall->CacheMCurr = NULL;
394 nwsmall->CacheMNext = NULL;
395 nwsmall->CacheMPrev = NULL;
396 }
397
398 //readmx.cpp
399 readmx.HeadingCount = 0;
400
401 //spfast.cpp
402 spfast.bGapScoreMatrixInit = false;
403 }
404
fillUidsVectors(int rowsCount)405 void MuscleContext::fillUidsVectors(int rowsCount) {
406 const int initialRowsCount = input_uIds.length();
407 tmp_uIds.reserve(initialRowsCount + rowsCount);
408 input_uIds.reserve(initialRowsCount + rowsCount);
409 for (unsigned i = initialRowsCount; i < initialRowsCount + rowsCount; i++) {
410 input_uIds.append(i);
411 tmp_uIds.append(i);
412 }
413 }
414
415 extern void FreeDPMemSPN();
416
~MuscleContext()417 MuscleContext::~MuscleContext() {
418 //params.g_ptrScoreMatrix - do not delete -> inited from static structures
419
420 delete ppscore.g_ptrPPScoreMSA1;
421 delete ppscore.g_ptrPPScoreMSA2;
422 //glbalignsp.cpp
423 if (glbalignsp.DPM.uLength > 0) {
424 for (unsigned i = 0; i < glbalignsp.DPM.uLength; ++i) {
425 delete[] glbalignsp.DPM.TraceBack[i];
426 delete[] glbalignsp.DPM.FreqsA[i];
427 delete[] glbalignsp.DPM.SortOrderA[i];
428 }
429 for (unsigned n = 0; n < 20; ++n) {
430 delete[] glbalignsp.DPM.ScoreMxB[n];
431 }
432
433 delete[] glbalignsp.DPM.MPrev;
434 delete[] glbalignsp.DPM.MCurr;
435 delete[] glbalignsp.DPM.MWork;
436 delete[] glbalignsp.DPM.DPrev;
437 delete[] glbalignsp.DPM.DCurr;
438 delete[] glbalignsp.DPM.DWork;
439 delete[] glbalignsp.DPM.uDeletePos;
440 delete[] glbalignsp.DPM.GapOpenA;
441 delete[] glbalignsp.DPM.GapOpenB;
442 delete[] glbalignsp.DPM.GapCloseA;
443 delete[] glbalignsp.DPM.GapCloseB;
444 delete[] glbalignsp.DPM.SortOrderA;
445 delete[] glbalignsp.DPM.FreqsA;
446 delete[] glbalignsp.DPM.ScoreMxB;
447 delete[] glbalignsp.DPM.TraceBack;
448 }
449
450
451 //glbalignspn.cpp
452 //FreeDPMemSPN();
453 const unsigned uOldLength = glbalignspn.DPM.uLength;
454 if (0 != uOldLength)
455 {
456 for (unsigned i = 0; i < uOldLength; ++i)
457 {
458 delete[] glbalignspn.DPM.TraceBack[i];
459 delete[] glbalignspn.DPM.FreqsA[i];
460 delete[] glbalignspn.DPM.SortOrderA[i];
461 }
462 for (unsigned n = 0; n < 4; ++n)
463 delete[] glbalignspn.DPM.ScoreMxB[n];
464
465 delete[] glbalignspn.DPM.MPrev;
466 delete[] glbalignspn.DPM.MCurr;
467 delete[] glbalignspn.DPM.MWork;
468 delete[] glbalignspn.DPM.DPrev;
469 delete[] glbalignspn.DPM.DCurr;
470 delete[] glbalignspn.DPM.DWork;
471 delete[] glbalignspn.DPM.uDeletePos;
472 delete[] glbalignspn.DPM.GapOpenA;
473 delete[] glbalignspn.DPM.GapOpenB;
474 delete[] glbalignspn.DPM.GapCloseA;
475 delete[] glbalignspn.DPM.GapCloseB;
476 delete[] glbalignspn.DPM.SortOrderA;
477 delete[] glbalignspn.DPM.FreqsA;
478 delete[] glbalignspn.DPM.ScoreMxB;
479 delete[] glbalignspn.DPM.TraceBack;
480 }
481 //glbalignle.cpp
482 if(glbalignle.DPM.uLength > 0) {
483 for (unsigned i = 0; i < glbalignle.DPM.uLength; ++i) {
484 delete[] glbalignle.DPM.TraceBack[i];
485 delete[] glbalignle.DPM.FreqsA[i];
486 delete[] glbalignle.DPM.SortOrderA[i];
487 }
488 for (unsigned n = 0; n < 20; ++n) {
489 delete[] glbalignle.DPM.ScoreMxB[n];
490 }
491
492 delete[] glbalignle.DPM.MPrev;
493 delete[] glbalignle.DPM.MCurr;
494 delete[] glbalignle.DPM.MWork;
495 delete[] glbalignle.DPM.DPrev;
496 delete[] glbalignle.DPM.DCurr;
497 delete[] glbalignle.DPM.DWork;
498 delete[] glbalignle.DPM.uDeletePos;
499 delete[] glbalignle.DPM.GapOpenA;
500 delete[] glbalignle.DPM.GapOpenB;
501 delete[] glbalignle.DPM.GapCloseA;
502 delete[] glbalignle.DPM.GapCloseB;
503 delete[] glbalignle.DPM.SortOrderA;
504 delete[] glbalignle.DPM.FreqsA;
505 delete[] glbalignle.DPM.ScoreMxB;
506 delete[] glbalignle.DPM.TraceBack;
507 #if OCC
508 delete[] glbalignle.DPM.OccA;
509 delete[] glbalignle.DPM.OccB;
510 #endif
511 }
512
513 //glbalignss.cpp
514 if (glbalignss.DPM.uLength > 0) {
515 for (unsigned i = 0; i < glbalignss.DPM.uLength; ++i) {
516 delete[] glbalignss.DPM.TraceBack[i];
517 }
518
519 delete[] glbalignss.DPM.MPrev;
520 delete[] glbalignss.DPM.MCurr;
521 delete[] glbalignss.DPM.MWork;
522 delete[] glbalignss.DPM.DPrev;
523 delete[] glbalignss.DPM.DCurr;
524 delete[] glbalignss.DPM.DWork;
525 delete[] glbalignss.DPM.MxRowA;
526 delete[] glbalignss.DPM.LettersB;
527 delete[] glbalignss.DPM.uDeletePos;
528 delete[] glbalignss.DPM.TraceBack;
529 }
530
531 //setnewhandler.cpp
532 free(setnewhandler.EmergencyReserve);
533
534 //scoregaps.cpp
535 scoregaps_struct::GAPINFO *GI = scoregaps.g_FreeList;
536 while(GI!=NULL) {
537 scoregaps_struct::GAPINFO *GI_del = GI;
538 GI=GI->Next;
539 delete GI_del;
540 }
541 delete[] scoregaps.g_Gaps;
542 delete[] scoregaps.g_ColDiff;
543
544 //msa2.cpp
545 delete[] msa2.g_MuscleWeights;
546 //created on stack only: delete msa2.g_ptrMuscleTree;
547
548 //mhack.cpp
549 delete[] mhack.M;
550
551 //nwsmall->cpp
552 for(int iThread=0;iThread<nThreads;iThread++) {
553 nwsmall_struct* nwsmall = &nwsmalls[iThread];
554 for (unsigned i = 0; i < nwsmall->uCachePrefixCountA; ++i) {
555 delete[] nwsmall->CacheTB[i];
556 }
557 delete[] nwsmall->CacheTB;
558 delete[] nwsmall->CacheDRow;
559 delete[] nwsmall->CacheMCurr;
560 delete[] nwsmall->CacheMNext;
561 delete[] nwsmall->CacheMPrev;
562 }
563 delete[] nwsmalls;
564
565 //options.cpp
566 for (int i=0; i<options.ValueOptCount; i++) {
567 delete options.ValueOpts[i].m_pstrValue;
568 }
569 delete[] options.ValueOpts;
570 delete[] options.FlagOpts;
571
572 //subfams.cpp
573 delete[] subfams.ShortestPathEstimate;
574 delete[] subfams.Predecessor;
575
576 //savebest.cpp
577 /* created on stack only
578 if(savebest.ptrBestMSA!=NULL) {
579 savebest.ptrBestMSA->Clear();
580 delete savebest.ptrBestMSA;
581 }*/
582
583 //validateids.cpp
584 #if _DEBUG
585 //created on stack only delete validateids.g_ptrMuscleSeqVect;
586 #endif
587
588 }
589
MuscleException(const char * _str)590 MuscleException::MuscleException(const char* _str) {
591 int len = strlen(_str);
592 assert(len < 4096);
593 memcpy(str, _str, len);
594 str[len] = '\0';
595 }
596
MuscleException()597 MuscleException::MuscleException() {
598 memset(str, '\0', 4096);
599 }
600