~ubuntu-branches/ubuntu/intrepid/primer3/intrepid

« back to all changes in this revision

Viewing changes to src/primer3.c

  • Committer: Bazaar Package Importer
  • Author(s): Charles Plessy
  • Date: 2006-09-28 20:18:54 UTC
  • Revision ID: james.westby@ubuntu.com-20060928201854-45pwapz5e3a6d684
Tags: upstream-1.0b
ImportĀ upstreamĀ versionĀ 1.0b

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/*
 
2
Copyright (c) 1996,1997,1998,1999,2000,2001,2004,2006
 
3
Whitehead Institute for Biomedical Research, Steve Rozen
 
4
(http://jura.wi.mit.edu/rozen), and Helen Skaletsky
 
5
All rights reserved.
 
6
 
 
7
Redistribution and use in source and binary forms, with or without
 
8
modification, are permitted provided that the following conditions are
 
9
met:
 
10
 
 
11
   * Redistributions of source code must retain the above copyright
 
12
notice, this list of conditions and the following disclaimer.
 
13
   * Redistributions in binary form must reproduce the above
 
14
copyright notice, this list of conditions and the following disclaimer
 
15
in the documentation and/or other materials provided with the
 
16
distribution.
 
17
   * Neither the names of the copyright holders nor contributors may
 
18
be used to endorse or promote products derived from this software
 
19
without specific prior written permission.
 
20
 
 
21
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 
22
"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 
23
LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
 
24
A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
 
25
OWNERS OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 
26
SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
 
27
LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
 
28
DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
 
29
THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 
30
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 
31
OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 
32
*/
 
33
 
 
34
#include <string.h>
 
35
#include <float.h>
 
36
#include "primer3.h"
 
37
#include "primer3_release.h"
 
38
 
 
39
#define MACRO_VALUE_AS_STRING(A) MACRO_STRING(A)
 
40
 
 
41
static int    check_intervals(const char *, const int,
 
42
                              interval_array_t, const int, seq_args *);
 
43
static char   dna_to_upper(char *, int);
 
44
static char   *strstr_nocase(char *, char *);
 
45
 
 
46
 
 
47
/* Default parameter values.  */
 
48
#define OPT_SIZE            20
 
49
#define MIN_SIZE             18
 
50
#define MAX_SIZE             27
 
51
 
 
52
#define OPT_TM             60.0
 
53
#define MIN_TM             57.0
 
54
#define MAX_TM             63.0
 
55
#define MAX_DIFF_TM       100.0
 
56
#define DEFAULT_OPT_GC_PERCENT PR_UNDEFINED_INT_OPT
 
57
#define MIN_GC             20.0
 
58
#define MAX_GC             80.0
 
59
#define SALT_CONC          50.0
 
60
#define DNA_CONC           50.0
 
61
#define NUM_NS_ACCEPTED       0
 
62
#define MAX_POLY_X            5
 
63
#define SELF_ANY            800
 
64
#define SELF_END            300
 
65
#define PAIR_COMPL_ANY      800
 
66
#define PAIR_COMPL_END      300
 
67
#define FILE_FLAG             0
 
68
#define EXPLAIN_FLAG          0
 
69
#define GC_CLAMP              0
 
70
#define LIBERAL_BASE          0
 
71
#define PICK_INTERNAL_OLIGO   0
 
72
#define PRIMER_TASK           0
 
73
#define INTERNAL_OLIGO_OPT_SIZE   20
 
74
#define INTERNAL_OLIGO_MIN_SIZE   18
 
75
#define INTERNAL_OLIGO_MAX_SIZE   27
 
76
#define INTERNAL_OLIGO_OPT_TM     60.0
 
77
#define INTERNAL_OLIGO_MIN_TM     57.0
 
78
#define INTERNAL_OLIGO_MAX_TM     63.0
 
79
#define INTERNAL_OLIGO_MIN_GC     20.0
 
80
#define INTERNAL_OLIGO_MAX_GC     80.0
 
81
#define INTERNAL_OLIGO_SALT_CONC         50.0
 
82
#define INTERNAL_OLIGO_DNA_CONC          50.0
 
83
#define INTERNAL_OLIGO_NUM_NS               0
 
84
#define INTERNAL_OLIGO_MAX_POLY_X           5 
 
85
#define INTERNAL_OLIGO_SELF_ANY          1200
 
86
#define INTERNAL_OLIGO_SELF_END          1200
 
87
#define INTERNAL_OLIGO_REPEAT_SIMILARITY 1200
 
88
#define REPEAT_SIMILARITY                1200
 
89
#define PAIR_REPEAT_SIMILARITY           2400
 
90
#define FIRST_BASE_INDEX            0
 
91
#define NUM_RETURN                  5
 
92
#define MIN_QUALITY                 0
 
93
#define QUALITY_RANGE_MIN           0
 
94
#define QUALITY_RANGE_MAX         100
 
95
#define DEFAULT_MAX_END_STABILITY    100.0
 
96
#define PRIMER_PRODUCT_OPT_SIZE      PR_UNDEFINED_INT_OPT
 
97
#define PRIMER_PRODUCT_OPT_TM        PR_UNDEFINED_DBL_OPT
 
98
#define MAX_TEMPLATE_MISPRIMING      PR_UNDEFINED_ALIGN_OPT
 
99
#define PAIR_MAX_TEMPLATE_MISPRIMING PR_UNDEFINED_ALIGN_OPT
 
100
#define IO_MAX_TEMPLATE_MISHYB       PR_UNDEFINED_ALIGN_OPT
 
101
 
 
102
#define LIB_AMBIGUITY_CODES_CONSENSUS 1
 
103
/*  For backward compatibility. It turns out that
 
104
    this _not_ what one normally wants, since many
 
105
    libraries contain strings of N, which then match
 
106
    every oligo (very bad).
 
107
*/
 
108
 
 
109
/* Weights for objective functions for oligos and pairs. */
 
