ThreeB 1.1
sampled_multispot.h
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 #ifndef SPOT_WITH_BACKGROUND_H
00021 #define SPOT_WITH_BACKGROUND_H
00022 
00023 #include <vector>
00024 #include <cvd/image_ref.h>
00025 #include <tr1/tuple>
00026 #include <TooN/TooN.h>
00027 
00028 #include "drift.h"
00029 
00030 typedef char State;
00031 
00032 namespace SampledMultispot
00033 {
00034 
00035 using namespace std;
00036 using namespace CVD;
00037 using namespace TooN;
00038 using namespace std::tr1;
00039 
00040 
00041 
00042 //The changes to SpotWithBackground to operate with/without drift and masking are
00043 //irritatingly pervasive. I can't think of any other way of doing it.
00044 #define SWBG_NAME SpotWithBackground
00045 #include "spot_with_background.hh"
00046 
00047 
00048 #define SWBG_NAME SpotWithBackgroundMasked
00049 #define SWBG_HAVE_MASK
00050 #include "spot_with_background.hh"
00051 
00052 //#define SWBG_NAME SpotWithBackgroundDrift
00053 //#define SWBG_HAVE_DRIFT
00054 //#include "spot_with_background.hh"
00055 
00056 //#undef SWBG_NAME
00057 //#define SWBG_NAME SpotWithBackgroundDriftMasked
00058 //#define SWBG_HAVE_DRIFT
00059 //#define SWBG_HAVE_MASK
00060 //#include "spot_with_background.hh"
00061 
00062 
00063 
00064 
00065 inline double intensity(double i)
00066 {
00067     return i;
00068 }
00069 
00070 inline double intensity(const pair<double, Vector<4> >& i)
00071 {
00072     return i.first;
00073 }
00074 
00075 
00076 //Add and remove a spot over the entire region
00077 template<class T>
00078 void remove_spot(vector<vector<double> >& current_sample_intensities, const vector<T>& spot_intensities, const vector<State>& spot_sample)
00079 {
00080     for(unsigned int frame=0; frame < current_sample_intensities.size(); frame++)
00081         if(spot_sample[frame] == 0) //Spot is on, so remove it
00082             for(unsigned int p=0; p < spot_intensities.size(); p++)
00083                 current_sample_intensities[frame][p] -= intensity(spot_intensities[p]);
00084 }
00085 
00086 template<class T>
00087 void add_spot(vector<vector<double> >& current_sample_intensities, const vector<T>& spot_intensities, const vector<State>& spot_sample)
00088 {
00089     for(unsigned int frame=0; frame < current_sample_intensities.size(); frame++)
00090         if(spot_sample[frame] == 0) //Spot is on, so add it
00091             for(unsigned int p=0; p < spot_intensities.size(); p++)
00092                 current_sample_intensities[frame][p] += intensity(spot_intensities[p]);
00093 }
00094 
00095 
00096 //Add and remove a spot only over a mask. 
00097 template<class T>
00098 void remove_spot(vector<vector<double> >& current_sample_intensities, const vector<T>& spot_intensities, const vector<State>& spot_sample, const vector<int>& mask)
00099 {
00100     for(unsigned int frame=0; frame < current_sample_intensities.size(); frame++)
00101         if(spot_sample[frame] == 0) //Spot is on, so remove it
00102             for(unsigned int p=0; p < mask.size(); p++)
00103                 current_sample_intensities[frame][mask[p]] -= intensity(spot_intensities[mask[p]]);
00104 }
00105 
00106 template<class T>
00107 void add_spot(vector<vector<double> >& current_sample_intensities, const vector<T>& spot_intensities, const vector<State>& spot_sample, const vector<int>& mask)
00108 {
00109     for(unsigned int frame=0; frame < current_sample_intensities.size(); frame++)
00110         if(spot_sample[frame] == 0) //Spot is on, so add it
00111             for(unsigned int p=0; p < mask.size(); p++)
00112                 current_sample_intensities[frame][mask[p]] += intensity(spot_intensities[mask[p]]);
00113 }
00114 
00115 
00116 //Add and remove a drifty spot only over a mask.
00117 template<class T>
00118 void remove_spot(vector<vector<double> >& current_sample_intensities, const vector<vector<T> > & spot_intensities, const vector<State>& spot_sample, const vector<int>& mask)
00119 {
00120     const int steps = spot_intensities.size();
00121     const int frames = current_sample_intensities.size();
00122 
00123     for(int frame=0; frame < frames; frame++)
00124     {
00125         int s = frame * steps / frames;
00126 
00127         if(spot_sample[frame] == 0) //Spot is on, so remove it
00128             for(unsigned int p=0; p < mask.size(); p++)
00129                 current_sample_intensities[frame][mask[p]] -= intensity(spot_intensities[s][mask[p]]);
00130     }
00131 }
00132 
00133 template<class T>
00134 void add_spot(vector<vector<double> >& current_sample_intensities, const vector<vector<T> >& spot_intensities, const vector<State>& spot_sample, const vector<int>& mask)
00135 {
00136     const int steps = spot_intensities.size();
00137     const int frames = current_sample_intensities.size();
00138 
00139     for(int frame=0; frame < frames; frame++)
00140     {
00141         int s = frame * steps / frames;
00142 
00143         if(spot_sample[frame] == 0) //Spot is on, so add it
00144             for(unsigned int p=0; p < mask.size(); p++)
00145                 current_sample_intensities[frame][mask[p]] += intensity(spot_intensities[s][mask[p]]);
00146     }
00147 }
00148 
00149 
00150 
00151 //Compute the spot intensity for a given spot at each pixel
00152 inline vector<double> compute_spot_intensity(const vector<ImageRef>& pixels, const Vector<4>& params)
00153 {
00154     vector<double> intensities(pixels.size());
00155 
00156     for(unsigned int i=0; i < pixels.size(); i++)
00157         intensities[i] = spot_shape(vec(pixels[i]), params);
00158 
00159     return intensities;
00160 }
00161 
00162 //Compute the spot intensity derivatives for a given spot at each pixel
00163 inline vector<pair<double, Vector<4> > > compute_spot_intensity_derivatives(const vector<ImageRef>& pixels, const Vector<4>& params)
00164 {
00165     vector<pair<double, Vector<4> > > derivatives(pixels.size());
00166 
00167     for(unsigned int i=0; i < pixels.size(); i++)
00168         derivatives[i] = spot_shape_diff_position(vec(pixels[i]), params);
00169     return derivatives;
00170 }
00171 
00172 inline vector<tuple<double, Vector<4>, Matrix<4> > > compute_spot_intensity_hessian(const vector<ImageRef>& pixels, const Vector<4>& params)
00173 {
00174     vector<tuple<double, Vector<4>, Matrix<4> > > hessian(pixels.size());
00175 
00176     for(unsigned int i=0; i < pixels.size(); i++)
00177         hessian[i] = spot_shape_hess_position(vec(pixels[i]), params);
00178     return hessian;
00179 }
00180 
00181 
00182 /**
00183 Create a sequence of integers. These can be used as observations
00184 in an observation class by forward_algorithm() and etc.
00185 @param n Length of sequence
00186 @ingroup gUtility
00187 */
00188 inline vector<int> sequence(int n)
00189 {
00190     vector<int> v;
00191     for(int i=0; i < n; i++)
00192         v.push_back(i);
00193     return v;
00194 }
00195 
00196 /*struct RndGrand48
00197 {
00198     double operator()()
00199     {
00200         return drand48();
00201     }
00202 };*/
00203 
00204 ///Draw samples from the spot states given the spots positions and some data.
00205 ///Variable naming matches that in FitSpots.
00206 ///@ingroup gStorm
00207 class GibbsSampler
00208 {
00209     const vector<vector<double> >& pixel_intensities;
00210     const vector<vector<double> >& spot_intensities;
00211     const vector<Vector<4> > spots;
00212     const Matrix<3> A; 
00213     const Vector<3> pi;
00214     const double base_variance;
00215     double variance;
00216 
00217     const int sample_iterations;
00218     const int num_frames, num_pixels;
00219     const vector<int> O;
00220   
00221     vector<vector<State> > current_sample;
00222     vector<vector<double> > current_sample_intensities;
00223 
00224     public:
00225     
00226     GibbsSampler(const vector<vector<double> >& pixel_intensities_,
00227                  const vector<vector<double> >& spot_intensities_,
00228                  const vector<Vector<4> >& spots_,
00229                  const Matrix<3> A_,
00230                  const Vector<3> pi_,
00231                  double variance_,
00232                  int sample_iterations_)
00233     :pixel_intensities(pixel_intensities_), //pixel_intensities: [frame][pixels]
00234      spot_intensities(spot_intensities_),   //spot_intensities: [spot][pixel]
00235      spots(spots_),
00236      A(A_),
00237      pi(pi_),
00238      base_variance(variance_),
00239      variance(variance_),
00240      sample_iterations(sample_iterations_),
00241      num_frames(pixel_intensities.size()),
00242      num_pixels(pixel_intensities[0].size()),
00243         //Observations vector. As usual for this application, the observations are just
00244         //numbered integers which refer to data held elsewhere.
00245      O(sequence(num_frames)),
00246         //Start all spots OFF, so the intensity is 0. OFF is 1 or 2, not 0!!!
00247         //sample_list: [sample][spot][frame]: list of samples drawn using Gibbs sampling
00248      current_sample(spots.size(), vector<State>(num_frames, 2)), //current sample [spot][frame]
00249         //pixel intensities assosciated with the current sample [frame][pixel]
00250      current_sample_intensities(num_frames, vector<double>(num_pixels))
00251     {
00252         //Check a bunch of stuff
00253         assert_same_size(pixel_intensities);
00254         assert_same_size(spot_intensities);
00255 
00256     }
00257 
00258     ///Update the noide variance. Used for adding thermal noise.
00259     ///@param v noise variance.
00260     void set_variance(double v)
00261     {
00262         variance = v;
00263     }
00264 
00265 
00266     ///Reset the gibbs sampler oro the initial state (all spots off)
00267     void reset()
00268     {
00269         vector<State> off(num_frames, 2);
00270         fill(current_sample.begin(), current_sample.end(), off);
00271 
00272         vector<double> black(num_pixels);
00273         fill(current_sample_intensities.begin(), current_sample_intensities.end(), black);
00274         variance = base_variance;
00275     }
00276     
00277     ///Get the next sample
00278     ///@param rng Random number generator
00279     template<class T> void next(T& rng)
00280     {
00281         for(int j=0; j < sample_iterations; j++)
00282             for(int k=0; k < (int) spots.size(); k++)
00283             {
00284                 //Subtract off the spot we're interested in.
00285                 remove_spot(current_sample_intensities, spot_intensities[k], current_sample[k]);
00286 
00287                 //Now current_sample_intensities is the image value for every spot in every frame,
00288                 //except the current spot, which is always set to off. This allows us to add it in 
00289                 //easily.
00290                 SpotWithBackground B(current_sample_intensities, spot_intensities[k], pixel_intensities, variance);
00291                 vector<array<double, 3> > delta = forward_algorithm_delta(A, pi, B, O);
00292                 current_sample[k] = backward_sampling<3,State, T>(A, delta, rng);
00293 
00294                 //Put the newly sampled spot in
00295                 add_spot(current_sample_intensities, spot_intensities[k], current_sample[k]);
00296             }
00297     }
00298     /*void next()
00299     {
00300         RngDrand48 rng;
00301         next(rng);
00302     }*/ 
00303 
00304     ///Retrieve the current sample
00305     const vector<vector<State> >& sample() const
00306     {
00307         return current_sample;
00308     }
00309     ///Retrieve the intensities for the current sample
00310     const vector<vector<double> >& sample_intensities() const
00311     {
00312         return current_sample_intensities;
00313     }
00314 
00315 };
00316 
00317 ///Gibbs sampling class which masks spots to reduce computation.
00318 ///
00319 ///This draws samples from, the spot states given the spots positions and some data. It is
00320 ///very similar to GibbsSampler, except that it only computes probabilities in a mask around each spot
00321 ///to save on computation. Variable naming matches that in FitSpots.
00322 ///@ingroup gStorm
00323 class GibbsSampler2
00324 {
00325     const vector<vector<double> >& pixel_intensities;
00326     const vector<vector<double> >& spot_intensities;
00327     const vector<Vector<4> > spots;
00328     const std::vector<std::vector<int> >& spot_pixels;
00329     const Matrix<3> A; 
00330     const Vector<3> pi;
00331     const double base_variance;
00332     double variance;
00333 
00334     const int sample_iterations;
00335     const int num_frames, num_pixels;
00336     const vector<int> O;
00337   
00338     vector<vector<State> > current_sample;
00339     vector<vector<double> > current_sample_intensities;
00340 
00341     vector<double> cutout_spot_intensities;
00342     vector<vector<double> > cutout_pixel_intensities;
00343     vector<vector<double> > cutout_current_sample_intensities;
00344 
00345     public:
00346     
00347     GibbsSampler2(const vector<vector<double> >& pixel_intensities_,
00348                  const vector<vector<double> >& spot_intensities_,
00349                  const vector<Vector<4> >& spots_,
00350                  const vector<vector<int> >& spot_pixels_,
00351                  const Matrix<3> A_,
00352                  const Vector<3> pi_,
00353                  double variance_,
00354                  int sample_iterations_)
00355     :pixel_intensities(pixel_intensities_), //pixel_intensities: [frame][pixels]
00356      spot_intensities(spot_intensities_),   //spot_intensities: [spot][pixel]
00357      spots(spots_),
00358      spot_pixels(spot_pixels_),
00359      A(A_),
00360      pi(pi_),
00361      base_variance(variance_),
00362      variance(variance_),
00363      sample_iterations(sample_iterations_),
00364      num_frames(pixel_intensities.size()),
00365      num_pixels(pixel_intensities[0].size()),
00366         //Observations vector. As usual for this application, the observations are just
00367         //numbered integers which refer to data held elsewhere.
00368      O(sequence(num_frames)),
00369         //Start all spots OFF, so the intensity is 0. OFF is 1 or 2, not 0!!!
00370         //sample_list: [sample][spot][frame]: list of samples drawn using Gibbs sampling
00371      current_sample(spots.size(), vector<State>(num_frames, 2)), //current sample [spot][frame]
00372         //pixel intensities assosciated with the current sample [frame][pixel]
00373      current_sample_intensities(num_frames, vector<double>(num_pixels)),
00374      cutout_pixel_intensities(num_frames),
00375      cutout_current_sample_intensities(num_frames)
00376     {
00377         //Check a bunch of stuff
00378         assert_same_size(pixel_intensities);
00379         assert_same_size(spot_intensities);
00380     }
00381 
00382     ///Update the noide variance. Used for adding thermal noise.
00383     ///@param v noise variance.
00384     void set_variance(double v)
00385     {
00386         variance = v;
00387     }
00388 
00389 
00390     ///Reset the gibbs sampler oro the initial state (all spots off)
00391     void reset()
00392     {
00393         vector<State> off(num_frames, 2);
00394         fill(current_sample.begin(), current_sample.end(), off);
00395 
00396         vector<double> black(num_pixels);
00397         fill(current_sample_intensities.begin(), current_sample_intensities.end(), black);
00398         variance = base_variance;
00399     }
00400     
00401     ///Get the next sample
00402     ///@param rng Random number generator
00403     template<class T> void next(T& rng)
00404     {
00405     
00406 //double remove=0;
00407 //double cut=0;
00408 //double swb=0;
00409 //double ff_masked=0;
00410 //double bs=0;
00411 //double add=0;
00412 //cvd_timer t;
00413     std::vector<array<double, 3> > delta3;
00414         for(int j=0; j < sample_iterations; j++)
00415             for(int k=0; k < (int) spots.size(); k++)
00416             {
00417 //t.reset();
00418                 //Subtract off the spot we're interested in.
00419                 remove_spot(current_sample_intensities, spot_intensities[k], current_sample[k], spot_pixels[k]);
00420 //remove+=t.reset();
00421 /*
00422                 //Cut out
00423                 //spot
00424                 cutout_spot_intensities.resize(spot_pixels[k].size());
00425                 for(unsigned int i=0; i < spot_pixels[k].size(); i++)
00426                     cutout_spot_intensities[i] = spot_intensities[k][spot_pixels[k][i]];
00427 
00428                 //others
00429                 for(int f=0; f < num_frames; f++)
00430                 {
00431                     cutout_current_sample_intensities[f].resize(spot_pixels[k].size());
00432                     cutout_pixel_intensities[f].resize(spot_pixels[k].size());
00433                     for(unsigned int i=0; i < spot_pixels[k].size();i++)
00434                     {
00435                         cutout_current_sample_intensities[f][i] = current_sample_intensities[f][spot_pixels[k][i]];
00436                         cutout_pixel_intensities[f][i] = pixel_intensities[f][spot_pixels[k][i]];
00437                     }
00438                 }*/
00439 //cut += t.reset();
00440                 //Now current_sample_intensities is the image value for every spot in every frame,
00441                 //except the current spot, which is always set to off. This allows us to add it in 
00442                 //easily.
00443 
00444 //              SpotWithBackground B(current_sample_intensities, spot_intensities[k], pixel_intensities, variance);
00445 //              vector<array<double, 3> > delta = forward_algorithm_delta(A, pi, B, O);
00446 
00447 //ff+=t.reset();
00448 
00449 //              SpotWithBackground B2(cutout_current_sample_intensities, cutout_spot_intensities, cutout_pixel_intensities, variance);
00450 //              std::vector<array<double, 3> > delta2 = forward_algorithm_delta(A, pi, B2, O);
00451 //ff_cut+=t.reset();
00452 
00453                 SpotWithBackgroundMasked B3(current_sample_intensities, spot_intensities[k], pixel_intensities, variance, spot_pixels[k]);
00454 //swb += t.reset();
00455                 forward_algorithm_delta2<3>(A, pi, B3, O, delta3);
00456 //f_masked+=t.reset();
00457                 /*for(unsigned int i=0; i < delta.size(); i++)
00458                 {
00459                     cout.precision(20);
00460                     cout.setf(cout.scientific);
00461                     std::cout << delta[i][0] << " " << delta[i][1] << " " <<delta[i][2] << std::endl;
00462                     std::cout << delta2[i][0] << " " << delta2[i][1] << " " <<delta2[i][2] << std::endl;
00463                     cout << endl;
00464                 }
00465                 std::exit(1);
00466 
00467 */
00468 
00469                 current_sample[k] = backward_sampling<3,State, T>(A, delta3, rng);
00470 //bs += t.reset();
00471                 //Put the newly sampled spot in
00472                 add_spot(current_sample_intensities, spot_intensities[k], current_sample[k], spot_pixels[k]);
00473 //add += t.reset();
00474             }
00475 //      cout << "remove=" <<remove << " cut=" << cut << " swb=" << swb<< " ff_mask=" << ff_masked << " bs=" <<bs << " add="<<add << endl;
00476     }
00477 /*  void next()
00478     {
00479         RngDrand48 rng;
00480         next(rng);
00481     }   
00482 */  
00483     ///Retrieve the current sample
00484     const vector<vector<State> >& sample() const
00485     {
00486         return current_sample;
00487     }
00488     ///Retrieve the intensities for the current sample
00489     const vector<vector<double> >& sample_intensities() const
00490     {
00491         return current_sample_intensities;
00492     }
00493 
00494 };
00495 
00496 #if 0
00497 ///Gibbs sampling class
00498 class GibbsSampler3
00499 {
00500     const vector<vector<double> >& pixel_intensities;
00501     const vector<vector<vector<double> > >& spot_intensities;
00502     const vector<Vector<4> > spots;
00503     const std::vector<std::vector<int> >& spot_pixels;
00504     const Matrix<3> A; 
00505     const Vector<3> pi;
00506     const double base_variance;
00507     double variance;
00508 
00509     const int sample_iterations;
00510     const int num_frames, num_pixels;
00511     const vector<int> O;
00512   
00513     vector<vector<State> > current_sample;
00514     vector<vector<double> > current_sample_intensities;
00515 
00516     public:
00517     
00518     GibbsSampler3(const vector<vector<double> >& pixel_intensities_,
00519                  const vector<vector<vector<double> > >& spot_intensities_,
00520                  const vector<Vector<4> >& spots_,
00521                  const vector<vector<int> >& spot_pixels_,
00522                  const Matrix<3> A_,
00523                  const Vector<3> pi_,
00524                  double variance_,
00525                  int sample_iterations_)
00526     :pixel_intensities(pixel_intensities_), //pixel_intensities: [frame][pixels]
00527      spot_intensities(spot_intensities_),   //spot_intensities: [spot][frame][pixel]
00528      spots(spots_),
00529      spot_pixels(spot_pixels_),
00530      A(A_),
00531      pi(pi_),
00532      base_variance(variance_),
00533      variance(variance_),
00534      sample_iterations(sample_iterations_),
00535      num_frames(pixel_intensities.size()),
00536      num_pixels(pixel_intensities[0].size()),
00537         //Observations vector. As usual for this application, the observations are just
00538         //numbered integers which refer to data held elsewhere.
00539      O(sequence(num_frames)),
00540         //Start all spots OFF, so the intensity is 0. OFF is 1 or 2, not 0!!!
00541         //sample_list: [sample][spot][frame]: list of samples drawn using Gibbs sampling
00542      current_sample(spots.size(), vector<State>(num_frames, 2)), //current sample [spot][frame]
00543         //pixel intensities assosciated with the current sample [frame][pixel]
00544      current_sample_intensities(num_frames, vector<double>(num_pixels))
00545     {
00546         //Check a bunch of stuff
00547         assert_same_size(pixel_intensities);
00548         assert_same_size(spot_intensities);
00549     }
00550     
00551     ///Update the noide variance. Used for adding thermal noise.
00552     ///@param v noise variance.
00553     void set_variance(double v)
00554     {
00555         variance = v;
00556     }
00557 
00558     ///Reset the gibbs sampler oro the initial state (all spots off)
00559     void reset()
00560     {
00561         vector<State> off(num_frames, 2);
00562         fill(current_sample.begin(), current_sample.end(), off);
00563 
00564         vector<double> black(num_pixels);
00565         fill(current_sample_intensities.begin(), current_sample_intensities.end(), black);
00566         variance = base_variance;
00567     }
00568     
00569     ///Get the next sample
00570     ///@param rng Random number generator
00571     template<class T> void next(T& rng)
00572     {
00573     
00574         std::vector<array<double, 3> > delta3;
00575         for(int j=0; j < sample_iterations; j++)
00576             for(int k=0; k < (int) spots.size(); k++)
00577             {
00578                 //Subtract off the spot we're interested in.
00579                 remove_spot(current_sample_intensities, spot_intensities[k], current_sample[k], spot_pixels[k]);
00580 
00581                 //Now current_sample_intensities is the image value for every spot in every frame,
00582                 //except the current spot, which is always set to off. This allows us to add it in 
00583                 //easily.
00584 
00585                 SpotWithBackgroundDriftMasked B3(current_sample_intensities, spot_intensities[k], pixel_intensities, variance, spot_pixels[k]);
00586                 forward_algorithm_delta2<3>(A, pi, B3, O, delta3);
00587 
00588                 current_sample[k] = backward_sampling<3,State, T>(A, delta3, rng);
00589                 //Put the newly sampled spot in
00590                 add_spot(current_sample_intensities, spot_intensities[k], current_sample[k], spot_pixels[k]);
00591             }
00592     }
00593 /*  void next()
00594     {
00595         RngDrand48 rng;
00596         next(rng);
00597     }   
00598 */
00599     ///Retrieve the current sample
00600     const vector<vector<State> >& sample() const
00601     {
00602         return current_sample;
00603     }
00604     ///Retrieve the intensities for the current sample
00605     const vector<vector<double> >& sample_intensities() const
00606     {
00607         return current_sample_intensities;
00608     }
00609 
00610 };
00611 #endif 
00612 
00613 }
00614 
00615 using SampledMultispot::SpotWithBackground;
00616 using SampledMultispot::remove_spot;
00617 using SampledMultispot::add_spot;
00618 using SampledMultispot::compute_spot_intensity;
00619 using SampledMultispot::compute_spot_intensity_hessian;
00620 using SampledMultispot::compute_spot_intensity_derivatives;
00621 using SampledMultispot::sequence;
00622 using SampledMultispot::GibbsSampler;
00623 using SampledMultispot::GibbsSampler2;
00624 
00625 #endif