1 /* $Id: blast_options_builder.cpp 617780 2020-10-06 16:24:16Z gouriano $
2 * ===========================================================================
3 *
4 * PUBLIC DOMAIN NOTICE
5 * National Center for Biotechnology Information
6 *
7 * This software/database is a "United States Government Work" under the
8 * terms of the United States Copyright Act. It was written as part of
9 * the author's official duties as a United States Government employee and
10 * thus cannot be copyrighted. This software/database is freely available
11 * to the public for use. The National Library of Medicine and the U.S.
12 * Government have not placed any restriction on its use or reproduction.
13 *
14 * Although all reasonable efforts have been taken to ensure the accuracy
15 * and reliability of the software and data, the NLM and the U.S.
16 * Government do not and cannot warrant the performance or results that
17 * may be obtained by using this software or data. The NLM and the U.S.
18 * Government disclaim all warranties, express or implied, including
19 * warranties of performance, merchantability or fitness for any particular
20 * purpose.
21 *
22 * Please cite the author in any work or product based on this material.
23 *
24 * ===========================================================================
25 *
26 * Author: Kevin Bealer
27 *
28 * ===========================================================================
29 */
30
31 /// @file blast_options_builder.cpp
32 /// Defines the CBlastOptionsBuilder class
33
34 #include <ncbi_pch.hpp>
35 #include <algo/blast/api/remote_blast.hpp>
36 #include <algo/blast/api/blast_options_builder.hpp>
37
38 /** @addtogroup AlgoBlast
39 *
40 * @{
41 */
42
43 BEGIN_NCBI_SCOPE
44 USING_SCOPE(objects);
BEGIN_SCOPE(blast)45 BEGIN_SCOPE(blast)
46
47 CBlastOptionsBuilder::
48 CBlastOptionsBuilder(const string & program,
49 const string & service,
50 CBlastOptions::EAPILocality locality,
51 bool ignore_unsupported_options)
52 : m_Program (program),
53 m_Service (service),
54 m_PerformCulling (false),
55 m_HspRangeMax (0),
56 m_QueryRange (TSeqRange::GetEmpty()),
57 m_Locality (locality),
58 m_IgnoreUnsupportedOptions(ignore_unsupported_options),
59 m_ForceMbIndex (false)
60 {
61 }
62
63 EProgram
ComputeProgram(const string & program,const string & service)64 CBlastOptionsBuilder::ComputeProgram(const string & program,
65 const string & service)
66 {
67 string p = program;
68 string s = service;
69
70 NStr::ToLower(p);
71 NStr::ToLower(s);
72
73 // a. is there a program for phiblast?
74 // b. others, like disco?
75
76 bool found = false;
77
78 if (p == "blastp") {
79 if (s == "rpsblast") {
80 p = "rpsblast";
81 found = true;
82 } else if (s == "psi") {
83 p = "psiblast";
84 found = true;
85 } else if (s == "phi") {
86 // phi is just treated as a blastp here
87 found = true;
88 } else if (s == "delta_blast") {
89 p = "deltablast";
90 found = true;
91 }
92 else if (s == "indexed") {
93 found = true;
94 }
95 } else if (p == "blastn") {
96 if (s == "megablast") {
97 p = "megablast";
98 found = true;
99 }
100 if (s == "vecscreen") {
101 p = "vecscreen";
102 found = true;
103 }
104 if ((s == "sra") || (s == "wgs") || (s=="indexed")) {
105 found = true;
106 }
107 } else if (p == "tblastn") {
108 if (s == "rpsblast") {
109 p = "rpstblastn";
110 found = true;
111 } else if (s == "psi") {
112 p = "psitblastn";
113 found = true;
114 }
115 if ((s == "sra") || (s == "wgs")) {
116 found = true;
117 }
118 } else if (p == "tblastx") {
119 found = true;
120 } else if (p == "blastx") {
121 if (s == "rpsblast") {
122 p = "rpstblastn";
123 found = true;
124 }
125 }
126
127 if (!found && !(s=="plain" || s=="multi_blast")) {
128 string msg = "Unsupported combination of program (";
129 msg += program;
130 msg += ") and service (";
131 msg += service;
132 msg += ").";
133
134 NCBI_THROW(CBlastException, eInvalidArgument, msg);
135 }
136
137 return ProgramNameToEnum(p);
138 }
139
140 void CBlastOptionsBuilder::
x_ProcessOneOption(CBlastOptionsHandle & opts,objects::CBlast4_parameter & p)141 x_ProcessOneOption(CBlastOptionsHandle & opts,
142 objects::CBlast4_parameter & p)
143 {
144 const CBlast4_value & v = p.GetValue();
145
146 // Note that this code does not attempt to detect or repair
147 // inconsistencies; since this request has already been processed
148 // by SplitD, the results are assumed to be correct, for now.
149 // This will remain so unless options validation code becomes
150 // available, in which case it could be used by this code. This
151 // could be considered as a potential "to-do" item.
152
153 if (! p.CanGetName() || p.GetName().empty()) {
154 NCBI_THROW(CBlastException,
155 eInvalidArgument,
156 "Option has no name.");
157 }
158
159 string nm = p.GetName();
160
161 bool found = true;
162
163 // This switch is not really necessary. I wanted to break things
164 // up for human consumption. But as long as I'm doing that, I may
165 // as well use a performance-friendly paragraph marker.
166
167 CBlastOptions & bo = opts.SetOptions();
168
169 switch(nm[0]) {
170 case 'B':
171 if (CBlast4Field::Get(eBlastOpt_BestHitScoreEdge).Match(p)) {
172 bo.SetBestHitScoreEdge(v.GetReal());
173 } else if (CBlast4Field::Get(eBlastOpt_BestHitOverhang).Match(p)) {
174 bo.SetBestHitOverhang(v.GetReal());
175 } else {
176 found = false;
177 }
178 break;
179
180 case 'C':
181 if (CBlast4Field::Get(eBlastOpt_CompositionBasedStats).Match(p)) {
182 ECompoAdjustModes adjmode = (ECompoAdjustModes) v.GetInteger();
183 bo.SetCompositionBasedStats(adjmode);
184 } else if (CBlast4Field::Get(eBlastOpt_Culling).Match(p)) {
185 m_PerformCulling = v.GetBoolean();
186 } else if (CBlast4Field::Get(eBlastOpt_CullingLimit).Match(p)) {
187 bo.SetCullingLimit(v.GetInteger());
188 } else if (CBlast4Field::Get(eBlastOpt_CutoffScore).Match(p)) {
189 opts.SetCutoffScore(v.GetInteger());
190 } else {
191 found = false;
192 }
193 break;
194
195 case 'D':
196 if (CBlast4Field::Get(eBlastOpt_DbGeneticCode).Match(p)) {
197 bo.SetDbGeneticCode(v.GetInteger());
198 } else if (CBlast4Field::Get(eBlastOpt_DbLength).Match(p)) {
199 opts.SetDbLength(v.GetBig_integer());
200 } else if (CBlast4Field::Get(eBlastOpt_DustFiltering).Match(p)) {
201 bo.SetDustFiltering(v.GetBoolean());
202 } else if (CBlast4Field::Get(eBlastOpt_DustFilteringLevel).Match(p)) {
203 bo.SetDustFilteringLevel(v.GetInteger());
204 } else if (CBlast4Field::Get(eBlastOpt_DustFilteringWindow).Match(p)) {
205 bo.SetDustFilteringWindow(v.GetInteger());
206 } else if (CBlast4Field::Get(eBlastOpt_DustFilteringLinker).Match(p)) {
207 bo.SetDustFilteringLinker(v.GetInteger());
208 } else if (CBlast4Field::Get(eBlastOpt_DbFilteringAlgorithmId).Match(p)) {
209 m_DbFilteringAlgorithmId = v.GetInteger();
210 } else if (CBlast4Field::Get(eBlastOpt_DbFilteringAlgorithmKey).Match(p)) {
211 m_DbFilteringAlgorithmKey = v.GetString();
212 } else if (CBlast4Field::Get(eBlastOpt_DomainInclusionThreshold).Match(p)) {
213 bo.SetDomainInclusionThreshold(v.GetReal());
214 } else {
215 found = false;
216 }
217 break;
218
219 case 'E':
220 if (CBlast4Field::Get(eBlastOpt_EffectiveSearchSpace).Match(p)) {
221 opts.SetEffectiveSearchSpace(v.GetBig_integer());
222 } else if (CBlast4Field::Get(eBlastOpt_EntrezQuery).Match(p)) {
223 m_EntrezQuery = v.GetString();
224 } else if (CBlast4Field::Get(eBlastOpt_EvalueThreshold).Match(p)
225 || p.GetName() == "EvalueThreshold") {
226 if (v.IsReal()) {
227 opts.SetEvalueThreshold(v.GetReal());
228 } else if (v.IsCutoff() && v.GetCutoff().IsE_value()) {
229 opts.SetEvalueThreshold(v.GetCutoff().GetE_value());
230 } else {
231 string msg = "EvalueThreshold has unsupported type.";
232 NCBI_THROW(CBlastException, eInvalidArgument, msg);
233 }
234 } else {
235 found = false;
236 }
237 break;
238
239 case 'F':
240 if (CBlast4Field::Get(eBlastOpt_FilterString).Match(p)) {
241 opts.SetFilterString(v.GetString().c_str(), true); /* NCBI_FAKE_WARNING */
242 } else if (CBlast4Field::Get(eBlastOpt_FinalDbSeq).Match(p)) {
243 m_FinalDbSeq = v.GetInteger();
244 } else if (CBlast4Field::Get(eBlastOpt_FirstDbSeq).Match(p)) {
245 m_FirstDbSeq = v.GetInteger();
246 } else if (CBlast4Field::Get(eBlastOpt_ForceMbIndex).Match(p)) {
247 m_ForceMbIndex = v.GetBoolean();
248 } else {
249 found = false;
250 }
251 break;
252
253 case 'G':
254 if (CBlast4Field::Get(eBlastOpt_GapExtensionCost).Match(p)) {
255 bo.SetGapExtensionCost(v.GetInteger());
256 } else if (CBlast4Field::Get(eBlastOpt_GapOpeningCost).Match(p)) {
257 bo.SetGapOpeningCost(v.GetInteger());
258 } else if (CBlast4Field::Get(eBlastOpt_GiList).Match(p)) {
259 #ifdef NCBI_STRICT_GI
260 const list<int>& int_list = v.GetInteger_list();
261 list<TGi> gi_list;
262 ITERATE ( list<int>, it, int_list ) {
263 gi_list.push_back(GI_FROM(int, *it));
264 }
265 m_GiList = gi_list;
266 #elif defined(NCBI_INT8_GI)
267 if(v.IsBig_integer_list()) {
268 m_GiList = v.GetBig_integer_list();
269 }
270 else {
271 const list<int>& int_list = v.GetInteger_list();
272 list<TGi> gi_list;
273 ITERATE ( list<int>, it, int_list ) {
274 gi_list.push_back(GI_FROM(int, *it));
275 }
276 m_GiList = gi_list;
277 }
278 #else
279 m_GiList = v.GetInteger_list();
280 #endif
281 } else if (CBlast4Field::Get(eBlastOpt_GapTracebackAlgorithm).Match(p)) {
282 bo.SetGapTracebackAlgorithm((EBlastTbackExt) v.GetInteger());
283 } else if (CBlast4Field::Get(eBlastOpt_GapTrigger).Match(p)) {
284 bo.SetGapTrigger(v.GetReal());
285 } else if (CBlast4Field::Get(eBlastOpt_GapXDropoff).Match(p)) {
286 bo.SetGapXDropoff(v.GetReal());
287 } else if (CBlast4Field::Get(eBlastOpt_GapXDropoffFinal).Match(p)) {
288 bo.SetGapXDropoffFinal(v.GetReal());
289 } else if (CBlast4Field::Get(eBlastOpt_GapExtnAlgorithm).Match(p)) {
290 bo.SetGapExtnAlgorithm(static_cast<EBlastPrelimGapExt>
291 (v.GetInteger()));
292 } else {
293 found = false;
294 }
295 break;
296
297 case 'H':
298 if (CBlast4Field::Get(eBlastOpt_HitlistSize).Match(p)) {
299 opts.SetHitlistSize(v.GetInteger());
300 } else if (CBlast4Field::Get(eBlastOpt_HspRangeMax).Match(p)) {
301 m_HspRangeMax = v.GetInteger();
302 } else {
303 found = false;
304 }
305 break;
306
307 case 'I':
308 if (CBlast4Field::Get(eBlastOpt_InclusionThreshold).Match(p)) {
309 bo.SetInclusionThreshold(v.GetReal());
310 } else if (CBlast4Field::Get(eBlastOpt_IgnoreMsaMaster).Match(p)) {
311 bo.SetIgnoreMsaMaster(v.GetBoolean());
312 } else {
313 found = false;
314 }
315 break;
316
317 case 'L':
318 if (CBlast4Field::Get(eBlastOpt_LCaseMask).Match(p))
319 {
320 if (!m_IgnoreQueryMasks)
321 {
322 _ASSERT(v.IsQuery_mask());
323 CRef<CBlast4_mask> refMask(new CBlast4_mask);
324 refMask->Assign(v.GetQuery_mask());
325
326 if (!m_QueryMasks.Have())
327 {
328 TMaskList listEmpty;
329 m_QueryMasks = listEmpty;
330 }
331 m_QueryMasks.GetRef().push_back(refMask);
332 }
333 } else if (CBlast4Field::Get(eBlastOpt_LongestIntronLength).Match(p)) {
334 bo.SetLongestIntronLength(v.GetInteger());
335 } else {
336 found = false;
337 }
338 break;
339
340 case 'M':
341 if (CBlast4Field::Get(eBlastOpt_MBTemplateLength).Match(p)) {
342 bo.SetMBTemplateLength(v.GetInteger());
343 } else if (CBlast4Field::Get(eBlastOpt_MBTemplateType).Match(p)) {
344 bo.SetMBTemplateType(v.GetInteger());
345 } else if (CBlast4Field::Get(eBlastOpt_MatchReward).Match(p)) {
346 bo.SetMatchReward(v.GetInteger());
347 } else if (CBlast4Field::Get(eBlastOpt_MatrixName).Match(p)) {
348 bo.SetMatrixName(v.GetString().c_str());
349 } else if (CBlast4Field::Get(eBlastOpt_MatrixTable).Match(p)) {
350 // This is no longer used.
351 } else if (CBlast4Field::Get(eBlastOpt_MismatchPenalty).Match(p)) {
352 bo.SetMismatchPenalty(v.GetInteger());
353 } else if (CBlast4Field::Get(eBlastOpt_MaskAtHash).Match(p)) {
354 bo.SetMaskAtHash(v.GetBoolean());
355 } else if (CBlast4Field::Get(eBlastOpt_MbIndexName).Match(p)) {
356 m_MbIndexName = v.GetString();
357 } else if (CBlast4Field::Get(eBlastOpt_MaxHspsPerSubject).Match(p)) {
358 bo.SetMaxHspsPerSubject(v.GetInteger());
359 }
360 else {
361 found = false;
362 }
363 break;
364
365 case 'O':
366 if (CBlast4Field::Get(eBlastOpt_OutOfFrameMode).Match(p)) {
367 bo.SetOutOfFrameMode(v.GetBoolean());
368 } else {
369 found = false;
370 }
371 break;
372
373 case 'N':
374 if (CBlast4Field::Get(eBlastOpt_NegativeGiList).Match(p)) {
375 #ifdef NCBI_STRICT_GI
376 const list<int>& int_list = v.GetInteger_list();
377 list<TGi> gi_list;
378 ITERATE ( list<int>, it, int_list ) {
379 gi_list.push_back(GI_FROM(int, *it));
380 }
381 m_NegativeGiList = gi_list;
382 #elif defined(NCBI_INT8_GI)
383 if (v.IsBig_integer_list()) {
384 m_NegativeGiList = v.GetBig_integer_list();
385 }
386 else {
387 const list<int>& int_list = v.GetInteger_list();
388 list<TGi> gi_list;
389 ITERATE ( list<int>, it, int_list ) {
390 gi_list.push_back(GI_FROM(int, *it));
391 }
392 m_NegativeGiList = gi_list;
393 }
394 #else
395 m_NegativeGiList = v.GetInteger_list();
396 #endif
397 }
398 else if (CBlast4Field::Get(eBlastOpt_NegativeTaxidList).Match(p)) {
399 #ifdef NCBI_STRICT_TAX_ID
400 const list<int>& int_list = v.GetInteger_list();
401 list<TTaxId> taxid_list;
402 ITERATE(list<int>, it, int_list) {
403 taxid_list.push_back(TAX_ID_FROM(int, *it));
404 }
405 m_NegativeTaxidList = taxid_list;
406 #else
407 m_NegativeTaxidList = v.GetInteger_list();
408 #endif
409 } else {
410 found = false;
411 }
412 break;
413
414 case 'P':
415 if (CBlast4Field::Get(eBlastOpt_PHIPattern).Match(p)) {
416 if (v.GetString() != "") {
417 bool is_na = !! Blast_QueryIsNucleotide(bo.GetProgramType());
418 bo.SetPHIPattern(v.GetString().c_str(), is_na);
419 }
420 } else if (CBlast4Field::Get(eBlastOpt_PercentIdentity).Match(p)) {
421 opts.SetPercentIdentity(v.GetReal());
422 } else if (CBlast4Field::Get(eBlastOpt_PseudoCount).Match(p)) {
423 bo.SetPseudoCount(v.GetInteger());
424 } else {
425 found = false;
426 }
427 break;
428
429 case 'Q':
430 if (CBlast4Field::Get(eBlastOpt_QueryGeneticCode).Match(p)) {
431 bo.SetQueryGeneticCode(v.GetInteger());
432 }else if (CBlast4Field::Get(eBlastOpt_QueryCovHspPerc).Match(p)) {
433 opts.SetQueryCovHspPerc(v.GetReal());
434 } else {
435 found = false;
436 }
437 break;
438
439 case 'R':
440 if (CBlast4Field::Get(eBlastOpt_RepeatFiltering).Match(p)) {
441 bo.SetRepeatFiltering(v.GetBoolean());
442 } else if (CBlast4Field::Get(eBlastOpt_RepeatFilteringDB).Match(p)) {
443 bo.SetRepeatFilteringDB(v.GetString().c_str());
444 } else if (CBlast4Field::Get(eBlastOpt_RequiredStart).Match(p)) {
445 m_QueryRange.SetFrom(v.GetInteger());
446 } else if (CBlast4Field::Get(eBlastOpt_RequiredEnd).Match(p)) {
447 m_QueryRange.SetToOpen(v.GetInteger());
448 } else {
449 found = false;
450 }
451 break;
452
453 case 'S':
454 if (CBlast4Field::Get(eBlastOpt_StrandOption).Match(p)) {
455 // These encodings use the same values.
456 ENa_strand strand = (ENa_strand) v.GetStrand_type();
457 bo.SetStrandOption(strand);
458 } else if (CBlast4Field::Get(eBlastOpt_SegFiltering).Match(p)) {
459 bo.SetSegFiltering(v.GetBoolean());
460 } else if (CBlast4Field::Get(eBlastOpt_SegFilteringWindow).Match(p)) {
461 bo.SetSegFilteringWindow(v.GetInteger());
462 } else if (CBlast4Field::Get(eBlastOpt_SegFilteringLocut).Match(p)) {
463 bo.SetSegFilteringLocut(v.GetReal());
464 } else if (CBlast4Field::Get(eBlastOpt_SegFilteringHicut).Match(p)) {
465 bo.SetSegFilteringHicut(v.GetReal());
466 } else if (CBlast4Field::Get(eBlastOpt_SumStatisticsMode).Match(p)) {
467 bo.SetSumStatisticsMode(v.GetBoolean());
468 } else if (CBlast4Field::Get(eBlastOpt_SmithWatermanMode).Match(p)) {
469 bo.SetSmithWatermanMode(v.GetBoolean());
470 } else if (CBlast4Field::Get(eBlastOpt_SubjectMaskingType).Match(p)) {
471 m_SubjectMaskingType = (ESubjectMaskingType) v.GetInteger();
472 } else {
473 found = false;
474 }
475 break;
476
477 case 'T':
478 if (CBlast4Field::Get(eBlastOpt_TaxidList).Match(p)) {
479 #ifdef NCBI_STRICT_TAX_ID
480 const list<int>& int_list = v.GetInteger_list();
481 list<TTaxId> taxid_list;
482 ITERATE(list<int>, it, int_list) {
483 taxid_list.push_back(TAX_ID_FROM(int, *it));
484 }
485 m_TaxidList = taxid_list;
486 #else
487 m_TaxidList = v.GetInteger_list();
488 #endif
489 } else {
490 found = false;
491 }
492 break;
493 case 'U':
494 if (CBlast4Field::Get(eBlastOpt_GappedMode).Match(p)) {
495 // Notes: (1) this is the inverse of the corresponding
496 // blast4 concept (2) blast4 always returns this option
497 // regardless of whether the value matches the default.
498
499 opts.SetGappedMode(! v.GetBoolean());
500 } else if (CBlast4Field::Get(eBlastOpt_UnifiedP).Match(p)) {
501 bo.SetUnifiedP(v.GetInteger());
502 } else if (CBlast4Field::Get(eBlastOpt_SubjectBestHit).Match(p)) {
503 bo.SetSubjectBestHit();
504 } else {
505 found = false;
506 }
507 break;
508
509 case 'W':
510 if (CBlast4Field::Get(eBlastOpt_WindowMaskerTaxId).Match(p)) {
511 opts.SetOptions().SetWindowMaskerTaxId(v.GetInteger());
512 } else if (CBlast4Field::Get(eBlastOpt_WindowSize).Match(p)) {
513 opts.SetWindowSize(v.GetInteger());
514 } else if (CBlast4Field::Get(eBlastOpt_WordSize).Match(p)) {
515 bo.SetWordSize(v.GetInteger());
516 } else if (CBlast4Field::Get(eBlastOpt_WordThreshold).Match(p)) {
517 bo.SetWordThreshold(v.GetInteger());
518 } else {
519 found = false;
520 }
521 break;
522
523 default:
524 found = false;
525 }
526
527 if (! found) {
528 if (m_IgnoreUnsupportedOptions)
529 return;
530
531 string msg = "Internal: Error processing option [";
532 msg += nm;
533 msg += "] type [";
534 msg += NStr::IntToString((int) v.Which());
535 msg += "].";
536
537 NCBI_THROW(CRemoteBlastException,
538 eServiceNotAvailable,
539 msg);
540 }
541 }
542
543 void
x_ProcessOptions(CBlastOptionsHandle & opts,const TValueList * L)544 CBlastOptionsBuilder::x_ProcessOptions(CBlastOptionsHandle & opts,
545 const TValueList * L)
546 {
547 if ( !L ) {
548 return;
549 }
550
551 ITERATE(TValueList, iter, *L) {
552 CBlast4_parameter & p = *const_cast<CBlast4_parameter*>(& **iter);
553 x_ProcessOneOption(opts, p);
554 }
555 }
556
x_ApplyInteractions(CBlastOptionsHandle & boh)557 void CBlastOptionsBuilder::x_ApplyInteractions(CBlastOptionsHandle & boh)
558 {
559 CBlastOptions & bo = boh.SetOptions();
560
561 if (m_PerformCulling) {
562 bo.SetCullingLimit(m_HspRangeMax);
563 }
564 if (m_ForceMbIndex) {
565 bo.SetUseIndex(true, m_MbIndexName, m_ForceMbIndex);
566 }
567 }
568
569 /// Finder class for matching CBlast4_parameter
570 struct SBlast4ParamFinder
571 {
SBlast4ParamFinderSBlast4ParamFinder572 SBlast4ParamFinder(EBlastOptIdx opt_idx)
573 : m_Target2Find(CBlast4Field::Get(opt_idx)) {}
operator ()SBlast4ParamFinder574 bool operator()(const CRef<CBlast4_parameter>& rhs) {
575 return rhs.NotEmpty() ? m_Target2Find.Match(*rhs) : false;
576 }
577 private:
578 CBlast4Field& m_Target2Find;
579 };
580
581 EProgram
AdjustProgram(const TValueList * L,EProgram program,const string & program_string)582 CBlastOptionsBuilder::AdjustProgram(const TValueList * L,
583 EProgram program,
584 const string & program_string)
585 {
586 if ( !L ) {
587 return program;
588 }
589
590 // PHI-BLAST pattern trumps all other options
591 if (find_if(L->begin(), L->end(),
592 SBlast4ParamFinder(eBlastOpt_PHIPattern)) != L->end()) {
593 switch(program) {
594 case ePHIBlastn:
595 case eBlastn:
596 _TRACE("Adjusting program to phiblastn");
597 return ePHIBlastn;
598
599 case ePHIBlastp:
600 case eBlastp:
601 _TRACE("Adjusting program to phiblastp");
602 return ePHIBlastp;
603
604 default:
605 {
606 string msg = "Incorrect combination of option (";
607 msg += CBlast4Field::GetName(eBlastOpt_PHIPattern);
608 msg += ") and program (";
609 msg += program_string;
610 msg += ")";
611
612 NCBI_THROW(CRemoteBlastException,
613 eServiceNotAvailable,
614 msg);
615 }
616 break;
617 }
618 }
619
620 ITERATE(TValueList, iter, *L) {
621 CBlast4_parameter & p = const_cast<CBlast4_parameter&>(**iter);
622 const CBlast4_value & v = p.GetValue();
623
624 if (CBlast4Field::Get(eBlastOpt_MBTemplateLength).Match(p)) {
625 if (v.GetInteger() != 0) {
626 _TRACE("Adjusting program to discontiguous Megablast");
627 return eDiscMegablast;
628 }
629 } else if (CBlast4Field::Get(eBlastOpt_Web_StepNumber).Match(p) ||
630 CBlast4Field::Get(eBlastOpt_Web_RunPsiBlast).Match(p) ||
631 CBlast4Field::Get(eBlastOpt_PseudoCount).Match(p) ||
632 CBlast4Field::Get(eBlastOpt_IgnoreMsaMaster).Match(p)
633 ) {
634 // FIXME: should we handle DELTA-BLAST here too?
635 _TRACE("Adjusting program to psiblast");
636 return ePSIBlast;
637 }
638 }
639
640 return program;
641 }
642
643 /// Convenience function to merge all lists into one object to facilitate
644 /// invoking AdjustProgram
645 static void
s_MergeCBlast4_parameters(const objects::CBlast4_parameters * aopts,const objects::CBlast4_parameters * popts,const objects::CBlast4_parameters * fopts,objects::CBlast4_parameters & retval)646 s_MergeCBlast4_parameters(const objects::CBlast4_parameters* aopts,
647 const objects::CBlast4_parameters* popts,
648 const objects::CBlast4_parameters* fopts,
649 objects::CBlast4_parameters& retval)
650 {
651 retval.Set().clear();
652 if (aopts) {
653 retval.Set().insert(retval.Set().end(), aopts->Get().begin(), aopts->Get().end());
654 }
655 if (popts) {
656 retval.Set().insert(retval.Set().end(), popts->Get().begin(), popts->Get().end());
657 }
658 if (fopts) {
659 retval.Set().insert(retval.Set().end(), fopts->Get().begin(), fopts->Get().end());
660 }
661 }
662
663 CRef<CBlastOptionsHandle> CBlastOptionsBuilder::
GetSearchOptions(const objects::CBlast4_parameters * aopts,const objects::CBlast4_parameters * popts,const objects::CBlast4_parameters * fopts,string * task_name)664 GetSearchOptions(const objects::CBlast4_parameters * aopts,
665 const objects::CBlast4_parameters * popts,
666 const objects::CBlast4_parameters* fopts,
667 string *task_name)
668 {
669 EProgram program = ComputeProgram(m_Program, m_Service);
670 objects::CBlast4_parameters all_params;
671 s_MergeCBlast4_parameters(aopts, popts, fopts, all_params);
672 program = AdjustProgram(&all_params.Get(), program, m_Program);
673
674 // Using eLocal allows more of the options to be returned to the user.
675
676 CRef<CBlastOptionsHandle>
677 cboh(CBlastOptionsFactory::Create(program, m_Locality));
678
679 if (task_name != NULL)
680 *task_name = EProgramToTaskName(program);
681
682 m_IgnoreQueryMasks = false;
683 x_ProcessOptions(*cboh, (aopts == NULL ? 0 : &aopts->Get()));
684
685 m_IgnoreQueryMasks = m_QueryMasks.Have();
686 x_ProcessOptions(*cboh, (popts == NULL ? 0 : &popts->Get()));
687
688 x_ApplyInteractions(*cboh);
689
690 return cboh;
691 }
692
HaveEntrezQuery()693 bool CBlastOptionsBuilder::HaveEntrezQuery()
694 {
695 return m_EntrezQuery.Have();
696 }
697
GetEntrezQuery()698 string CBlastOptionsBuilder::GetEntrezQuery()
699 {
700 return m_EntrezQuery.Get();
701 }
702
HaveFirstDbSeq()703 bool CBlastOptionsBuilder::HaveFirstDbSeq()
704 {
705 return m_FirstDbSeq.Have();
706 }
707
GetFirstDbSeq()708 int CBlastOptionsBuilder::GetFirstDbSeq()
709 {
710 return m_FirstDbSeq.Get();
711 }
712
HaveFinalDbSeq()713 bool CBlastOptionsBuilder::HaveFinalDbSeq()
714 {
715 return m_FinalDbSeq.Have();
716 }
717
GetFinalDbSeq()718 int CBlastOptionsBuilder::GetFinalDbSeq()
719 {
720 return m_FinalDbSeq.Get();
721 }
722
HaveGiList()723 bool CBlastOptionsBuilder::HaveGiList()
724 {
725 return m_GiList.Have();
726 }
727
GetGiList()728 list<TGi> CBlastOptionsBuilder::GetGiList()
729 {
730 return m_GiList.Get();
731 }
732
HasDbFilteringAlgorithmId()733 bool CBlastOptionsBuilder::HasDbFilteringAlgorithmId()
734 {
735 return m_DbFilteringAlgorithmId.Have();
736 }
737
GetDbFilteringAlgorithmId()738 int CBlastOptionsBuilder::GetDbFilteringAlgorithmId()
739 {
740 return m_DbFilteringAlgorithmId.Get();
741 }
742
HasDbFilteringAlgorithmKey()743 bool CBlastOptionsBuilder::HasDbFilteringAlgorithmKey()
744 {
745 return m_DbFilteringAlgorithmKey.Have();
746 }
747
GetDbFilteringAlgorithmKey()748 string CBlastOptionsBuilder::GetDbFilteringAlgorithmKey()
749 {
750 return m_DbFilteringAlgorithmKey.Get();
751 }
752
HasSubjectMaskingType()753 bool CBlastOptionsBuilder::HasSubjectMaskingType()
754 {
755 return m_SubjectMaskingType.Have();
756 }
757
GetSubjectMaskingType()758 ESubjectMaskingType CBlastOptionsBuilder::GetSubjectMaskingType()
759 {
760 return m_SubjectMaskingType.Get();
761 }
762
763
HaveNegativeGiList()764 bool CBlastOptionsBuilder::HaveNegativeGiList()
765 {
766 return m_NegativeGiList.Have();
767 }
768
GetNegativeGiList()769 list<TGi> CBlastOptionsBuilder::GetNegativeGiList()
770 {
771 return m_NegativeGiList.Get();
772 }
773
HaveQueryMasks()774 bool CBlastOptionsBuilder::HaveQueryMasks()
775 {
776 return m_QueryMasks.Have();
777 }
778
779 CBlastOptionsBuilder::TMaskList
GetQueryMasks()780 CBlastOptionsBuilder::GetQueryMasks()
781 {
782 return m_QueryMasks.Get();
783 }
784
SetIgnoreUnsupportedOptions(bool ignore)785 void CBlastOptionsBuilder::SetIgnoreUnsupportedOptions(bool ignore)
786 {
787 m_IgnoreUnsupportedOptions = ignore;
788 }
789
HaveTaxidList()790 bool CBlastOptionsBuilder::HaveTaxidList()
791 {
792 return m_TaxidList.Have();
793 }
794
GetTaxidList()795 list<TTaxId> CBlastOptionsBuilder::GetTaxidList()
796 {
797 return m_TaxidList.Get();
798 }
799
HaveNegativeTaxidList()800 bool CBlastOptionsBuilder::HaveNegativeTaxidList()
801 {
802 return m_NegativeTaxidList.Have();
803 }
804
GetNegativeTaxidList()805 list<TTaxId> CBlastOptionsBuilder::GetNegativeTaxidList()
806 {
807 return m_NegativeTaxidList.Get();
808 }
809
810
811 END_SCOPE(blast)
812 END_NCBI_SCOPE
813
814 /* @} */
815