Classes | |
struct | ConjugateGradientOnly< Size, Precision > |
Class for performing optimization with Conjugate Gradient, where only the derivatives are available. More... | |
class | NullGraphics |
Graphics class which does absoloutely nothing. More... | |
class | DataForMCMC |
Closure hoding the data required do use GibbsSampler2 See FitSpots for naming of variables. More... | |
struct | SampledBackgroundData |
Closure holding image data generated using samples drawn from the model. More... | |
struct | SpotNegProbabilityDiffWithSampledBackground |
Compute the derivative of the negative log probability with respect to the parameters of one spot, given some samples of the other spots. More... | |
class | FreeEnergyHessian |
Class for computing the Hessian of the negative free energy. More... | |
class | FitSpots |
Mega class which actually does the meat of the spot fitting. More... | |
class | SampledMultispot::GibbsSampler |
Draw samples from the spot states given the spots positions and some data. More... | |
class | SampledMultispot::GibbsSampler2 |
Gibbs sampling class which masks spots to reduce computation. More... | |
Functions | |
auto_ptr< FitSpotsGraphics > | null_graphics () |
void | get_spot_pixels (const vector< ImageRef > &pixels, const Vector< 4 > &spot, vector< int > &out) |
StateParameters | parse_log_file (istream &in) |
double | brightness_motion_limit (double mu, double sigma, bool not_one) |
void | fit_spots_new (const vector< Image< float > > &ims, StateParameters &p, ofstream &save_spots, FitSpotsGraphics &gr) |
template<class B > | |
double | spot_shape_s (const TooN::Vector< 2 > &x, const TooN::Vector< 4, double, B > &phi) |
template<class B > | |
std::pair< double, TooN::Vector< 4 > > | spot_shape_diff_position (const TooN::Vector< 2 > &x, const TooN::Vector< 4, double, B > &phi) |
template<class B > | |
std::tr1::tuple< double, TooN::Vector < 4 >, TooN::Matrix< 4 > > | spot_shape_hess_position (const TooN::Vector< 2 > &x, const TooN::Vector< 4, double, B > &phi) |
template<class B > | |
std::tr1::tuple< double, TooN::Vector < 2 >, TooN::Matrix< 2 > > | spot_shape_hess (const TooN::Vector< 2 > &x, const TooN::Vector< 4, double, B > &phi) |
template<class B > | |
std::pair< double, TooN::Vector< 2 > > | spot_shape_diff (const TooN::Vector< 2 > &x, const TooN::Vector< 4, double, B > &phi) |
template<class B > | |
double | spot_shape (const TooN::Vector< 2 > &x, const TooN::Vector< 4, double, B > &phi) |
template<class Base > | |
std::tr1::tuple< double, TooN::Vector < 2 >, TooN::Matrix< 2 > > | log_probability_spot_hess (const CVD::SubImage< float > &im, double variance, const TooN::Vector< 4, double, Base > &spot_parameters) |
template<class Base > | |
std::pair< double, TooN::Vector< 2 > > | log_probability_spot_diff (const CVD::SubImage< float > &im, double variance, const TooN::Vector< 4, double, Base > &spot_parameters) |
template<class Base > | |
double | log_probability_spot (const CVD::SubImage< float > &im, double variance, const TooN::Vector< 4, double, Base > &spot_parameters) |
double | log_normal_std (double mu, double sigma) |
double | log_normal_mode (double mu, double sigma) |
double | log_log_normal (double x, double mu, double sigma) |
double | diff_log_log_normal (double x, double mu, double sigma) |
double | hess_log_log_normal (double x, double mu, double sigma) |
auto_ptr<FitSpotsGraphics> null_graphics | ( | ) |
Factory function to generate an instance of NullGraphics.
Definition at line 86 of file multispot5.cc.
Referenced by mmain().
00087 { 00088 return auto_ptr<FitSpotsGraphics>(new NullGraphics); 00089 }
void get_spot_pixels | ( | const vector< ImageRef > & | pixels, | |
const Vector< 4 > & | spot, | |||
vector< int > & | out | |||
) |
Which pixels belong to a given spot? Find the indices of those pixels.
Definition at line 884 of file multispot5.cc.
Referenced by fit_spots_historic(), FitSpots::optimize_each_spot_in_turn_for_several_passes(), and FitSpots::try_modifying_model().
00885 { 00886 //Go out to three sigma 00887 00888 vector<ImageRef> pix = getDisc(spot[1]*6 + 1); 00889 out.resize(0); 00890 ImageRef offset = ir_rounded(spot.slice<2,2>()); 00891 for(unsigned int j=0; j < pix.size(); j++) 00892 { 00893 int pos = lower_bound(pixels.begin(), pixels.end(), pix[j] + offset) - pixels.begin(); 00894 if(pos != (int)pixels.size() && pixels[pos] == pix[j] + offset) 00895 out.push_back(pos); 00896 } 00897 00898 if(out.size() == 0) 00899 { 00900 cout << "********************************\n"; 00901 cout << "********************************\n"; 00902 cout << "********************************\n"; 00903 cout << "********************************\n"; 00904 cout << "********************************\n"; 00905 cout << "Oe noes!11one\n"; 00906 cout << pix.size() << endl; 00907 } 00908 }
StateParameters parse_log_file | ( | istream & | in | ) |
Parser for multispot 5 log files.
Log files are mostly line oriented and contain various records
The main records are:
Iteraton: #ITNUM MAIN: <list of spot parameters>
Pass: #PASSNUM MT19937 <random number generator state> PASS#PASSNUM: <list of spot parameters> ENDCHECKPOINT
Note that MAIN is redundant since it will be the same as the following PASS 1 (or the first pass computed if restoring from a checkpoint).
Data should only be considered valid after ENDCHECKPOINT has been read
Iteration is written once per iteration, not once per pass. (FIXME)
Which moron invented this file format???
Note that the file format hasn't beren fixed, so that the output can easily be compared to the output of the historic version which is known to be good.
in | Stream to parse file from |
Definition at line 984 of file multispot5.cc.
References StateParameters::iteration, StateParameters::pass, StateParameters::pixels, StateParameters::rng, split(), StateParameters::spots, and xtoa().
Referenced by mmain().
00985 { 00986 //A line read from the file 00987 string line; 00988 00989 //State lines known to be OK 00990 string rngline, passline, iterationline; 00991 bool state_ok=0; 00992 00993 //State lines read in, with flags of goodness 00994 string new_rngline, new_passline, new_iterationline; 00995 bool new_rngline_ok=0, new_passline_ok=0, new_iterationline_ok=0; 00996 00997 unsigned int lineno=0; 00998 bool doing_gvars = 0; 00999 01000 vector<ImageRef> pixels; 01001 01002 while(!in.eof()) 01003 { 01004 getline(in, line); 01005 if(in.fail()) 01006 break; 01007 01008 lineno++; 01009 01010 if(line == "ENDGVARLIST") 01011 { 01012 if(!doing_gvars) 01013 throw LogFileParseError("Spurious end of GVars"); 01014 doing_gvars = 0; 01015 } 01016 else if(doing_gvars) 01017 { 01018 GUI.ParseLine(line); 01019 } 01020 else if(line == "BEGINGVARLIST") 01021 { 01022 doing_gvars = 1; 01023 } 01024 if(line.substr(0, 11) == "Iteration: ") 01025 { 01026 new_iterationline = line; 01027 new_iterationline_ok = true; 01028 } 01029 else if(line.substr(0, 4) == "PASS") 01030 { 01031 new_passline = line; 01032 if(new_passline_ok) 01033 throw LogFileParseError("Duplicate PASS on line " + xtoa(lineno)); 01034 new_passline_ok = true; 01035 } 01036 else if(line.substr(0, 8) == "MT19937 ") 01037 { 01038 new_rngline = line; 01039 if(new_rngline_ok) 01040 throw LogFileParseError("Duplicate MT19937 on line " + xtoa(lineno)); 01041 01042 new_rngline_ok = true; 01043 01044 } 01045 else if(line == "ENDCHECKPOINT") 01046 { 01047 if(new_passline_ok && new_rngline_ok && new_iterationline_ok) 01048 { 01049 iterationline = new_iterationline; 01050 rngline = new_rngline; 01051 passline = new_passline; 01052 } 01053 else 01054 throw LogFileParseError("Reached checkpoint with missing data: " 01055 "it=" + xtoa(new_iterationline_ok) + 01056 " pa=" + xtoa(new_passline_ok) + 01057 " rg=" + xtoa(new_rngline_ok) + " on line " + xtoa(lineno)); 01058 01059 //Don't reset iteration since it only appears once for each entire 01060 //set of passes. 01061 new_rngline_ok = 0; 01062 new_passline_ok = 0; 01063 01064 state_ok = true; 01065 } 01066 else if(line.substr(0, 7) == "PIXELS ") 01067 { 01068 vector<string> l = split(line); 01069 if( (l.size() - 1)%2 == 0) 01070 { 01071 int n = (l.size()-1)/2; 01072 pixels.resize(n); 01073 for(int i=0; i < n; i++) 01074 { 01075 pixels[i].x = atox<int>(l[i*2+1+0], "pixels"); 01076 pixels[i].y = atox<int>(l[i*2+1+1], "pixels"); 01077 } 01078 } 01079 else 01080 throw LogFileParseError("Bad PIXELS line"); 01081 } 01082 } 01083 01084 if(!state_ok) 01085 throw LogFileParseError("No state found"); 01086 01087 if(pixels.size() == 0) 01088 throw LogFileParseError("No pixels, or pixels is empty"); 01089 01090 //Now parse the lines 01091 StateParameters p; 01092 vector<string> l; 01093 01094 //Parse the iterations 01095 l = split(iterationline); 01096 p.iteration = atox<int>(l[1], "iteration"); 01097 01098 //Parse the random number generator 01099 p.rng =shared_ptr<MT19937>(new MT19937); 01100 { 01101 istringstream rng_s(rngline); 01102 try{ 01103 p.rng->read(rng_s); 01104 } 01105 catch(MT19937::ParseError p) 01106 { 01107 throw LogFileParseError("Error parsing MT19937"); 01108 } 01109 } 01110 01111 //Parse PASS and the listing of spots 01112 l = split(passline); 01113 if( (l.size() - 1 ) % 4 == 0) 01114 { 01115 p.pass = atox<int>(l[0].substr(4), "pass"); 01116 01117 for(unsigned int i=0; i < (l.size()-1)/4; i++) 01118 { 01119 cout << l[i*4+1+0] << endl; 01120 cout << l[i*4+1+1] << endl; 01121 cout << l[i*4+1+2] << endl; 01122 cout << l[i*4+1+3] << endl; 01123 p.spots.push_back(makeVector( 01124 atox<double>(l[i*4+1+0], "spot"), 01125 atox<double>(l[i*4+1+1], "spot"), 01126 atox<double>(l[i*4+1+2], "spot"), 01127 atox<double>(l[i*4+1+3], "spot"))); 01128 } 01129 01130 } 01131 else 01132 throw LogFileParseError("Wrong number of elements in PASS line"); 01133 01134 //Set up the pixels (oh god the pixels) 01135 p.pixels = pixels; 01136 01137 return p; 01138 }
double brightness_motion_limit | ( | double | mu, | |
double | sigma, | |||
bool | not_one | |||
) |
How far should steps in brightness be limited to?
Definition at line 1332 of file multispot5.cc.
References log_normal_std().
01333 { 01334 if(not_one) 01335 return log_normal_std(mu, sigma); 01336 else 01337 return 1; 01338 }
void fit_spots_new | ( | const vector< Image< float > > & | ims, | |
StateParameters & | p, | |||
ofstream & | save_spots, | |||
FitSpotsGraphics & | gr | |||
) |
Wrapper function for using FitSpots.
Definition at line 2836 of file multispot5.cc.
References FitSpots::run().
Referenced by mmain().
02837 { 02838 FitSpots fit(ims, gr, p, save_spots); 02839 fit.run(); 02840 }
double spot_shape_s | ( | const TooN::Vector< 2 > & | x, | |
const TooN::Vector< 4, double, B > & | phi | |||
) | [inline] |
See spot_shape().
x | ![]() | |
phi | ![]() |
Definition at line 36 of file storm.h.
Referenced by spot_shape(), spot_shape_diff(), spot_shape_diff_position(), spot_shape_hess(), and spot_shape_hess_position().
std::pair<double, TooN::Vector<4> > spot_shape_diff_position | ( | const TooN::Vector< 2 > & | x, | |
const TooN::Vector< 4, double, B > & | phi | |||
) | [inline] |
Compute the spot shape and its derivative with respect to posision.
See also spot_shape()
x | ![]() | |
phi | ![]() |
Definition at line 46 of file storm.h.
References spot_shape_s(), and sq().
Referenced by SampledMultispot::compute_spot_intensity_derivatives().
00047 { 00048 using namespace TooN; 00049 00050 double s = spot_shape_s(x, phi); 00051 double r_2_pi = sqrt(2*M_PI); 00052 00053 double prob = exp(s) * phi[0]/(phi[1]*r_2_pi); 00054 00055 Vector<4> deriv = (exp(s) / (phi[1]*r_2_pi)) * 00056 makeVector(1, 00057 -phi[0] * (1 + 2*s)/phi[1], 00058 (x[0] - phi[2])*(phi[0]/sq(phi[1])), 00059 (x[1] - phi[3])*(phi[0]/sq(phi[1]))); 00060 return std::make_pair(prob, deriv); 00061 }
std::tr1::tuple<double, TooN::Vector<4>, TooN::Matrix<4> > spot_shape_hess_position | ( | const TooN::Vector< 2 > & | x, | |
const TooN::Vector< 4, double, B > & | phi | |||
) | [inline] |
Compute the spot shape and its Hessian with respect to posision.
See also spot_shape()
x | ![]() | |
phi | ![]() |
Definition at line 68 of file storm.h.
References spot_shape_s(), and sq().
Referenced by SampledMultispot::compute_spot_intensity_hessian().
00069 { 00070 using namespace TooN; 00071 using namespace std::tr1; 00072 00073 double s = spot_shape_s(x, phi); 00074 double r_2_pi = sqrt(2*M_PI); 00075 00076 double es = exp(s); 00077 00078 double prob = es * phi[0]/(phi[1]*r_2_pi); 00079 00080 Vector<4> deriv = (es / (phi[1]*r_2_pi)) * 00081 makeVector(1, 00082 -phi[0] * (1 + 2*s)/phi[1], 00083 (x[0] - phi[2])*(phi[0]/sq(phi[1])), 00084 (x[1] - phi[3])*(phi[0]/sq(phi[1]))); 00085 00086 Matrix<4> hess; 00087 hess[0][0] = 0; 00088 00089 hess[0][1] = -es*(1+2*s) / (phi[1] * phi[1] * r_2_pi); 00090 hess[1][0] = hess[0][1]; 00091 00092 hess[0][2] = es * (x[0] - phi[2]) / (pow(phi[1], 3)*r_2_pi); 00093 hess[2][0] = es * (x[0] - phi[2]) / (pow(phi[1], 3)*r_2_pi); 00094 00095 hess[0][3] = es * (x[1] - phi[3]) / (pow(phi[1], 3)*r_2_pi); 00096 hess[3][0] = es * (x[1] - phi[3]) / (pow(phi[1], 3)*r_2_pi); 00097 00098 hess[1][1] = 2*phi[0]*es*(1 + 5*s + 2*s*s) / ( pow(phi[1], 3) * r_2_pi); 00099 00100 hess[1][2] = -phi[0] * es * (3 + 2*s) * (x[0] - phi[2]) / (pow(phi[1], 4) * r_2_pi); 00101 hess[1][3] = -phi[0] * es * (3 + 2*s) * (x[1] - phi[3]) / (pow(phi[1], 4) * r_2_pi); 00102 00103 hess[2][1] = hess[1][2]; 00104 hess[3][1] = hess[1][3]; 00105 00106 hess[2][2] = phi[0] * es * (sq(x[0] - phi[2]) - sq(phi[1])) / (r_2_pi * pow(phi[1], 5)); 00107 hess[3][3] = phi[0] * es * (sq(x[1] - phi[3]) - sq(phi[1])) / (r_2_pi * pow(phi[1], 5)); 00108 00109 hess[2][3] = phi[0] * es * (x[0] - phi[2])*(x[1] - phi[3]) / (r_2_pi * pow(phi[1], 5)); 00110 hess[3][2] = hess[2][3]; 00111 00112 00113 return make_tuple(prob, deriv, hess); 00114 }
std::tr1::tuple<double, TooN::Vector<2>, TooN::Matrix<2> > spot_shape_hess | ( | const TooN::Vector< 2 > & | x, | |
const TooN::Vector< 4, double, B > & | phi | |||
) | [inline] |
Value of the spot, given the parameters and input location.
The spot is described by the following formula:
where
This describes a generic blobby spot function of a variable size. The light output can be tuned by varying , and the level of blur can be changed independently by varying
. The derivative is:
And the hessian is:
x | ![]() | |
phi | ![]() |
Definition at line 143 of file storm.h.
References spot_shape_s().
Referenced by log_probability_spot_hess().
00144 { 00145 double s = spot_shape_s(x, phi); 00146 double r_2_pi = sqrt(2*M_PI); 00147 00148 double prob = exp(s) * phi[0]/(phi[1]*r_2_pi); 00149 TooN::Vector<2> deriv = (exp(s) / (phi[1]*r_2_pi)) * TooN::makeVector(1, -phi[0] * (1 + 2*s)/phi[1]); 00150 TooN::Matrix<2> hess; 00151 00152 hess[0][0] = 0; 00153 hess[0][1] = -exp(s)*(1+2*s) / (phi[1] * phi[1] * r_2_pi); 00154 hess[1][0] = hess[0][1]; 00155 hess[1][1] = 2*phi[0]*exp(s)*(1 + 5*s + 2*s*s) / ( pow(phi[1], 3) * r_2_pi); 00156 00157 return std::tr1::make_tuple(prob, deriv, hess); 00158 }
std::pair<double, TooN::Vector<2> > spot_shape_diff | ( | const TooN::Vector< 2 > & | x, | |
const TooN::Vector< 4, double, B > & | phi | |||
) | [inline] |
x | ![]() | |
phi | ![]() |
Definition at line 165 of file storm.h.
References spot_shape_s().
Referenced by log_probability_spot_diff().
00166 { 00167 double s = spot_shape_s(x, phi); 00168 double r_2_pi = sqrt(2*M_PI); 00169 00170 double prob = exp(s) * phi[0]/(phi[1]*r_2_pi); 00171 TooN::Vector<2> deriv = (exp(s) / (phi[1]*r_2_pi)) * TooN::makeVector(1, -phi[0] * (1 + 2*s)/phi[1]); 00172 return std::make_pair(prob, deriv); 00173 }
double spot_shape | ( | const TooN::Vector< 2 > & | x, | |
const TooN::Vector< 4, double, B > & | phi | |||
) | [inline] |
x | ![]() | |
phi | ![]() |
Definition at line 181 of file storm.h.
References spot_shape_s().
Referenced by SampledMultispot::SpotWithBackgroundMasked::hess_log(), and log_probability_spot().
00182 { 00183 double s = spot_shape_s(x, phi); 00184 double r_2_pi = sqrt(2*M_PI); 00185 00186 // FIXME FIXME FIXME and don't forget to fix the HESSIAN AND DERIVATIVE 00187 // Should be: 1/(2 pi s^2) for two dimensions 00188 // vvvvvvvvvvvvv http://lol.i.trollyou.com/ 00189 double prob = exp(s) * phi[0]/(phi[1]*r_2_pi); 00190 00191 00192 return prob; 00193 }
std::tr1::tuple<double, TooN::Vector<2>, TooN::Matrix<2> > log_probability_spot_hess | ( | const CVD::SubImage< float > & | im, | |
double | variance, | |||
const TooN::Vector< 4, double, Base > & | spot_parameters | |||
) | [inline] |
Find the log probability of an image patch, assuming zero base-line mean and the given variance.
This function makes use of the spot shape. It is assumed that the centre pixel of the image is at 0,0. Since the noise is Gaussian:
where I is the image, and N is the number of pixels. See also log_probability_no_spot and (spot_shape). The derivatives are:
im | Image | |
variance | ![]() | |
spot_parameters | ![]() |
Definition at line 236 of file storm.h.
References spot_shape_hess(), and sq().
00237 { 00238 using namespace TooN; 00239 using namespace std::tr1; 00240 00241 //-1 because if the image is 3x3, ie 0,1,2 then 1,1 is the centre. 00242 //If it is 2x2, ie 0,1 then .5,.5 is the centre 00243 Vector<2> centre = makeVector((im.size().x-1) / 2.0, (im.size().y-1) / 2.0); 00244 00245 double logprob_part=0; 00246 Vector<2> diff = Zeros; 00247 Matrix<2> hess = Zeros; 00248 for(int y=0; y < im.size().y; y++) 00249 for(int x=0; x < im.size().x; x++) 00250 { 00251 Vector<2> d = TooN::makeVector(x, y) - centre; 00252 00253 double mu; 00254 Vector<2> diff_mu; 00255 Matrix<2> hess_mu; 00256 tie(mu, diff_mu, hess_mu) = spot_shape_hess(d, spot_parameters); 00257 00258 double e = im[y][x] - mu; 00259 00260 logprob_part += -sq(e); 00261 diff += diff_mu * e; 00262 hess += e * hess_mu - diff_mu.as_col() * diff_mu.as_row(); 00263 } 00264 return make_tuple( logprob_part / (2*variance) - im.size().area() * log(2*M_PI*variance)/2, 00265 diff / variance, 00266 hess / variance); 00267 }
std::pair<double, TooN::Vector<2> > log_probability_spot_diff | ( | const CVD::SubImage< float > & | im, | |
double | variance, | |||
const TooN::Vector< 4, double, Base > & | spot_parameters | |||
) | [inline] |
See log_probability_spot_hess.
im | Image | |
variance | ![]() | |
spot_parameters | ![]() |
Definition at line 276 of file storm.h.
References spot_shape_diff(), and sq().
00277 { 00278 using namespace TooN; 00279 using namespace std::tr1; 00280 using namespace std; 00281 //-1 because if the image is 3x3, ie 0,1,2 then 1,1 is the centre. 00282 //If it is 2x2, ie 0,1 then .5,.5 is the centre 00283 Vector<2> centre = makeVector((im.size().x-1) / 2.0, (im.size().y-1) / 2.0); 00284 00285 double logprob_part=0; 00286 Vector<2> diff = Zeros; 00287 for(int y=0; y < im.size().y; y++) 00288 for(int x=0; x < im.size().x; x++) 00289 { 00290 Vector<2> d = makeVector(x, y) - centre; 00291 00292 double mu; 00293 Vector<2> diff_mu; 00294 tie(mu, diff_mu) = spot_shape_diff(d, spot_parameters); 00295 00296 double e = im[y][x] - mu; 00297 00298 logprob_part += -sq(e); 00299 diff += diff_mu * e; 00300 } 00301 return make_pair(logprob_part / (2*variance) - im.size().area() * log(2*M_PI*variance)/2, diff / variance); 00302 }
double log_probability_spot | ( | const CVD::SubImage< float > & | im, | |
double | variance, | |||
const TooN::Vector< 4, double, Base > & | spot_parameters | |||
) | [inline] |
See log_probability_spot_hess.
im | Image | |
variance | ![]() | |
spot_parameters | ![]() |
Definition at line 311 of file storm.h.
References spot_shape(), and sq().
00312 { 00313 //-1 because if the image is 3x3, ie 0,1,2 then 1,1 is the centre. 00314 //If it is 2x2, ie 0,1 then .5,.5 is the centre 00315 TooN::Vector<2> centre = TooN::makeVector((im.size().x-1) / 2.0, (im.size().y-1) / 2.0); 00316 00317 double logprob_part=0; 00318 for(int y=0; y < im.size().y; y++) 00319 for(int x=0; x < im.size().x; x++) 00320 { 00321 TooN::Vector<2> d = TooN::makeVector(x, y) - centre; 00322 00323 double mu = spot_shape(d, spot_parameters); 00324 00325 double e = im[y][x] - mu; 00326 00327 logprob_part += -sq(e); 00328 } 00329 return logprob_part / (2*variance) - im.size().area() * log(2*M_PI*variance)/2; 00330 }
double log_normal_std | ( | double | mu, | |
double | sigma | |||
) | [inline] |
Compute the standard deviation of a log-normal distribution.
See log_normal().
sigma | ![]() | |
mu | ![]() |
Definition at line 343 of file storm.h.
References sq().
Referenced by brightness_motion_limit().
double log_normal_mode | ( | double | mu, | |
double | sigma | |||
) | [inline] |
Compute the mode of a log-normal distribution.
See log_normal().
sigma | ![]() | |
mu | ![]() |
Definition at line 359 of file storm.h.
Referenced by fit_spots_historic(), generate_state_parameters_ye_olde(), and FitSpots::try_modifying_model().
double log_log_normal | ( | double | x, | |
double | mu, | |||
double | sigma | |||
) | [inline] |
Log-normal distribution.
This is given by:
x | x | |
mu | ![]() | |
sigma | ![]() |
Definition at line 374 of file storm.h.
Referenced by NegativeFreeEnergy::compute_with_mask(), and NegativeFreeEnergy::operator()().
double diff_log_log_normal | ( | double | x, | |
double | mu, | |||
double | sigma | |||
) | [inline] |
Derivative of the log of the log-normal distribution:
.
x | x | |
mu | ![]() | |
sigma | ![]() |
Definition at line 388 of file storm.h.
Referenced by FreeEnergyHessian::hessian(), SpotNegProbabilityDiffWithSampledBackground::operator()(), sampled_background_spot_hessian2(), and sampled_background_spot_hessian_ffbs().
double hess_log_log_normal | ( | double | x, | |
double | mu, | |||
double | sigma | |||
) | [inline] |
Second derivative of the log of the log-normal distribution:
.
x | x | |
mu | ![]() | |
sigma | ![]() |
Definition at line 403 of file storm.h.
Referenced by FreeEnergyHessian::hessian(), sampled_background_spot_hessian2(), sampled_background_spot_hessian_FAKE(), and sampled_background_spot_hessian_ffbs().