From a52f31c95b65f4e97b4ee56ea1b1e5ea51991860 Mon Sep 17 00:00:00 2001 From: Zhuoran Han Date: Fri, 10 Aug 2018 16:40:17 -0500 Subject: [PATCH] Add cpp math functions to calc normalization arrays Add math funcs to reconstruct FET --- src/math_functions.cpp | 52 +++++++++++++++++++++++++++++++++++++++--- src/math_functions.h | 46 ++++++++++++++++++++++++++++++++++++- 2 files changed, 94 insertions(+), 4 deletions(-) diff --git a/src/math_functions.cpp b/src/math_functions.cpp index 8a36f7b922..b5c718d5d0 100644 --- a/src/math_functions.cpp +++ b/src/math_functions.cpp @@ -103,17 +103,36 @@ void calc_pn_c(int n, double x, double pnx[]) { } -double evaluate_legendre_c(int n, const double data[], double x) { +double evaluate_legendre_c(int n, const double data[], double x, int flag) { double pnx[n + 1]; double val = 0.0; calc_pn_c(n, x, pnx); - for (int l = 0; l <= n; l++) { - val += (l + 0.5) * data[l] * pnx[l]; + + if (flag == 1) { + for (int l = 0; l <= n; l++) { + val += (l + 0.5) * data[l] * pnx[l]; + } + } else if (flag == 2) { + for (int l = 0; l <= n; l++) { + val += data[l] * pnx[l]; + } } + return val; } +void calc_norm_pn_c(int l, double zmin, double zmax, double norm_vec[]) { + // =========================================================================== + // Set up the normalization vector for FET Pn reconstruction + double d = zmax - zmin; + + for (int i = 0; i <= l; i++) { + norm_vec[i] = (2 * i + 1)/d; + } +} + + void calc_rn_c(int n, const double uvw[3], double rn[]){ // rn[] is assumed to have already been allocated to the correct size @@ -615,6 +634,33 @@ void calc_zn_rad_c(int n, double rho, double zn_rad[]) { } +void calc_norm_zn_rad_c(int n, double r, double norm_vec[]) { + // =========================================================================== + // Set up the normalization vector for FET Zn Radial only reconstruction + int index; + + for (int i = 0; i <= n; i+=2) { + index = int(i/2); + norm_vec[index] = (i + 1)/(PI * r * r); + } +} + + +double evaluate_zernike_rad_c(int n, const double data[], double r) { + double rnr[n / 2 + 1]; + double val = 0.0; + calc_zn_rad_c(n, r, rnr); + int index; + + for (int i = 0; i <= n; i+ = 2) { + index = int(i/2); + val += data[index] * rnr[index]; + } + + return val; +} + + void rotate_angle_c(double uvw[3], double mu, double* phi) { // Copy original directional cosines double u0 = uvw[0]; // original cosine in x direction diff --git a/src/math_functions.h b/src/math_functions.h index 7ef39c4461..30e210e500 100644 --- a/src/math_functions.h +++ b/src/math_functions.h @@ -53,11 +53,27 @@ extern "C" void calc_pn_c(int n, double x, double pnx[]); //! @param data The polynomial expansion coefficient data; without the (2l+1)/2 //! factor. //! @param x The value to evaluate at; x is expected to be within [-1,1] +//! @param flag A flag of which method to calculate Legendre polynomials. 1 for +//! default (2l+1)/2; 2 for custom normalized coefficients //! @return The requested Legendre polynomials of order 0 to n (inclusive) //! evaluated at x //============================================================================== -extern "C" double evaluate_legendre_c(int n, const double data[], double x); +extern "C" double evaluate_legendre_c(int n, const double data[], double x, int flag = 1); + +//============================================================================== +//! Evaluate the normalization vector when reconstructing functional expansion +//! from Legendre polynomials. +//! +//! The normalization is (2 * l + 1)/(zmax - zmin) +//! +//! @param l The maximum order requested +//! @param zmin The minimum location of 1-d domain +//! @param zmax The maximum location of 1-d domain +//! @param norm_vec The requested normalization of order 0 to n (inclusive) +//============================================================================== + +extern "C" void calc_norm_pn_c(int l, double zmin, double zmax, double norm_vec[]); //============================================================================== //! Calculate the n-th order real spherical harmonics for a given angle (in @@ -69,6 +85,7 @@ extern "C" double evaluate_legendre_c(int n, const double data[], double x); //! evaluated at uvw. //============================================================================== + extern "C" void calc_rn_c(int n, const double uvw[3], double rn[]); //============================================================================== @@ -109,6 +126,33 @@ extern "C" void calc_zn_c(int n, double rho, double phi, double zn[]); extern "C" void calc_zn_rad_c(int n, double rho, double zn_rad[]); +//============================================================================== +//! Evaluate the normalization vector when reconstructing functional expansion +//! from radial only Zernike polynomials. +//! +//! The normalization is (2 * l + 1)/(zmax - zmin) +//! +//! @param n The maximum order requested +//! @param r The radius of the disk +//! @param norm_vec The requested normalization of even orders from 0 to n +//============================================================================== + +extern "C" void calc_norm_zn_rad_c(int n, double r, double norm_vec[]); + +//============================================================================== +//! Find the value of f(x) given a set of even order Zernike coefficients and +//! the value of x if there is only radial dependency. +//! +//! @param n The maximum order of the expansion +//! @param data The overall polynomial expansion coefficient data; +//! normalization should be considered by users +//! @param r The value to evaluate at; r is expected to be within [0,1] +//! @return The requested Zernike polynomials of order 0 to n (inclusive) +//! evaluated at r +//============================================================================== + +extern "C" double evaluate_zernike_rad_c(int n, const double data[], double r); + //============================================================================== //! Rotate the direction cosines through a polar angle whose cosine is mu and //! through an azimuthal angle sampled uniformly.