models/cdm: applied fixes from AlTarFramework#11 submitted by @gracebato

altar/bayesian/Metropolis.py

@@ -51,7 +51,7 @@ def initialize(self, application):
         rng = application.rng.rng
         # set up the distribution for building the sample multiplicities; use a strictly
         # positive distribution to avoid generating candidates with zero displacement
-        self.uniform = altar.pdf.uniform_pos(support=(0,1), rng=rng)
+        self.uniform = altar.pdf.uniform_pos(rng=rng)
         # set up the distribution for the random walk displacement vectors
         self.uninormal = altar.pdf.ugaussian(rng=rng)
 

models/cdm/cdm/Fast.py

@@ -89,8 +89,10 @@ def dataLikelihood(self, model, step):
         for sample in range(samples):
             # get the residuals
             residuals = predicted.getRow(sample)
-            # compute the norm, and normalize it
-            llk = normalization - norm.eval(v=residuals, sigma_inv=cd_inv) / 2
+            # compute the norm
+            nrm = norm.eval(v=residuals, sigma_inv=cd_inv)
+            # and normalize it
+            llk = normalization - nrm**2 / 2
             # store it
             dataLLK[sample] = llk
 

models/cdm/cdm/Native.py

@@ -55,13 +55,13 @@ def dataLikelihood(self, model, step):
         xIdx = model.xIdx
         yIdx = model.yIdx
         dIdx = model.dIdx
+        openingIdx = model.openingIdx
         aXIdx = model.aXIdx
         aYIdx = model.aYIdx
         aZIdx = model.aZIdx
         omegaXIdx = model.omegaXIdx
         omegaYIdx = model.omegaYIdx
         omegaZIdx = model.omegaZIdx
-        openingIdx = model.openingIdx
         offsetIdx = model.offsetIdx
 
         # get the observations
@@ -91,9 +91,9 @@ def dataLikelihood(self, model, step):
             omegaZ = parameters[omegaZIdx]
 
             # 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,
                          ax=aX, ay=aY, az=aZ, omegaX=omegaX, omegaY=omegaY, omegaZ=omegaZ,
-                         opening=opening, v=model.nu)
+                         v=model.nu)
             # compute the expected displacement
             u = cdm.displacements(locations=locations, los=los)
 
@@ -105,9 +105,9 @@ def dataLikelihood(self, model, step):
                 u[obs] -= parameters[offsetIdx + oid[obs]]
 
             # compute the norm of the displacements
-            residual = norm.eval(v=u, sigma_inv=cd_inv)
+            nrm = 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 - nrm**2 /2
 
         # all done
         return self

models/cdm/cdm/Source.py

@@ -78,9 +78,10 @@ def displacements(self, locations, los):
         # allocate space for the result
         u = altar.vector(shape=len(locations))
         # compute the displacements
-        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(X=Xf, Y=Yf, X0=x_src, Y0=y_src, depth=d_src,
+                          ax=ax_src, ay=ay_src, az=az_src,
+                          omegaX=omegaX_src, omegaY=omegaY_src, omegaZ=omegaZ_src,
+                          opening=opening, nu=v)
         # go through each observation location
         for idx, (ux,uy,uz) in enumerate(zip(ue, un, uv)):
             # project the expected displacement along LOS and store
@@ -91,8 +92,9 @@ def displacements(self, locations, los):
 
 
     # 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,
+                 ax=ax, ay=ay, az=az, omegaX=omegaX, omegaY=omegaY, omegaZ=omegaZ,
+                 opening=opening, v=v, **kwds):
         # chain up
         super().__init__(**kwds)
         # store the location

models/cdm/cdm/libcdm.py

@@ -35,7 +35,7 @@ def norm(v):
     return numpy.sqrt(v.dot(v))
 
 
-def CDM(X, Y, X0, Y0, depth, omegaX, omegaY, omegaZ, ax, ay, az, opening, nu):
+def CDM(X, Y, X0, Y0, depth, ax, ay, az, omegaX, omegaY, omegaZ, opening, nu):
     """
     CDM
     calculates the surface displacements and potency associated with a CDM

models/cdm/lib/libcdm/cdm.cc

@@ -364,9 +364,9 @@ inline static
 altar::models::cdm::mat_t
 altar::models::cdm::
 operator*(const mat_t & m1, const mat_t & m2) {
-    return { m1[0]*m2[0] + m1[2]*m2[3] + m1[2]*m2[6],
-             m1[0]*m2[1] + m1[2]*m2[4] + m1[2]*m2[7],
-             m1[0]*m2[2] + m1[2]*m2[5] + m1[2]*m2[8],
+    return { m1[0]*m2[0] + m1[1]*m2[3] + m1[2]*m2[6],
+             m1[0]*m2[1] + m1[1]*m2[4] + m1[2]*m2[7],
+             m1[0]*m2[2] + m1[1]*m2[5] + m1[2]*m2[8],
 
              m1[3]*m2[0] + m1[4]*m2[3] + m1[5]*m2[6],
              m1[3]*m2[1] + m1[4]*m2[4] + m1[5]*m2[7],
@@ -411,24 +411,24 @@ inline static
 altar::models::cdm::vec_t
 altar::models::cdm::
 xform(const mat_t & m, const vec_t & v) {
-    return {m[0]*v[0] + m[1]*v[2] + m[2]*v[3],
-            m[3]*v[0] + m[4]*v[2] + m[5]*v[3],
-            m[6]*v[0] + m[7]*v[2] + m[8]*v[3]};
+    return {m[0]*v[0] + m[1]*v[1] + m[2]*v[2],
+            m[3]*v[0] + m[4]*v[1] + m[5]*v[2],
+            m[6]*v[0] + m[7]*v[1] + m[8]*v[2]};
 };
 
 // trig
 inline
 double
 altar::models::cdm::
 sin(double omega) {
-    return std::sin(omega * 180/pi);
+    return std::sin(omega * pi/180);
 };
 
 inline
 double
 altar::models::cdm::
 cos(double omega) {
-    return std::cos(omega * 180/pi);
+    return std::cos(omega * pi/180);
 };
 
 // end-of-file