Skip to content

Conversation

@gracebato
Copy link

  1. proper calculation of LLK (e.g. LLK = normalization - normeval**2.0 / 2.0 )
  2. avoiding 0 values in the random generator within Metropolis.py
  3. dV is now treated as a non-Jeffrey parameter to allow negative source (deflation).
  4. LOS matrix for the EW component should have a negative sign.

… 2) avoiding 0 values in the random generator within Metropolis.py, 3) dV is treated as a non-Jeffrey parameter.
@lijun99
Copy link
Contributor

lijun99 commented Mar 24, 2019

The framework changes were also included in my previous pull request.

  1. proper calculation of LLK (e.g. LLK = normalization - normeval**2.0 / 2.0 )
    I have proposed to use ddot (v*v) instead of dnrm2 (|v|) in L2.py:eval. Maybe a compromise is that we define a new function evalsquare for (v*v) so that users can pick either eval (LLK = -0.5*norm*norm) or evalsquare (LLK = -0.5*normsq) in their models. I will test which one is more efficient, but judging from openblas x86_64 kernels, ddot is better optimized with streaming SIMD extensions.

  2. avoiding 0 values in the random generator within Metropolis.py
    gsl offers a positive uniform distribution implementation, which is now included in pyre. Use
    self.uniform = altar.pdf.uniform_pos(rng=rng) instead. Adding a range requires rescaling for each generated random number, which might not be efficient.

@aivazis aivazis self-assigned this Apr 11, 2019
libcdm.layout(source,
model.xIdx, model.dIdx,
model.openingIdx, model.aXIdx, model.omegaXIdx,
model.xIdx, model.dIdx, model.openingIdx, model.aXIdx, model.omegaXIdx,
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Modified function variable referencing/order for consistency.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

this looks like merging two lines into one; am i missing something

# compute the norm, and normalize it
llk = normalization - norm.eval(v=residuals, sigma_inv=cd_inv) / 2
normeval = norm.eval(v=residuals, sigma_inv=cd_inv)
llk = normalization - normeval**2.0 / 2.0
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Corrected the calculation of LLK.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

good catch!

xIdx = model.xIdx
yIdx = model.yIdx
dIdx = model.dIdx
openingIdx = model.openingIdx
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Nothing significant. Just reordered the variable for consistency.


# make a source using the sample parameters
cdm = source(x=x, y=y, d=d,
cdm = source(x=x, y=y, d=d, opening=opening,
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Modified function variable referencing/order for consistency.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

with python pass-by-name, the order of arguments is not important

normeval = norm.eval(v=u, sigma_inv=cd_inv)
# normalize and store it as the data log likelihood
dataLLK[sample] = normalization - residual/2
dataLLK[sample] = normalization - normeval**2.0 /2
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Corrected the calculation of the LLK.

ue, un, uv = CDM(Xf, Yf, x_src, y_src, d_src,
omegaX_src, omegaY_src, omegaZ_src, ax_src, ay_src, az_src,
opening, v)
ue, un, uv = CDM(Xf, Yf, x_src, y_src, d_src, opening,
Copy link
Author

@gracebato gracebato Apr 25, 2019

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Lines 81-83: Modified function variable referencing/order for consistency.

# meta-methods
def __init__(self, x=x, y=y, d=d, omegaX=omegaX, omegaY=omegaY, omegaZ=omegaZ,
ax=ax, ay=ay, az=az, opening=opening, v=v, **kwds):
def __init__(self, x=x, y=y, d=d, opening=opening, ax=ax, ay=ay, az=az,
Copy link
Author

@gracebato gracebato Apr 25, 2019

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Lines 95-96: Same here--modified function variable referencing/order for consistency.



def CDM(X, Y, X0, Y0, depth, omegaX, omegaY, omegaZ, ax, ay, az, opening, nu):
def CDM(X, Y, X0, Y0, depth, opening, ax, ay, az, omegaX, omegaY, omegaZ, nu):
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Also here--modified function variable referencing/order for consistency.

import itertools
# build a set of points on a grid
stations = itertools.product(range(-5,6), range(-5,6))
import numpy as np
Copy link
Author

@gracebato gracebato Apr 25, 2019

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Line 171-177: I think for the synthetic case, this is better to constrain the number of points.

# make a source
source = altar.models.cdm.source(
x=self.x, y=self.y, d=self.d,
x=self.x, y=self.y, d=self.d, opening=self.opening,
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Lines 99-102: Modified function variable referencing/order for consistency.

for idx, observation in enumerate(data):
# set the covariance to a fraction of the "observed" displacement
correlation[idx,idx] = 1.0 #.01 * observation.d
correlation[idx,idx] = 0.0001 #1.0 #.01 * observation.d
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

1.0 m is quite large for an observation error (variance).

// compute the displacements
cdm(_locations,
xSrc, ySrc, dSrc,
openingSrc,
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Modified function variable referencing/order for consistency.

auto uy = gsl_matrix_get(disp, loc, 1);
auto ud = gsl_matrix_get(disp, loc, 2);

// project; don't forget {ud} is positive into the ground
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is not a depth value. It is the vertical component.

// project; don't forget {ud} is positive into the ground
auto u = ux*nx + uy*ny - ud*nz;
// project
auto u = ux*nx + uy*ny + ud*nz;
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is not a depth value and should be '+' sign.

double opening,
double aX, double aY, double aZ,
double omegaX, double omegaY, double omegaZ,
double opening,
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Modified function variable referencing/order for consistency.

namespace cdm {
void cdm(const gsl_matrix * locations,
double X0, double Y0, double depth,
double opening,
Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Modified function variable referencing/order for consistency.

@gracebato gracebato changed the title Changes made to make mogi run properly Changes made to make mogi and cdm run properly Apr 25, 2019
gracebato pushed a commit to gracebato/altar that referenced this pull request Jul 15, 2019
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants