source: trunk/SEQ_QUALITY/SQ_main.cxx

Last change on this file was 19425, checked in by westram, 17 months ago
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 12.7 KB
Line 
1//  ==================================================================== //
2//                                                                       //
3//    File      : SQ_main.cxx                                            //
4//    Purpose   : Entrypoint to Seq. Quality analysis; calls functions   //
5//                                                                       //
6//                                                                       //
7//  Coded by Juergen Huber in July 2003 - February 2004                  //
8//  Coded by Kai Bader (baderk@in.tum.de) in 2007 - 2008                 //
9//  Copyright Department of Microbiology (Technical University Munich)   //
10//                                                                       //
11//  Visit our web site at: http://www.arb-home.de/                       //
12//                                                                       //
13//  ==================================================================== //
14
15#include "seq_quality.h"
16#include "SQ_functions.h"
17
18#include <awt_filter.hxx>
19#include <awt_sel_boxes.hxx>
20
21#include <aw_awars.hxx>
22#include <aw_msg.hxx>
23#include <aw_root.hxx>
24
25#include <arb_progress.h>
26#include <TreeNode.h>
27#include <arb_global_defs.h>
28#include <awt_config_manager.hxx>
29
30// --------------------------------------------------------------------------------
31
32#define AWAR_SQ_PERM "seq_quality/"     // saved in properties
33#define AWAR_SQ_TEMP "tmp/seq_quality/" // not saved in properties
34
35#define AWAR_SQ_WEIGHT_BASES     AWAR_SQ_PERM "weight_bases"
36#define AWAR_SQ_WEIGHT_DEVIATION AWAR_SQ_PERM "weight_deviation"
37#define AWAR_SQ_WEIGHT_HELIX     AWAR_SQ_PERM "weight_helix"
38#define AWAR_SQ_WEIGHT_CONSENSUS AWAR_SQ_PERM "weight_consensus"
39#define AWAR_SQ_WEIGHT_IUPAC     AWAR_SQ_PERM "weight_iupac"
40#define AWAR_SQ_WEIGHT_GC        AWAR_SQ_PERM "weight_gc"
41
42#define AWAR_SQ_MARK_ONLY_FLAG AWAR_SQ_PERM "mark_only_flag"
43#define AWAR_SQ_MARK_FLAG      AWAR_SQ_PERM "mark_flag"
44#define AWAR_SQ_MARK_BELOW     AWAR_SQ_PERM "mark_below"
45#define AWAR_SQ_REEVALUATE     AWAR_SQ_PERM "reevaluate"
46#define AWAR_SQ_FILTER_NAME    AWAR_SQ_TEMP "filter/name"
47
48void SQ_create_awars(AW_root *aw_root, AW_default aw_def) {
49    aw_root->awar_int(AWAR_SQ_WEIGHT_BASES,     5,  aw_def);
50    aw_root->awar_int(AWAR_SQ_WEIGHT_DEVIATION, 15, aw_def);
51    aw_root->awar_int(AWAR_SQ_WEIGHT_HELIX,     15, aw_def);
52    aw_root->awar_int(AWAR_SQ_WEIGHT_CONSENSUS, 50, aw_def);
53    aw_root->awar_int(AWAR_SQ_WEIGHT_IUPAC,     5,  aw_def);
54    aw_root->awar_int(AWAR_SQ_WEIGHT_GC,        10, aw_def);
55    aw_root->awar_int(AWAR_SQ_MARK_ONLY_FLAG,   0,  aw_def);
56    aw_root->awar_int(AWAR_SQ_MARK_FLAG,        1,  aw_def);
57    aw_root->awar_int(AWAR_SQ_MARK_BELOW,       40, aw_def);
58    aw_root->awar_int(AWAR_SQ_REEVALUATE,       0,  aw_def);
59
60    awt_create_filter_awars(aw_root, aw_def, AWAR_SQ_FILTER_NAME, AWAR_DEFAULT_ALIGNMENT);
61}
62
63// --------------------------------------------------------------------------------
64
65
66inline size_t count_nodes(TreeNode *node) {
67    // calculate number of nodes in tree
68    return GBT_count_leafs(node)*2-1;
69}
70
71static void sq_calc_seq_quality_cb(AW_window *aww, adfiltercbstruct *acbs, GBDATA *gb_main) {
72    AW_root        *aw_root     = aww->get_root();
73    GB_ERROR        error       = NULp;
74    TreeNode       *tree        = NULp;
75    bool            marked_only = (aw_root->awar(AWAR_SQ_MARK_ONLY_FLAG)->read_int() > 0);
76    arb_progress    main_progress("Calculating sequence quality");
77    GB_transaction  ta(gb_main);
78
79    if (!error) {
80        char *ali = GBT_get_default_alignment(gb_main);
81        if (ali) {
82            GB_alignment_type type = GBT_get_alignment_type(gb_main, ali);
83            if (type == GB_AT_UNKNOWN) {
84                error = GB_await_error();
85            }
86            else if (type != GB_AT_RNA && type != GB_AT_DNA) {
87                error = "does only support RNA/DNA sequence data";
88            }
89            if (error) {
90                error = GBS_global_string("Invalid alignment '%s' selected (%s)", ali, error);
91            }
92        }
93        else {
94            error = "No default alignment selected";
95        }
96        free(ali);
97    }
98
99    if (!error) {
100        char *treename = aw_root->awar(AWAR_TREE)->read_string();
101
102        if (treename && strcmp(treename, NO_TREE_SELECTED) != 0) {
103            if (!error) {
104                tree = GBT_read_tree(gb_main, treename, new SimpleRoot);
105                if (!tree) error = GB_await_error();
106                else {
107                    error = GBT_link_tree(tree, gb_main, false, NULp, NULp);
108                    if (!error) {
109                        GBT_TreeRemoveType mode = marked_only ? GBT_KEEP_MARKED : GBT_REMOVE_ZOMBIES;
110                        tree = GBT_remove_leafs(tree, mode, NULp, NULp, NULp);
111                        if (!tree || tree->is_leaf()) {
112                            error = GBS_global_string("Tree contains less than 2 species after removing zombies%s",
113                                                      marked_only ? " and non-marked" : "");
114                        }
115                    }
116                }
117            }
118        }
119        free(treename);
120    }
121
122    // if tree == 0 -> do basic quality calculations that are possible without tree information
123    // otherwise    -> use all groups found in tree and compare sequences against the groups they are contained in
124
125    if (!error) {
126        struct SQ_weights weights;
127
128        weights.bases             = aw_root->awar(AWAR_SQ_WEIGHT_BASES)->read_int();
129        weights.diff_from_average = aw_root->awar(AWAR_SQ_WEIGHT_DEVIATION)->read_int();
130        weights.helix             = aw_root->awar(AWAR_SQ_WEIGHT_HELIX)->read_int();
131        weights.consensus         = aw_root->awar(AWAR_SQ_WEIGHT_CONSENSUS)->read_int();
132        weights.iupac             = aw_root->awar(AWAR_SQ_WEIGHT_IUPAC)->read_int();
133        weights.gc                = aw_root->awar(AWAR_SQ_WEIGHT_GC)->read_int();
134
135        int mark_flag  = aw_root->awar(AWAR_SQ_MARK_FLAG)->read_int();
136        int mark_below = aw_root->awar(AWAR_SQ_MARK_BELOW)->read_int();
137        int reevaluate = aw_root->awar(AWAR_SQ_REEVALUATE)->read_int();
138
139        // Load and use Sequence-Filter
140        AP_filter *filter = awt_get_filter(acbs);
141        error             = awt_invalid_filter(filter);
142
143        /*
144          SQ_evaluate() generates the final estimation for the quality of an alignment.
145          It takes the values from the different containers, which are generated by the other functions, weights them
146          and calculates a final value. The final value is stored in "value_of_evaluation" (see options).
147          With the values stored in "weights" one can customize how important a value stored in a container becomes
148          for the final result.
149        */
150
151        if (!error) {
152            if (!tree) {
153                if (reevaluate) {
154                    error = SQ_mark_species(gb_main, mark_below, marked_only);
155                }
156                else {
157                    arb_progress  progress(GBT_get_species_count(gb_main)*2);
158                    SQ_GroupData *globalData = new SQ_GroupData_RNA;
159
160                    progress.subtitle("pass1");
161                    error = SQ_pass1_no_tree(globalData, gb_main, filter, progress);
162                    if (!error) {
163                        progress.subtitle("pass2");
164                        error = SQ_pass2_no_tree(globalData, gb_main, filter, progress);
165                        if (!error) {
166                            error = SQ_evaluate(gb_main, weights, marked_only);
167                            if (mark_flag && !error) {
168                                error = SQ_mark_species(gb_main, mark_below, marked_only);
169                            }
170                        }
171                    }
172                    if (error) progress.done();
173                    delete globalData;
174                }
175            }
176            else {
177                SQ_TREE_ERROR check = SQ_check_tree_structure(tree);
178                if (check != NONE) {
179                    switch (check) {
180                        case ZOMBIE:
181                            error = "Found one or more zombies in the tree.\n"
182                                "Please remove them or use another tree before running the quality check tool.";
183                            break;
184                        case MISSING_NODE:
185                            error = "Missing node(s) or unusable tree structure.\n"
186                                "Please fix the tree before running the quality check tool.";
187                            break;
188                        default:
189                            error = "An error occurred while traversing the tree.\n"
190                                "Please fix the tree before running the quality check tool.";
191                            break;
192                    }
193                }
194                else if (reevaluate) {
195                    error = SQ_mark_species(gb_main, mark_below, marked_only);
196                }
197                else {
198                    arb_progress progress(count_nodes(tree)*2);
199                    SQ_GroupData *globalData = new SQ_GroupData_RNA;
200
201                    progress.subtitle("pass1");
202                    error = SQ_pass1_on_tree(tree, gb_main, globalData, filter, progress);
203                    if (!error) {
204                        progress.subtitle("pass2");
205                        error = SQ_pass2_on_tree(tree, gb_main, globalData, filter, progress);
206                        if (!error) {
207                            error = SQ_evaluate(gb_main, weights, marked_only);
208                            if (mark_flag && !error) {
209                                error = SQ_mark_species(gb_main, mark_below, marked_only);
210                            }
211                        }
212                    }
213                    if (error) progress.done();
214                    delete globalData;
215                }
216            }
217        }
218        awt_destroy_filter(filter);
219    }
220
221    error = ta.close(error);
222    if (error) aw_message(error);
223
224    SQ_clear_group_dictionary();
225    UNCOVERED();
226    destroy(tree);
227}
228
229static void sq_remove_quality_entries_cb(AW_window*, GBDATA *gb_main) {
230    GB_ERROR error = SQ_remove_quality_entries(gb_main);
231    aw_message_if(error);
232}
233
234static AWT_config_mapping_def seq_quality_config_mapping[] = {
235    { AWAR_SQ_WEIGHT_BASES,     "wbases" },
236    { AWAR_SQ_WEIGHT_DEVIATION, "wdeviation" },
237    { AWAR_SQ_WEIGHT_HELIX,     "whelix" },
238    { AWAR_SQ_WEIGHT_CONSENSUS, "wconsens" },
239    { AWAR_SQ_WEIGHT_IUPAC,     "wiupac" },
240    { AWAR_SQ_WEIGHT_GC,        "wgc" },
241    { AWAR_SQ_MARK_ONLY_FLAG,   "onlymarked" },
242    { AWAR_SQ_MARK_FLAG,        "markbad" },
243    { AWAR_SQ_MARK_BELOW,       "markbelow" },
244    { AWAR_SQ_REEVALUATE,       "reeval" },
245
246    { NULp, NULp }
247};
248
249AW_window *SQ_create_seq_quality_window(AW_root *aw_root, GBDATA *gb_main) {
250    // create window for sequence quality calculation (called only once)
251
252    AW_window_simple *aws = new AW_window_simple;
253
254    aws->init(aw_root, "CALC_SEQ_QUALITY", "CALCULATE SEQUENCE QUALITY");
255    aws->load_xfig("seq_quality.fig");
256
257    aws->at("close");
258    aws->callback(AW_POPDOWN);
259    aws->create_button("CLOSE", "CLOSE", "C");
260
261    aws->at("help");
262    aws->callback(makeHelpCallback("seq_quality.hlp"));
263    aws->create_button("HELP", "HELP", "H");
264
265    aws->at("base");
266    aws->create_input_field(AWAR_SQ_WEIGHT_BASES, 3);
267
268    aws->at("deviation");
269    aws->create_input_field(AWAR_SQ_WEIGHT_DEVIATION, 3);
270
271    aws->at("no_helices");
272    aws->create_input_field(AWAR_SQ_WEIGHT_HELIX, 3);
273
274    aws->at("consensus");
275    aws->create_input_field(AWAR_SQ_WEIGHT_CONSENSUS, 3);
276
277    aws->at("iupac");
278    aws->create_input_field(AWAR_SQ_WEIGHT_IUPAC, 3);
279
280    aws->at("gc_proportion");
281    aws->create_input_field(AWAR_SQ_WEIGHT_GC, 3);
282
283    aws->at("monly");
284    aws->create_toggle(AWAR_SQ_MARK_ONLY_FLAG);
285
286    aws->at("mark");
287    aws->create_toggle(AWAR_SQ_MARK_FLAG);
288
289    aws->at("mark_below");
290    aws->create_input_field(AWAR_SQ_MARK_BELOW, 3);
291
292    aws->at("tree");
293    awt_create_TREE_selection_list(gb_main, aws, AWAR_TREE);
294
295    aws->at("filter");
296    adfiltercbstruct *adfilter = awt_create_select_filter(aws->get_root(), gb_main, AWAR_SQ_FILTER_NAME);
297    aws->callback(makeCreateWindowCallback(awt_create_select_filter_win, adfilter));
298    aws->create_button("SELECT_FILTER", AWAR_SQ_FILTER_NAME);
299
300    aws->at("go");
301    aws->callback(makeWindowCallback(sq_calc_seq_quality_cb, adfilter, gb_main));
302    aws->highlight();
303    aws->create_button("GO", "GO", "G");
304
305    aws->at("config");
306    AWT_insert_config_manager(aws, AW_ROOT_DEFAULT, "seq_quality", seq_quality_config_mapping);
307
308    aws->at("reevaluate");
309    aws->label("Re-Evaluate only");
310    aws->create_toggle(AWAR_SQ_REEVALUATE);
311
312    aws->at("remove");
313    aws->callback(makeWindowCallback(sq_remove_quality_entries_cb, gb_main));
314    aws->create_button("Remove", "Remove", "R");
315
316    return aws;
317}
Note: See TracBrowser for help on using the repository browser.