110
#define PRIMER_WT_TM_GT          1
 
111
#define PRIMER_WT_TM_LT          1
 
112
#define PRIMER_WT_SIZE_LT        1
 
113
#define PRIMER_WT_SIZE_GT        1
 
114
#define PRIMER_WT_GC_PERCENT_LT  0
 
115
#define PRIMER_WT_GC_PERCENT_GT  0
 
116
#define PRIMER_WT_COMPL_ANY      0
 
117
#define PRIMER_WT_COMPL_END      0
 
118
#define PRIMER_WT_NUM_NS         0
 
119
#define PRIMER_WT_REP_SIM        0
 
120
#define PRIMER_WT_SEQ_QUAL       0
 
121
#define PRIMER_WT_END_QUAL       0
 
122
#define PRIMER_WT_POS_PENALTY    1
 
123
#define PRIMER_WT_END_STABILITY  0
 
124
 
 
125
#define IO_WT_TM_GT          1
 
126
#define IO_WT_TM_LT          1
 
127
#define IO_WT_SIZE_LT        1
 
128
#define IO_WT_SIZE_GT        1
 
129
#define IO_WT_GC_PERCENT_LT  0
 
130
#define IO_WT_GC_PERCENT_GT  0
 
131
#define IO_WT_COMPL_ANY      0
 
132
#define IO_WT_COMPL_END      0
 
133
#define IO_WT_NUM_NS         0
 
134
#define IO_WT_REP_SIM        0
 
135
#define IO_WT_SEQ_QUAL       0
 
136
#define IO_WT_END_QUAL       0
 
137
 
 
138
#define PAIR_WT_PRIMER_PENALTY      1
 
139
#define PAIR_WT_IO_PENALTY          0
 
140
#define PAIR_WT_DIFF_TM             0
 
141
#define PAIR_WT_COMPL_ANY           0
 
142
#define PAIR_WT_COMPL_END           0
 
143
#define PAIR_WT_REP_SIM             0
 
144
#define PAIR_WT_PRODUCT_TM_LT       0
 
145
#define PAIR_WT_PRODUCT_TM_GT       0
 
146
#define PAIR_WT_PRODUCT_SIZE_LT     0
 
147
#define PAIR_WT_PRODUCT_SIZE_GT     0
 
148
 
 
149
void
 
150
pr_set_default_global_args(a)
 
151
    primer_args *a;
 
152
{
 
153
    memset(a, 0, sizeof(*a));  
 
154
    a->primer_opt_size  = OPT_SIZE;
 
155
    a->primer_min_size  = MIN_SIZE;
 
156
    a->primer_max_size  = MAX_SIZE;
 
157
    a->opt_tm           = OPT_TM;
 
158
    a->min_tm           = MIN_TM;
 
159
    a->max_tm           = MAX_TM;
 
160
    a->max_diff_tm      = MAX_DIFF_TM;
 
161
    a->min_gc           = MIN_GC;
 
162
    a->opt_gc_content   = DEFAULT_OPT_GC_PERCENT;
 
163
    a->max_gc           = MAX_GC;
 
164
    a->salt_conc        = SALT_CONC;
 
165
    a->dna_conc         = DNA_CONC;
 
166
    a->num_ns_accepted  = NUM_NS_ACCEPTED;
 
167
    a->self_any         = SELF_ANY;
 
168
    a->self_end         = SELF_END;
 
169
    a->pair_compl_any   = PAIR_COMPL_ANY;
 
170
    a->pair_compl_end   = PAIR_COMPL_END;
 
171
    a->file_flag        = FILE_FLAG;
 
172
    a->explain_flag     = EXPLAIN_FLAG;
 
173
    a->gc_clamp         = GC_CLAMP;
 
174
    a->max_poly_x       = MAX_POLY_X;
 
175
    a->liberal_base      = LIBERAL_BASE;
 
176
    a->primer_task       = PRIMER_TASK;
 
177
    a->first_base_index  = FIRST_BASE_INDEX;
 
178
    a->num_return        = NUM_RETURN;
 
179
    a->pr_min[0]         = 100;
 
180
    a->pr_max[0]         = 300;
 
181
    a->num_intervals     = 1;
 
182
    a->repeat_compl      = REPEAT_SIMILARITY;
 
183
    a->pair_repeat_compl = PAIR_REPEAT_SIMILARITY;
 
184
    a->min_quality       = MIN_QUALITY;
 
185
    a->min_end_quality   = MIN_QUALITY;
 
186
    a->quality_range_min = QUALITY_RANGE_MIN;
 
187
    a->quality_range_max = QUALITY_RANGE_MAX;
 
188
    a->outside_penalty   = PR_DEFAULT_OUTSIDE_PENALTY;
 
189
    a->inside_penalty    = PR_DEFAULT_INSIDE_PENALTY;
 
190
    a->max_end_stability = DEFAULT_MAX_END_STABILITY;
 
191
    a->product_max_tm    = PR_DEFAULT_PRODUCT_MAX_TM;
 
192
    a->product_min_tm    = PR_DEFAULT_PRODUCT_MIN_TM;
 
193
    a->product_opt_tm    = PRIMER_PRODUCT_OPT_TM;
 
194
    a->product_opt_size  = PRIMER_PRODUCT_OPT_SIZE;
 
195
    a->max_template_mispriming
 
196
                          = MAX_TEMPLATE_MISPRIMING;
 
197
    a->pair_max_template_mispriming
 
198
                          = PAIR_MAX_TEMPLATE_MISPRIMING;
 
199
 
 
200
    a->io_primer_opt_size = INTERNAL_OLIGO_OPT_SIZE;
 
201
    a->io_primer_min_size = INTERNAL_OLIGO_MIN_SIZE;
 
202
    a->io_primer_max_size = INTERNAL_OLIGO_MAX_SIZE;
 
203
    a->io_opt_tm          = INTERNAL_OLIGO_OPT_TM;
 
204
    a->io_min_tm          = INTERNAL_OLIGO_MIN_TM;
 
205
    a->io_max_tm          = INTERNAL_OLIGO_MAX_TM;
 
206
    a->io_min_gc          = INTERNAL_OLIGO_MIN_GC;
 
207
    a->io_opt_gc_content  = DEFAULT_OPT_GC_PERCENT;
 
208
    a->io_max_gc          = INTERNAL_OLIGO_MAX_GC;
 
209
    a->io_max_poly_x      = INTERNAL_OLIGO_MAX_POLY_X;
 
210
    a->io_salt_conc       = INTERNAL_OLIGO_SALT_CONC;
 
211
    a->io_dna_conc        = INTERNAL_OLIGO_DNA_CONC;
 
212
    a->io_num_ns_accepted = INTERNAL_OLIGO_NUM_NS;
 
213
    a->io_self_any        = INTERNAL_OLIGO_SELF_ANY;
 
214
    a->io_self_end        = INTERNAL_OLIGO_SELF_END;
 
215
    a->io_repeat_compl    = INTERNAL_OLIGO_REPEAT_SIMILARITY;
 
216
    a->io_min_quality     = MIN_QUALITY;
 
217
    a->io_min_end_quality = MIN_QUALITY;
 
218
    a->io_max_template_mishyb
 
219
                          = IO_MAX_TEMPLATE_MISHYB;
 
220
 
 
221
    a->primer_weights.temp_gt       = PRIMER_WT_TM_GT;
 
222
    a->primer_weights.temp_lt       = PRIMER_WT_TM_LT;
 
223
    a->primer_weights.length_gt     = PRIMER_WT_SIZE_GT;
 
224
    a->primer_weights.length_lt     = PRIMER_WT_SIZE_LT;
 
225
    a->primer_weights.gc_content_gt = PRIMER_WT_GC_PERCENT_GT;
 
226
    a->primer_weights.gc_content_lt = PRIMER_WT_GC_PERCENT_LT;
 
227
    a->primer_weights.compl_any     = PRIMER_WT_COMPL_ANY;
 
228
    a->primer_weights.compl_end     = PRIMER_WT_COMPL_END;
 
229
    a->primer_weights.num_ns        = PRIMER_WT_NUM_NS;
 
230
    a->primer_weights.repeat_sim    = PRIMER_WT_REP_SIM;
 
231
    a->primer_weights.seq_quality   = PRIMER_WT_SEQ_QUAL;
 
232
    a->primer_weights.end_quality   = PRIMER_WT_END_QUAL;
 
233
    a->primer_weights.pos_penalty   = PRIMER_WT_POS_PENALTY;
 
234
    a->primer_weights.end_stability = PRIMER_WT_END_STABILITY;
 
235
 
 
236
    a->io_weights.temp_gt     = IO_WT_TM_GT;
 
237
    a->io_weights.temp_lt     = IO_WT_TM_LT;
 
238
    a->io_weights.length_gt   = IO_WT_SIZE_GT;
 
239
    a->io_weights.length_lt   = IO_WT_SIZE_LT;
 
240
    a->io_weights.gc_content_gt = IO_WT_GC_PERCENT_GT;
 
241
    a->io_weights.gc_content_lt = IO_WT_GC_PERCENT_LT;
 
242
    a->io_weights.compl_any   = IO_WT_COMPL_ANY;
 
243
    a->io_weights.compl_end   = IO_WT_COMPL_END;
 
244
    a->io_weights.num_ns      = IO_WT_NUM_NS;
 
245
    a->io_weights.repeat_sim  = IO_WT_REP_SIM;
 
246
    a->io_weights.seq_quality = IO_WT_SEQ_QUAL;
 
247
    a->io_weights.end_quality = IO_WT_END_QUAL;
 
248
 
 
249
    a->pr_pair_weights.primer_quality  = PAIR_WT_PRIMER_PENALTY;
 
250
    a->pr_pair_weights.io_quality      = PAIR_WT_IO_PENALTY;
 
251
    a->pr_pair_weights.diff_tm         = PAIR_WT_DIFF_TM;
 
252
    a->pr_pair_weights.compl_any       = PAIR_WT_COMPL_ANY;
 
253
    a->pr_pair_weights.compl_end       = PAIR_WT_COMPL_END;
 
254
    a->pr_pair_weights.repeat_sim      = PAIR_WT_REP_SIM;
 
255
    a->pr_pair_weights.product_tm_lt   = PAIR_WT_PRODUCT_TM_LT;
 
256
    a->pr_pair_weights.product_tm_gt   = PAIR_WT_PRODUCT_TM_GT;
 
257
    a->pr_pair_weights.product_size_lt = PAIR_WT_PRODUCT_SIZE_LT;
 
258
    a->pr_pair_weights.product_size_gt = PAIR_WT_PRODUCT_SIZE_GT;
 
259
 
 
260
    /* a->short_match = 0; not used */
 
261
 
 
262
    a->lib_ambiguity_codes_consensus   = LIB_AMBIGUITY_CODES_CONSENSUS;
 
263
}
 
264
 
 
265
/*
 
266
 * Return 1 on error, 0 on success.  Set sa->trimmed_seq and possibly modify
 
267
 * sa->tar.  Upcase and check all bases in sa->trimmed_seq.
 
268
 */
 
269
int
 
270
_pr_data_control(pa, sa)
 
271
    primer_args *pa;
 
272
    seq_args *sa;
 
