ThreeB 1.1

- Main Page
- Related Pages
- Modules
- Namespaces
- Classes
- Files

- File List
- File Members

multispot5_headless.cc
Go to the documentation of this file.

    00001 /*
    00002     This file is part of B-cubed.
    00003
    00004     Copyright (C) 2009, 2010, 2011, Edward Rosten and Susan Cox
    00005
    00006     B-cubed is free software; you can redistribute it and/or
    00007     modify it under the terms of the GNU Lesser General Public
    00008     License as published by the Free Software Foundation; either
    00009     version 3.0 of the License, or (at your option) any later version.
    00010
    00011     This library is distributed in the hope that it will be useful,
    00012     but WITHOUT ANY WARRANTY; without even the implied warranty of
    00013     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
    00014     Lesser General Public License for more details.
    00015
    00016     You should have received a copy of the GNU General Public License
    00017     along with this program.  If not, see <http://www.gnu.org/licenses/>
    00018 */
    00019
    00020 #include <tag/printf.h>
    00021 #undef make_tuple
    00022
    00023 #include <tr1/tuple>
    00024 #include <algorithm>
    00025 #include <climits>
    00026 #include <iomanip>
    00027 #include <map>
    00028 #include <cvd/image_io.h>
    00029 #include <cvd/image_convert.h>
    00030 #include <cvd/morphology.h>
    00031 #include <cvd/connected_components.h>
    00032 #include <cvd/draw.h>
    00033 #include <cvd/vector_image_ref.h>
    00034
    00035 #include <gvars3/instances.h>
    00036
    00037 #include "storm_imagery.h"
    00038 #include "multispot5.h"
    00039 #include "multispot5_place_choice.h"
    00040 #include "utility.h"
    00041
    00042 using namespace std;
    00043 using namespace std::tr1;
    00044 using namespace CVD;
    00045 using namespace GVars3;
    00046 using namespace TooN;
    00047
    00048
    00049
    00050 vector<vector<ImageRef> > get_regions(const SubImage<double>& log_ratios)
    00051 {
    00052     gvar3<double> radius("radius", 0, 1);
    00053
    00054     //Set the liklihood ratio threshold/spot density prior
    00055     //same thing.
    00056     double threshold = GV3::get<double>("threshold", 0, -1);
    00057     int edge = GV3::get<int>("edge", 0, -1);
    00058
    00059
    00060     //Threshold image
    00061     Image<byte> thresholded(log_ratios.size(), 0);
    00062     for(int r=0; r < thresholded.size().y; r++)
    00063         for(int c=0; c < min(thresholded.size().x, edge); c++)
    00064             thresholded[r][c] = 255 * (log_ratios[r][c] > threshold);
    00065
    00066     //Dilate
    00067     Image<byte> dilated = morphology(thresholded, getDisc(*radius), Morphology::BinaryDilate<byte>());
    00068
    00069     transform(dilated.begin(), dilated.end(), dilated.begin(), bind1st(multiplies<int>(), 255));
    00070
    00071     //Connected components of dilated image
    00072     vector<ImageRef> fg;
    00073     for(int r=0; r < thresholded.size().y; r++)
    00074         for(int c=0; c < min(thresholded.size().x, edge); c++)
    00075             if(dilated[r][c])
    00076                 fg.push_back(ImageRef(c, r));
    00077
    00078     vector<vector<ImageRef> > regions;
    00079     connected_components(fg, regions);
    00080
    00081     return regions;
    00082 }
    00083
    00084 void mmain(int argc, char** argv)
    00085 {
    00086     GUI.LoadFile("multispot5.cfg");
    00087     int lastarg = GUI.parseArguments(argc, argv);
    00088     if(lastarg >= argc)
    00089     {
    00090         cerr << "Specify the images to load\n";
    00091         exit(1);
    00092     }
    00093     vector<string> files(argv + lastarg, argv + argc);
    00094
    00095     //Save this now since the de-checkpointing code will kl0bber it
    00096     //when it reloads the gvars
    00097     string save_spots_file = GV3::get<string>("save_spots", "", -1);
    00098
    00099     string checkpoint_file = GV3::get<string>("load_checkpoint", "", 1);
    00100
    00101     if(checkpoint_file != "")
    00102     {
    00103         //Load and de-checkpointing
    00104         ifstream chk;
    00105         open_or_die(chk, checkpoint_file);
    00106
    00107         StateParameters p;
    00108
    00109         try{
    00110             p = parse_log_file(chk);
    00111         }
    00112         catch(LogFileParseError e)
    00113         {
    00114             cerr << "SI TEH FUX0R11ONEone!oneleven: " << e.what << endl;
    00115             exit(1);
    00116         }
    00117
    00118         vector<Image<float> > ims = load_and_normalize_images(files);
    00119
    00120         //Restore kl0bbered variable
    00121         GV3::get<string>("save_spots") = save_spots_file;
    00122
    00123         ofstream save_spots;
    00124         open_or_die(save_spots, save_spots_file);
    00125
    00126         fit_spots_new(ims, p, save_spots, *null_graphics());
    00127
    00128     }
    00129
    00130     vector<Image<float> > ims = load_and_normalize_images(files);
    00131
    00132     //Load the log_ratios image.
    00133     //We will use this as a starting point for searching for spots.
    00134     Image<double> log_ratios;
    00135     try
    00136     {
    00137         log_ratios = img_load(GV3::get<string>("log_ratios", "", -1));
    00138     }
    00139     catch(Exceptions::All e)
    00140     {
    00141         cerr << "Error loading " << GV3::get<string>("log_ratios", "") << ": " << e.what << endl;
    00142         exit(1);
    00143     }
    00144
    00145
    00146     gvar3<int> cluster_to_show("cluster_to_show", 0, -1);
    00147     gvar3<int> use_largest("use_largest", 0, 1);
    00148
    00149     vector<vector<ImageRef> > regions;
    00150
    00151     regions = get_regions(log_ratios);
    00152     if(regions.size() == 0)
    00153     {
    00154         cerr << "There are no regions!\n";
    00155
    00156         ofstream save_spots;
    00157         open_or_die(save_spots, save_spots_file);
    00158         save_spots << "NOREGIONS\n";
    00159
    00160         exit(1);
    00161     }
    00162
    00163     if(*use_largest && !regions.empty())
    00164     {
    00165         *cluster_to_show=0;
    00166         for(unsigned int i=1; i < regions.size(); i++)
    00167             if(regions[i].size() > regions[*cluster_to_show].size())
    00168                 *cluster_to_show = i;
    00169
    00170     }
    00171     else
    00172         *cluster_to_show = max(min(*cluster_to_show, (int)regions.size() - 1), 0);
    00173
    00174
    00175     auto_ptr<FitSpotsGraphics> gr = null_graphics();
    00176     place_and_fit_spots(ims, regions[*cluster_to_show], log_ratios, save_spots_file, *gr);
    00177 }
    00178
    00179
    00180 int main(int argc, char** argv)
    00181 {
    00182     try{
    00183         mmain(argc, argv);
    00184     }
    00185     catch(Exceptions::All e)
    00186     {
    00187         cerr << "Fatal error: " << e.what << endl;
    00188     }
    00189 }