273
{
 
274
    static char s1[MAX_PRIMER_LENGTH+1];
 
275
    int i, pr_min;
 
276
    int seq_len = strlen(sa->sequence);
 
277
    char offending_char = '\0';
 
278
 
 
279
    if (pa->io_max_template_mishyb >= 0)
 
280
      pr_append_new_chunk(&pa->glob_err,
 
281
                          "PRIMER_INTERNAL_OLIGO_MAX_TEMPLATE_MISHYB is not supported");
 
282
 
 
283
    if (pa->primer_min_size < 1)
 
284
      pr_append_new_chunk(&pa->glob_err, "PRIMER_MIN_SIZE must be >= 1");
 
285
 
 
286
    if (pa->primer_max_size > MAX_PRIMER_LENGTH) {
 
287
      pr_append_new_chunk(&pa->glob_err,
 
288
                          "PRIMER_MAX_SIZE exceeds built-in maximum of ");
 
289
      pr_append(&pa->glob_err, MACRO_VALUE_AS_STRING(MAX_PRIMER_LENGTH));
 
290
      return 1;
 
291
    }
 
292
 
 
293
    if (pa->primer_opt_size > pa->primer_max_size) {
 
294
        pr_append_new_chunk(&pa->glob_err,
 
295
                            "PRIMER_{OPT,DEFAULT}_SIZE > PRIMER_MAX_SIZE");
 
296
        return 1;
 
297
    }
 
298
 
 
299
    if (pa->primer_opt_size < pa->primer_min_size) {
 
300
        pr_append_new_chunk(&pa->glob_err,
 
301
                            "PRIMER_{OPT,DEFAULT}_SIZE < PRIMER_MIN_SIZE");
 
302
        return 1;
 
303
    }
 
304
 
 
305
    if (pa->io_primer_max_size > MAX_PRIMER_LENGTH) {
 
306
        pr_append_new_chunk(&pa->glob_err,
 
307
                  "PRIMER_INTERNAL_OLIGO_MAX_SIZE exceeds built-in maximum");
 
308
        return 1;
 
309
    }
 
310
 
 
311
    if (pa->io_primer_opt_size > pa->io_primer_max_size) {
 
312
        pr_append_new_chunk(&pa->glob_err,
 
313
                  "PRIMER_INTERNAL_OLIGO_{OPT,DEFAULT}_SIZE > MAX_SIZE");
 
314
        return 1;
 
315
    }
 
316
 
 
317
    if (pa->io_primer_opt_size < pa->io_primer_min_size) {
 
318
        pr_append_new_chunk(&pa->glob_err,
 
319
                  "PRIMER_INTERNAL_OLIGO_{OPT,DEFAULT}_SIZE < MIN_SIZE");
 
320
        return 1;
 
321
    }
 
322
 
 
323
    if (pa->gc_clamp > pa->primer_min_size) {
 
324
        pr_append_new_chunk(&pa->glob_err,
 
325
                            "PRIMER_GC_CLAMP > PRIMER_MIN_SIZE");
 
326
        return 1;
 
327
    }
 
328
 
 
329
    if (NULL == sa->sequence_name && pa->file_flag) {
 
330
        pr_append_new_chunk(&sa->error,
 
331
                            "Need PRIMER_SEQUENCE_ID if PRIMER_FILE_FLAG != 0");
 
332
        return 1;
 
333
    }
 
334
 
 
335
    if (0 == pa->num_intervals) {
 
336
        pr_append_new_chunk(&pa->glob_err,
 
337
                            "Empty value for PRIMER_PRODUCT_SIZE_RANGE");
 
338
        return 1;
 
339
    }
 
340
    for (i = 0; i < pa->num_intervals; i++) {
 
341
        if (pa->pr_min[i] > pa->pr_max[i] || pa->pr_min[i] < 0) {
 
342
            pr_append_new_chunk(&pa->glob_err,
 
343
                                "Illegal element in PRIMER_PRODUCT_SIZE_RANGE");
 
344
            return 1;
 
345
        }
 
346
    }
 
347
 
 
348
    pr_min = INT_MAX;
 
349
    for(i=0;i<pa->num_intervals;i++)
 
350
        if(pa->pr_min[i]<pr_min) pr_min=pa->pr_min[i];
 
351
 
 
352
    if (pa->primer_max_size > pr_min) {
 
353
        pr_append_new_chunk(&pa->glob_err,
 
354
                            "PRIMER_MAX_SIZE > min PRIMER_PRODUCT_SIZE_RANGE");
 
355
        return 1;
 
356
    }
 
357
 
 
358
    if ((pick_pcr_primers_and_hyb_probe == pa->primer_task 
 
359
         || pick_hyb_probe_only == pa->primer_task)
 
360
        && pa->io_primer_max_size > pr_min) {
 
361
        pr_append_new_chunk(&pa->glob_err,
 
362
                 "PRIMER_INTERNAL_OLIGO_MAX_SIZE > min PRIMER_PRODUCT_SIZE_RANGE");
 
363
        return 1;
 
364
    }
 
365
 
 
366
    if (pa->num_return < 1) {
 
367
        pr_append_new_chunk(&pa->glob_err,
 
368
                            "PRIMER_NUM_RETURN < 1");
 
369
        return 1;
 
370
    }
 
371
 
 
372
    if (sa->incl_l >= INT_MAX) {
 
373
        pr_append_new_chunk(&sa->error, "Value for INCLUDED_REGION too large");
 
374
        return 1;
 
375
    }
 
376
 
 
377
    if (sa->incl_s < 0 || sa->incl_l < 0 
 
378
        || sa->incl_s + sa->incl_l > seq_len) {
 
379
        pr_append_new_chunk(&sa->error, "Illegal value for INCLUDED_REGION");
 
380
        return 1;
 
381
    }
 
382
    
 
383
    if (sa->incl_l < pr_min && pa->primer_task != pick_hyb_probe_only
 
384
        && pa->primer_task != pick_left_only
 
385
        && pa->primer_task != pick_right_only) {
 
386
        pr_append_new_chunk(&sa->error,
 
387
           "INCLUDED_REGION length < min PRIMER_PRODUCT_SIZE_RANGE");
 
388
        return 1;
 
389
    }
 
390
 
 
391
    if (pa->max_end_stability < 0) {
 
392
        pr_append_new_chunk(&sa->error,
 
393
                            "PRIMER_MAX_END_STABILITY must be non-negative");
 
394
        return 1;
 
395
    }
 
396
 
 
397
    if (!PR_START_CODON_POS_IS_NULL(sa)) {
 
398
      if (!PR_POSITION_PENALTY_IS_NULL(pa)) {
 
399
        pr_append_new_chunk(&sa->error,
 
400
           "Cannot accept both PRIMER_START_CODON_POSITION and non-default ");
 
401
        pr_append(&sa->error,
 
402
           "arguments for PRIMER_INSIDE_PENALTY or PRIMER_OUTSIDE_PENALTY");
 
403
      }
 
404
      if (sa->start_codon_pos  > (sa->incl_s + sa->incl_l - 3))
 
405
        pr_append_new_chunk(&sa->error,
 
406
           "Start codon position not contained in INCLUDED_REGION");
 
407
      else {
 
408
        if (sa->start_codon_pos >= 0
 
409
            && ((sa->sequence[sa->start_codon_pos] != 'A'
 
410
                 && sa->sequence[sa->start_codon_pos] != 'a')
 
411
                || (sa->sequence[sa->start_codon_pos + 1] != 'T'
 
412
                    && sa->sequence[sa->start_codon_pos + 1] != 't')
 
413
                || (sa->sequence[sa->start_codon_pos + 2] != 'G'
 
414
                    && sa->sequence[sa->start_codon_pos + 2] != 'g')))
 
415
          pr_append_new_chunk(&sa->error,
 
416
                              "No start codon at PRIMER_START_CODON_POSITION");
 
417
      }
 
418
    }
 
419
 
 
420
    sa->trimmed_seq = pr_safe_malloc(sa->incl_l + 1);
 
421
    _pr_substr(sa->sequence, sa->incl_s, sa->incl_l, sa->trimmed_seq);
 
422
 
 
423
    sa->upcased_seq = pr_safe_malloc(strlen(sa->sequence) + 1);
 
424
    strcpy(sa->upcased_seq, sa->sequence);
 
425
    if ((offending_char = dna_to_upper(sa->upcased_seq, 1))) {
 
426
      offending_char = '\0';
 
427
      /* NEW */
 
428
      /* TODO add warning or error (depending on liberal base)
 
429
         here. */
 
430
    }
 
431
    sa->upcased_seq_r = pr_safe_malloc(strlen(sa->sequence) + 1);
 
432
    _pr_reverse_complement(sa->upcased_seq, sa->upcased_seq_r);
 
433
 
 
434
    if (check_intervals("TARGET", sa->num_targets, sa->tar, seq_len, sa)
 
435
        == 1) return 1;
 
436
    sa->start_codon_pos -= sa->incl_s;
 
437
 
 
438
    if (check_intervals("EXCLUDED_REGION", sa->num_excl, sa->excl,
 
439
                        seq_len, sa)
 
440
        == 1) return 1;
 
441
 
 
442
    if (check_intervals("PRIMER_INTERNAL_OLIGO_EXCLUDED_REGION",
 
443
                        sa->num_internal_excl, sa->excl_internal,
 
444
                        seq_len, sa)
 
445
        == 1) return 1;
 
446
 
 
447
    if (NULL != sa->quality){
 
448
        if(pa->min_quality != 0 && pa->min_quality < pa->quality_range_min) {
 
449
           pr_append_new_chunk(&pa->glob_err,
 
450
               "PRIMER_MIN_QUALITY < PRIMER_QUALITY_RANGE_MIN");
 
451
           return 1;
 
452
        }
 
453
        if(pa->min_quality != 0 && pa->min_quality > pa->quality_range_max) {
 
454
           pr_append_new_chunk(&pa->glob_err,
 
455
               "PRIMER_MIN_QUALITY > PRIMER_QUALITY_RANGE_MAX");
 
456
           return 1;
 
457
        }
 
458
        if(pa->io_min_quality != 0 && pa->io_min_quality <pa->quality_range_min) {
 
459
           pr_append_new_chunk(&pa->glob_err,
 
460
            "PRIMER_INTERNAL_OLIGO_MIN_QUALITY < PRIMER_QUALITY_RANGE_MIN");
 
461
           return 1;
 
462
        }
 
463
        if(pa->io_min_quality != 0 && pa->io_min_quality > pa->quality_range_max) {
 
464
           pr_append_new_chunk(&pa->glob_err,
 
465
             "PRIMER_INTERNAL_OLIGO_MIN_QUALITY > PRIMER_QUALITY_RANGE_MAX");
 
466
           return 1;
 
467
        }
 
468
        for(i=0; i< seq_len; i++) {
 
469
           if(sa->quality[i] < pa->quality_range_min ||
 
470
                 sa->quality[i] > pa->quality_range_max) {
 
471
             pr_append_new_chunk(&sa->error,
 
472
                "Sequence quality score out of range");
 
473
             return 1;
 
474
           }
 
475
        }
 
476
    }
 
477
    else if (pa->primer_weights.seq_quality || pa->io_weights.seq_quality) {
 
478
         pr_append_new_chunk(&sa->error,
 
479
              "Sequence quality is part of objective function but sequence quality is not defined");
 
480
         return 1;
 
481
    }
 
482
 
 
483
    if ((offending_char = dna_to_upper(sa->trimmed_seq, 0))) {
 
484
      if (pa->liberal_base) {
 
485
        pr_append_new_chunk(&sa->warning,
 
486
                            "Unrecognized base in input sequence");
 
487
      }
 
488
      else {
 
489
        pr_append_new_chunk(&sa->error, "Unrecognized base in input sequence");
 
490
        return 1;
 
491
      }
 
492
    }
 
493
 
 
494
    if(pa->opt_tm < pa->min_tm || pa->opt_tm > pa->max_tm) {
 
495
         pr_append_new_chunk(&pa->glob_err,
 
496
                             "Optimum primer Tm lower than minimum or higher than maximum");
 
497
         return 1;
 
498
    }
 
499
    if(pa->io_opt_tm < pa->io_min_tm || pa->io_opt_tm > pa->io_max_tm) {
 
500
         pr_append_new_chunk(&pa->glob_err,
 
501
                             "Illegal values for PRIMER_INTERNAL_OLIGO_TM");
 
502
         return 1;
 
503
    }
 
504
    if(pa->min_gc>pa->max_gc||pa->min_gc>100||pa->max_gc<0){
 
505
         pr_append_new_chunk(&pa->glob_err,
 
506
                             "Illegal value for PRIMER_MAX_GC and PRIMER_MIN_GC");
 
507
         return 1;
 
508
    }
 
509
    if(pa->io_min_gc>pa->io_max_gc||pa->io_min_gc>100||pa->io_max_gc<0){
 
510
         pr_append_new_chunk(&pa->glob_err,
 
511
                             "Illegal value for PRIMER_INTERNAL_OLIGO_GC");
 
512
         return 1;
 
513
    }
 
514
    if(pa->num_ns_accepted<0){
 
515
         pr_append_new_chunk(&pa->glob_err,
 
516
                             "Illegal value for PRIMER_NUM_NS_ACCEPTED");
 
517
         return 1;
 
518
    }
 
519
    if(pa->io_num_ns_accepted<0){ 
 
520
         pr_append_new_chunk(&pa->glob_err,
 
521
                             "Illegal value for PRIMER_INTERNAL_OLIGO_NUM_NS");
 
522
         return 1;
 
523
    }
 
524
    if(pa->self_any<0||pa->self_end<0
 
525
                       ||pa->pair_compl_any<0||pa->pair_compl_end<0){
 
526
         pr_append_new_chunk(&pa->glob_err,
 
527
             "Illegal value for primer complementarity restrictions");
 
528
         return 1;
 
529
    }
 
530
    if(pa->io_self_any<0||pa->io_self_end<0){
 
531
          pr_append_new_chunk(&pa->glob_err,
 
532
              "Illegal value for internal oligo complementarity restrictions");
 
533
          return 1;
 
534
    }
 
535
    if(pa->salt_conc<=0||pa->dna_conc<=0){
 
536
          pr_append_new_chunk(&pa->glob_err,
 
537
              "Illegal value for primer salt or dna concentration");
 
538
          return 1;
 
539
    }
 
540
    if(pa->io_salt_conc<=0||pa->io_dna_conc<=0){
 
541
          pr_append_new_chunk(&pa->glob_err,
 
542
              "Illegal value for internal oligo salt or dna concentration");
 
543
          return 1;
 
544
    }
 
545
    if (!_PR_DEFAULT_POSITION_PENALTIES(pa) && sa->num_targets > 1) {
 
546
      pr_append_new_chunk(&sa->error,
 
547
                          "Non-default inside penalty or outside penalty ");
 
548
      pr_append(&sa->error,
 
549
                "is valid only when number of targets <= 1");
 
550
    }
 
551
    if (!_PR_DEFAULT_POSITION_PENALTIES(pa) && 0 == sa->num_targets) {
 
552
      pr_append_new_chunk(&sa->warning,
 
553
                          "Non-default inside penalty or outside penalty ");
 
554
      pr_append(&sa->warning,
 
555
                "has no effect when number of targets is 0");
 
556
    }
 
557
    if (pa->primer_task != pick_pcr_primers_and_hyb_probe 
 
558
        && pa->primer_task != pick_hyb_probe_only
 
559
        && sa->internal_input) {
 
560
      pr_append_new_chunk(&sa->error,
 
561
                          "Not specified to pick internal oligos");
 
562
      pr_append(&sa->error,
 
563
                " but a specific internal oligo is provided");
 
564
    }
 
565
    if (sa->internal_input) {
 
566
      if (strlen(sa->internal_input) > pa->io_primer_max_size)
 
567
        pr_append_new_chunk(&sa->error, "Specified internal oligo too long");
 
568
 
 
569
      if (strlen(sa->internal_input) < pa->io_primer_min_size)
 
570
        pr_append_new_chunk(&sa->error, "Specified internal oligo too short");
 
571
 
 
572
      if (!strstr_nocase(sa->sequence, sa->internal_input))
 
573
        pr_append_new_chunk(&sa->error,
 
574
                            "Specified internal oligo not in sequence");
 
575
      else if (!strstr_nocase(sa->trimmed_seq, sa->internal_input))
 
576
        pr_append_new_chunk(&sa->error,
 
577
                            "Specified internal oligo not in Included Region");
 
578
    }
 
579
    if (sa->left_input) {
 
580
      if (strlen(sa->left_input) > pa->primer_max_size)
 
581
        pr_append_new_chunk(&sa->error, "Specified left primer too long");
 
582
      if (strlen(sa->left_input) < pa->primer_min_size)
 
583
        pr_append_new_chunk(&sa->error, "Specified left primer too short");
 
584
      if (!strstr_nocase(sa->sequence, sa->left_input))
 
585
        pr_append_new_chunk(&sa->error,
 
586
                            "Specified left primer not in sequence");
 
587
      else if (!strstr_nocase(sa->trimmed_seq, sa->left_input))
 
588
        pr_append_new_chunk(&sa->error,
 
589
                            "Specified left primer not in Included Region");
 
590
    }
 
591
    if (sa->right_input) {
 
592
      if (strlen(sa->right_input) < pa->primer_min_size)
 
593
        pr_append_new_chunk(&sa->error, "Specified right primer too short");
 
594
      if (strlen(sa->right_input) > pa->primer_max_size) {
 
595
        pr_append_new_chunk(&sa->error, "Specified right primer too long");
 
596
       } else { /* We do not want to overflow s1. */
 
597
        _pr_reverse_complement(sa->right_input,s1);
 
598
        if (!strstr_nocase(sa->sequence, s1))
 
599
          pr_append_new_chunk(&sa->error,
 
600
                              "Specified right primer not in sequence");
 
601
        else if (!strstr_nocase(sa->trimmed_seq, s1))
 
602
          pr_append_new_chunk(&sa->error,
 
603
                              "Specified right primer not in Included Region");
 
604
      }
 
605
    }
 
606
 
 
607
    if ((pa->pr_pair_weights.product_tm_lt || 
 
608
         pa->pr_pair_weights.product_tm_gt)
 
609
        && pa->product_opt_tm == PR_UNDEFINED_DBL_OPT) {
 
610
        pr_append_new_chunk(&pa->glob_err, 
 
611
           "Product temperature is part of objective function while optimum temperature is not defined");
 
612
        return 1;
 
613
     }
 
614
        
 
615
    if((pa->pr_pair_weights.product_size_lt ||
 
616
        pa->pr_pair_weights.product_size_gt) 
 
617
       && pa->product_opt_size == PR_UNDEFINED_INT_OPT){
 
618
       pr_append_new_chunk(&pa->glob_err,
 
619
          "Product size is part of objective function while optimum size is not defined");
 
620
       return 1;
 
621
    }
 
622
 
 
623
    if ((pa->primer_weights.gc_content_lt || 
 
624
         pa->primer_weights.gc_content_gt)
 
625
        && pa->opt_gc_content == DEFAULT_OPT_GC_PERCENT) {
 
626
        pr_append_new_chunk(&pa->glob_err, 
 
627
           "Primer GC content is part of objective function while optimum gc_content is not defined");
 
628
        return 1;
 
629
     }
 
630
        
 
631
    if ((pa->io_weights.gc_content_lt || 
 
632
         pa->io_weights.gc_content_gt)
 
633
        && pa->io_opt_gc_content == DEFAULT_OPT_GC_PERCENT) {
 
634
        pr_append_new_chunk(&pa->glob_err, 
 
635
           "Hyb probe GC content is part of objective function while optimum gc_content is not defined");
 
636
        return 1;
 
637
     }
 
638
        
 
639
    if ((pa->primer_task != pick_pcr_primers_and_hyb_probe 
 
640
         && pa->primer_task != pick_hyb_probe_only ) &&
 
641
                        (pa->pr_pair_weights.io_quality)) {
 
642
       pr_append_new_chunk(&pa->glob_err,
 
643
          "Internal oligo quality is part of objective function while internal oligo choice is not required");
 
644
       return 1;
 
645
    }
 
646
 
 
647
    if (pa->primer_weights.repeat_sim && (!pa->repeat_lib.seq_num)) {
 
648
       pr_append_new_chunk(&pa->glob_err,
 
649
          "Mispriming score is part of objective function, but mispriming library is not defined");
 
650
       return 1;
 
651
    }
 
652
 
 
653
    if(pa->io_weights.repeat_sim && (!pa->io_mishyb_library.seq_num)){
 
654
      pr_append_new_chunk(&pa->glob_err,
 
655
      "Internal oligo mispriming score is part of objective function while mishyb library is not defined");
 
656
      return 1;
 
657
    }
 
658
 
 
659
    if(pa->pr_pair_weights.repeat_sim && (!pa->repeat_lib.seq_num)){
 
660
      pr_append_new_chunk(&pa->glob_err,
 
661
        "Mispriming score is part of objective function, but mispriming library is not defined");
 
662
      return 1;
 
663
    }
 
664
 
 
665
    if(pa->pr_pair_weights.io_quality 
 
666
        && pa->primer_task != pick_pcr_primers_and_hyb_probe ) {
 
667
          pr_append_new_chunk(&pa->glob_err,
 
668
           "Internal oligo quality is part of objective function while internal oligo choice is not required");
 
669
        return 1;
 
670
    }
 
671
 
 
672
    return (NULL == sa->error.data && NULL == pa->glob_err.data) ? 0 : 1;
 
673
}
 
674
 
 
675
/* 
 
676
 * Check intervals, and add to sa->error.  Update the start of each interval to
 
677
 * be relative to the start of the included region.
 
678
 */ 
 
679
static int
 
680
check_intervals(tag_name, num_intervals, intervals, seq_len, sa)
 
681
    const char *tag_name;
 
682
    const int num_intervals;
 
683
    interval_array_t intervals;
 
684
    const int seq_len;
 
685
    seq_args *sa;
 
686
{
 
687
    int i;
 
688
    int outside_warning_issued = 0;
 
689
    for (i=0; i < num_intervals; i++) {
 
690
        if (intervals[i][0] + intervals[i][1] > seq_len) {
 
691
            pr_append_new_chunk(&sa->error, tag_name);
 
692
            pr_append(&sa->error, " beyond end of sequence");
 
693
            return 1;
 
694
        }
 
695
        /* Cause the interval start to be relative to the included region. */
 
696
        intervals[i][0] -= sa->incl_s;
 
697
        /* Check that intervals are within the included region. */
 
698
        if (intervals[i][0] < 0
 
699
            || intervals[i][0] + intervals[i][1] > sa->incl_l) {
 
700
            if (!outside_warning_issued) {
 
701
                pr_append_new_chunk(&sa->warning, tag_name);
 
702
                pr_append(&sa->warning,
 
703
                          " outside of INCLUDED_REGION");
 
704
                outside_warning_issued = 1;
 
705
            }
 
706
        }
 
707
        if (intervals[i][1] < 0) {
 
708
            pr_append_new_chunk(&sa->error, "Negative ");
 
709
            pr_append(&sa->error, tag_name);
 
710
            pr_append(&sa->error, " length");
 
711
            return 1;
 
712
        }
 
713
    }
 
714
    return 0;
 
715
}
 
716
 
 
717
static char * strstr_nocase(s1, s2)
 
718
char *s1, *s2;
 
719
{
 
720
   int  n1, n2;
 
721
   char *p, q, *tmp;
 
722
 
 
723
   if(s1 == NULL || s2 == NULL) return NULL;
 
724
   n1 = strlen(s1); n2 = strlen(s2);
 
725
   if(n1 < n2) return NULL;
 
726
 
 
727
   tmp = pr_safe_malloc(n1 + 1);
 
728
   strcpy(tmp, s1);
 
729
 
 
730
   q = *tmp; p = tmp;
 
731
   while(q != '\0' && q != '\n'){
 
732
      q = *(p + n2);
 
733
      *(p + n2) = '\0';
 
734
      if(strcmp_nocase(p, s2)){
 
735
         *(p + n2) = q; p++; continue;
 
736
      }
 
737
      else {free(tmp); return p;}
 
738
   }
 
739
   free(tmp); return NULL;
 
740
}
 
741
 
 
742
/* Upcase a DNA string, s, in place.  If amibiguity_code_ok is false the
 
743
   string can consist of acgtnACGTN.  If it is true then the IUB/IUPAC
 
744
   ambiguity codes are are allowed.  Return the first unrecognized letter if
 
745
   any is seen (and turn it to 'N' in s).  Otherwise return '\0'.  */
 
746
static char
 
747
dna_to_upper(s, ambiguity_code_ok)
 
748
    char * s;
 
749
    int ambiguity_code_ok;
 
750
{
 
751
  char *p = s;
 
752
  int unrecognized_base = '\0';
 
753
  while (*p) {
 
754
    switch (*p) {
 
755
      case 'a': case 'A': *p='A'; break;
 
756
      case 'c': case 'C': *p='C'; break;
 
757
      case 'g': case 'G': *p='G'; break;
 
758
      case 't': case 'T': *p='T'; break;
 
759
      case 'n': case 'N': *p='N'; break;
 
760
      default: 
 
761
        if (ambiguity_code_ok) {
 
762
          switch (*p) {
 
763
          case 'r': case 'R': *p = 'R'; break;
 
764
          case 'y': case 'Y': *p = 'Y'; break;
 
765
          case 'm': case 'M': *p = 'M'; break;
 
766
          case 'w': case 'W': *p = 'W'; break;
 
767
          case 's': case 'S': *p = 'S'; break;
 
768
          case 'k': case 'K': *p = 'K'; break;
 
769
          case 'd': case 'D': *p = 'D'; break;
 
770
          case 'h': case 'H': *p = 'H'; break;
 
771
          case 'v': case 'V': *p = 'V'; break;
 
772
          case 'b': case 'B': *p = 'B'; break;
 
773
          }
 
774
        } else {
 
775
          if (!unrecognized_base) unrecognized_base = *p;
 
776
          *p = 'N';
 
777
        }
 
778
        break;
 
779
      }
 
780
    p++;
 
781
  }
 
782
  return unrecognized_base;
 
783
}
 
784
 
 
785
int
 
786
_pr_need_template_mispriming(pa)
 
787
  const primer_args *pa;
 
788
{
 
789
  return 
 
790
    pa->max_template_mispriming >= 0
 
791
    || pa->primer_weights.template_mispriming > 0.0
 
792
    || _pr_need_pair_template_mispriming(pa);
 
793
}
 
794
 
 
795
int
 
796
_pr_need_pair_template_mispriming(pa)
 
797
  const primer_args *pa;
 
798
{
 
799
  return 
 
800
    pa->pair_max_template_mispriming >= 0
 
801
    || pa->pr_pair_weights.template_mispriming > 0.0;
 
802
}
 
803
 
 
804
/* Takes substring of seq starting from n with length m and puts it to s.    */
 
805
void
 
806
_pr_substr(const char *seq, int n, int m, char *s)
 
807
{
 
808
        int i;
 
809
        for(i=n;i<n+m;i++)s[i-n]=seq[i];
 
810
        s[m]='\0';
 
811
}
 
812
 
 
813
/* Reverse and complement the sequence seq and put the result in s. */ 
 
814
void
 
815
_pr_reverse_complement(const char *seq, char *s)
 
816
{
 
817
    const char *p = seq;
 
818
    char *q = s;
 
819
 
 
820
    while (*p != '\0') p++;
 
821
    p--;
 
822
    while (p >= seq) {
 
823
        switch (*p)
 
824
        {
 
825
        case 'A': *q='T'; break;
 
826
        case 'C': *q='G'; break;
 
827
        case 'G': *q='C'; break;
 
828
        case 'T': *q='A'; break;
 
829
        case 'U': *q='A'; break;
 
830
 
 
831
        case 'B': *q='V'; break;
 
832
        case 'D': *q='H'; break;
 
833
        case 'H': *q='D'; break;
 
834
        case 'V': *q='B'; break;
 
835
        case 'R': *q='Y'; break;
 
836
        case 'Y': *q='R'; break;
 
837
        case 'K': *q='M'; break;
 
838
        case 'M': *q='K'; break;
 
839
        case 'S': *q='S'; break;
 
840
        case 'W': *q='W'; break;
 
841
 
 
842
        case 'N': *q='N'; break;
 
843
 
 
844
        case 'a': *q='t'; break;
 
845
        case 'c': *q='g'; break;
 
846
        case 'g': *q='c'; break;
 
847
        case 't': *q='a'; break;
 
848
        case 'u': *q='a'; break;
 
849
 
 
850
        case 'b': *q='v'; break;
 
851
        case 'd': *q='h'; break;
 
852
        case 'h': *q='d'; break;
 
853
        case 'v': *q='b'; break;
 
854
        case 'r': *q='y'; break;
 
855
        case 'y': *q='r'; break;
 
856
        case 'k': *q='m'; break;
 
857
        case 'm': *q='k'; break;
 
858
        case 's': *q='s'; break;
 
859
        case 'w': *q='w'; break;
 
860
 
 
861
        case 'n': *q='n'; break;
 
862
        }
 
863
        p--; q++;
 
864
    }
 
865
    *q = '\0';
 
866
}