/** @file  psConstants.h
 *
 *  This file will hold definitions of various constants as well as common
 *  macros used throughout psLib.
 *
 *  @author GLG, MHPCC
 *
 *  @version $Revision: 1.73 $ $Name: not supported by cvs2svn $
 *  @date $Date: 2005-06-08 23:40:45 $
 *
 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
 *
 *  XXX: Add parenthesis around all arguments so that these macros can be
 *       called with complex expressions.
 *
 *  XXX: All functions which use the PS_CHECK macros must be scrutinized so
 *  that we ensure that an argument which is expected to be output is
 *  psFree'ed before reurning NULL.
 *
 *  XXX: The macros have a name similar to PS_CHECK_CONDITION() and generally
 *  throw a psError if the CONDITION is true.  However, some throw the error
 *  if the CONDITION is false.  This should be consistant.
 *
 *  XXX: rename all these with form PS_ASSERT_CONDITION().
 *
 */

#include<math.h> // for M_PI

/*****************************************************************************
These constants are used by various functions in the psLib.
 *****************************************************************************/
#define PS_DETERMINE_BRACKET_STEP_SIZE 0.10
#define PS_MAX_LMM_ITERATIONS 100
#define PS_MAX_MINIMIZE_ITERATIONS 100
#define PS_LEFT_SPLINE_DERIV 0.0
#define PS_RIGHT_SPLINE_DERIV 0.0
/*****************************************************************************
These are common mathimatical constants used by various functions in the psLib.
 *****************************************************************************/
#ifndef M_PI
#define M_PI   3.1415926535897932384626433832795029  /* pi */
#define M_PI_2 1.5707963267948966192313216916397514  /* pi/2 */
#define M_PI_4 0.7853981633974483096156608458198757  /* pi/4 */
#define M_1_PI 0.3183098861837906715377675267450287  /* 1/pi */
#define M_2_PI 0.6366197723675813430755350534900574  /* 2/pi */
#endif // #ifndef M_PI

#define DEG_TO_RAD(DEGREES) ((DEGREES) * M_PI / 180.0)
#define MIN_TO_RAD(MINUTES) ((MINUTES) * M_PI / (180.0 * 60.0))
#define SEC_TO_RAD(SECONDS) ((SECONDS) * M_PI / (180.0 * 60.0 * 60.0))
#define RAD_TO_DEG(RADIANS) ((RADIANS) * 180.0 / M_PI)
#define RAD_TO_MIN(RADIANS) ((RADIANS) * 180.0 * 60.0 / M_PI)
#define RAD_TO_SEC(RADIANS) ((RADIANS) * 180.0 * 60.0 * 60.0 / M_PI)

/*****************************************************************************
 
*****************************************************************************/

#define PS_ASSERT_INT_UNEQUAL(NAME1, NAME2, RVAL) \
if ((NAME1) == (NAME2)) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: %s and %s are equal.", \
            #NAME1, #NAME2); \
    return(RVAL); \
}

#define PS_ASSERT_INT_EQUAL(NAME1, NAME2, RVAL) \
if ((NAME1) != (NAME2)) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: %s and %s are not equal.", \
            #NAME1, #NAME2); \
    return(RVAL); \
}

#define PS_ASSERT_INT_NONNEGATIVE(NAME, RVAL) \
if ((NAME) < 0) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: %s is less than 0.", #NAME); \
    return(RVAL); \
}

#define PS_ASSERT_INT_POSITIVE(NAME, RVAL) \
if ((NAME) < 1) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: %s is 0 or less.", #NAME); \
    return(RVAL); \
}

#define PS_ASSERT_INT_ZERO(NAME, RVAL) \
if ((NAME) != 0) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: %s is 0.", #NAME); \
    return(RVAL); \
}

#define PS_ASSERT_INT_NONZERO(NAME, RVAL) \
if ((NAME) == 0) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: %s is 0.", #NAME); \
    return(RVAL); \
}

// XXX: Where did these int casts come from?
#define PS_ASSERT_INT_WITHIN_RANGE(NAME, LOWER, UPPER, RVAL) \
if ((int)(NAME) < LOWER || (int)(NAME) > UPPER) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: %s, %d, is out of range.  Must be between %d and %d.", \
            #NAME,(int)NAME,LOWER,UPPER); \
    return RVAL; \
}

#define PS_ASSERT_INT_LARGER_THAN(NAME1, NAME2, RVAL) \
if (!((NAME1) > (NAME2))) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: !(%s > %s) (%d %d).", \
            #NAME1, #NAME2, NAME1, NAME2); \
    return(RVAL); \
}

#define PS_ASSERT_INT_LARGER_THAN_OR_EQUAL(NAME1, NAME2, RVAL) \
if (!((NAME1) >= (NAME2))) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: !(%s >= %s) (%d %d).", \
            #NAME1, #NAME2, NAME1, NAME2); \
    return(RVAL); \
}
#define PS_ASSERT_FLOAT_LARGER_THAN(NAME1, NAME2, RVAL) \
if (!((NAME1) > (NAME2))) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: !(%s > %s) (%f %f).", \
            #NAME1, #NAME2, NAME1, NAME2); \
    return(RVAL); \
}

#define PS_ASSERT_FLOAT_LARGER_THAN_OR_EQUAL(NAME1, NAME2, RVAL) \
if (!((NAME1) >= (NAME2))) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: !(%s >= %s) (%f %f).", \
            #NAME1, #NAME2, NAME1, NAME2); \
    return(RVAL); \
}

// Produce an error if (NAME1 > NAME2)
// XXX: Get rid of this, use above macros.
#define PS_INT_COMPARE(NAME1, NAME2, RVAL) \
if ((NAME1) > (NAME2)) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: (%s > %s) (%d %d).", \
            #NAME1, #NAME2, NAME1, NAME2); \
    return(RVAL); \
}

// Produce an error if ((NAME1 > NAME2)
// XXX: Get rid of this, use above macros.
#define PS_FLOAT_COMPARE(NAME1, NAME2, RVAL) \
if ((NAME1) > (NAME2)) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: (%s > %s) (%f %f)", \
            #NAME1, #NAME2, NAME1, NAME2); \
    return(RVAL); \
}

#define PS_ASSERT_FLOAT_NON_EQUAL(NAME1, NAME2, RVAL) \
if (fabs((NAME2) - (NAME1)) < FLT_EPSILON) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: %s and %s are equal.", \
            #NAME1, #NAME2); \
    return(RVAL); \
}

#define PS_ASSERT_FLOAT_EQUAL(NAME1, NAME2, RVAL) \
if (fabs((NAME2) - (NAME1)) > FLT_EPSILON) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: %s and %s are not equal.", \
            #NAME1, #NAME2); \
    return(RVAL); \
}

// Return an error if the arg is lies outside the supplied range.
#define PS_ASSERT_FLOAT_WITHIN_RANGE(NAME, LOWER, UPPER, RVAL) \
if ((NAME) < (LOWER) || (NAME) > (UPPER)) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: %s, %f, is out of range.  Must be between %f and %f.", \
            #NAME, NAME, LOWER, UPPER); \
    return RVAL; \
}

// Return an error if the arg lies outside the supplied range
#define PS_ASSERT_LONG_WITHIN_RANGE(NAME, LOWER, UPPER, RVAL) \
if ((NAME) < (LOWER) || (NAME) > (UPPER)) { \
    psError(PS_ERR_BAD_PARAMETER_VALUE, true, \
            "Error: %s, %lld, is out of range.", \
            #NAME, NAME, LOWER, UPPER); \
    return RVAL; \
}

/*****************************************************************************
Macros which take a generic psLib type and determine if it is NULL, or has
the wrong type.
*****************************************************************************/
#define PS_ASSERT_PTR_NON_NULL(NAME, RVAL) PS_ASSERT_GENERAL_PTR_NON_NULL(NAME, return RVAL)
#define PS_ASSERT_GENERAL_PTR_NON_NULL(NAME, CLEANUP) \
if ((NAME) == NULL) { \
    psError(PS_ERR_BAD_PARAMETER_NULL, true, \
            "Unallowable operation: %s is NULL.", \
            #NAME); \
    CLEANUP; \
}

#define PS_ASSERT_PTR_TYPE(NAME, TYPE, RVAL) \
if ((NAME)->type.type != TYPE) { \
    psError(PS_ERR_BAD_PARAMETER_TYPE, true, \
            "Unallowable operation: %s has incorrect type.", \
            #NAME); \
    return(RVAL); \
}

#define PS_ASSERT_PTR_DIMEN(NAME, DIMEN, RVAL) PS_ASSERT_GENERAL_PTR_DIMEN(NAME, DIMEN, return RVAL)
#define PS_ASSERT_GENERAL_PTR_DIMEN(NAME, DIMEN, CLEANUP) \
if ((NAME)->type.dimen != DIMEN) { \
    psError(PS_ERR_BAD_PARAMETER_TYPE, true, \
            "Unallowable operation: %s has incorrect dimensionality.", \
            #NAME); \
    CLEANUP; \
}

#define PS_ASSERT_PTR_DIMEN_GENERAL_NOT(NAME, DIMEN, CLEANUP) \
if ((NAME)->type.dimen == DIMEN) { \
    psError(PS_ERR_BAD_PARAMETER_TYPE, true, \
            "Unallowable operation: %s has incorrect dimensionality.", \
            #NAME); \
    CLEANUP; \
}


#define PS_ASSERT_PTRS_SIZE_EQUAL(PTR1, PTR2, RVAL) \
if (PTR1->n != PTR2->n) { \
    psError(PS_ERR_BAD_PARAMETER_SIZE, true, \
            "ptr %s has size %d, ptr %s has size %d.", \
            #PTR1, PTR1->n, #PTR2, PTR2->n); \
    return(RVAL); \
}

#define PS_ASSERT_PTR_TYPE_EQUAL(PTR1, PTR2, RVAL) PS_ASSERT_PTRS_TYPE_EQUAL_GENERAL(PTR1, PTR2, return RVAL)
#define PS_ASSERT_PTRS_TYPE_EQUAL_GENERAL(PTR1, PTR2, CLEANUP) \
if (PTR1->type.type != PTR2->type.type) { \
    psError(PS_ERR_BAD_PARAMETER_TYPE, true, \
            "ptr %s has type %d, ptr %s has type %d.", \
            #PTR1, PTR1->type.type, #PTR2, PTR2->type.type); \
    CLEANUP; \
}


/*****************************************************************************
    PS_VECTOR macros:
 *****************************************************************************/
#define PS_ASSERT_VECTOR_NON_NULL(NAME, RVAL) PS_ASSERT_GENERAL_VECTOR_NON_NULL(NAME, return RVAL)
#define PS_ASSERT_GENERAL_VECTOR_NON_NULL(NAME, CLEANUP) \
if ((NAME) == NULL || (NAME)->data.U8 == NULL) { \
    psError(PS_ERR_BAD_PARAMETER_NULL, true, \
            "Unallowable operation: psVector %s or its data is NULL.", \
            #NAME); \
    CLEANUP; \
} \

#define PS_ASSERT_VECTOR_NON_EMPTY(NAME, RVAL) PS_ASSERT_GENERAL_VECTOR_NON_EMPTY(NAME, return RVAL)
#define PS_ASSERT_GENERAL_VECTOR_NON_EMPTY(NAME, CLEANUP) \
if ((NAME)->n < 1) { \
    psError(PS_ERR_BAD_PARAMETER_SIZE, true, \
            "Unallowable operation: psVector %s has no elements.", \
            #NAME); \
    CLEANUP; \
} \

#define PS_ASSERT_VECTOR_TYPE_F32_OR_F64(NAME, RVAL) \
if (((NAME)->type.type != PS_TYPE_F32) && ((NAME)->type.type != PS_TYPE_F64)) { \
    psError(PS_ERR_BAD_PARAMETER_TYPE, true, \
            "psVector %s: bad type(%d)", \
            #NAME, NAME->type.type); \
    return(RVAL); \
} \

#define PS_ASSERT_VECTOR_TYPE_S16_S32_F32(NAME, RVAL) \
if (((NAME)->type.type != PS_TYPE_S16) && ((NAME)->type.type != PS_TYPE_S32) && ((NAME)->type.type != PS_TYPE_F32)) { \
    psError(PS_ERR_BAD_PARAMETER_TYPE, true, \
            "psVector %s: bad type(%d)", \
            #NAME, NAME->type.type); \
    return(RVAL); \
} \

#define PS_ASSERT_VECTOR_TYPE(NAME, TYPE, RVAL) \
if ((NAME)->type.type != TYPE) { \
    psError(PS_ERR_BAD_PARAMETER_TYPE, true, \
            "Unallowable operation: psVector %s has incorrect type.", \
            #NAME); \
    return(RVAL); \
}

#define PS_ASSERT_VECTORS_SIZE_EQUAL(VEC1, VEC2, RVAL) \
if (VEC1->n != VEC2->n) { \
    psError(PS_ERR_BAD_PARAMETER_SIZE, true, \
            "psVector %s has size %d, psVector %s has size %d.", \
            #VEC1, VEC1->n, #VEC2, VEC2->n); \
    return(RVAL); \
}

#define PS_ASSERT_VECTOR_SIZE(VEC, SIZE, RVAL) \
if (VEC->n != SIZE) { \
    psError(PS_ERR_BAD_PARAMETER_SIZE, true, \
            "psVector %s has size %d, should be %d." \
            #VEC, VEC->n, SIZE); \
    return(RVAL); \
}

#define PS_ASSERT_VECTOR_TYPE_EQUAL(VEC1, VEC2, RVAL) \
if (VEC1->type.type != VEC2->type.type) { \
    psError(PS_ERR_BAD_PARAMETER_SIZE, true, \
            "psVector %s has size %d, psVector %s has size %d.", \
            #VEC1, VEC1->type.type, #VEC2, VEC2->type.type); \
    return(RVAL); \
}

#define PS_VECTOR_F64_TO_F32(X64, X32) \
psVector *X32 = psVectorAlloc(X64->n, PS_TYPE_F32); \
for (int i=0;i<X64->n;i++) { \
    X32->data.F32[i] = (float) X64->data.F64[i]; \
} \

#define PS_VECTOR_F32_TO_F64(X32, X64) \
psVector *X64 = psVectorAlloc(X32->n, PS_TYPE_F64); \
for (int i=0;i<X32->n;i++) { \
    X64->data.F64[i] = (float) X32->data.F32[i]; \
} \

#define PS_VECTOR_PRINT_F32(NAME) \
for (int my_i=0;my_i<(NAME)->n;my_i++) { \
    printf("%s->data.F32[%d] is %f\n", #NAME, my_i, (NAME)->data.F32[my_i]); \
} \
printf("\n"); \


#define PS_VECTOR_CONVERT_F64_TO_F32_STATIC(OLD, NEW_PTR32, NEW_STATIC32) \
if (OLD->type.type == PS_TYPE_F32) { \
    NEW_PTR32 = (psVector *) OLD; \
} else if (OLD->type.type == PS_TYPE_F64) { \
    NEW_STATIC32 = psVectorRecycle(NEW_STATIC32, OLD->n, PS_TYPE_F32); \
    p_psMemSetPersistent(NEW_STATIC32, true); \
    p_psMemSetPersistent(NEW_STATIC32->data.U8, true); \
    for (i=0; i < OLD->n ; i++) { \
        NEW_STATIC32->data.F32[i] = (float) OLD->data.F64[i]; \
    } \
    NEW_PTR32 = NEW_STATIC32; \
} \

#define PS_VECTOR_CONVERT_F32_TO_F64_STATIC(OLD, NEW_PTR64, NEW_STATIC64) \
if (OLD->type.type == PS_TYPE_F64) { \
    NEW_PTR64 = (psVector *) OLD; \
} else if (OLD->type.type == PS_TYPE_F32) { \
    NEW_STATIC64 = psVectorRecycle(NEW_STATIC64, OLD->n, PS_TYPE_F64); \
    p_psMemSetPersistent(NEW_STATIC64, true); \
    p_psMemSetPersistent(NEW_STATIC64->data.U8, true); \
    for (i=0; i < OLD->n ; i++) { \
        NEW_STATIC64->data.F64[i] = (double) OLD->data.F32[i]; \
    } \
    NEW_PTR64 = NEW_STATIC64; \
} \

#define PS_VECTOR_GEN_YERR_STATIC_F32(VEC, N) \
VEC = psVectorRecycle(VEC, N, PS_TYPE_F32); \
p_psMemSetPersistent(VEC, true); \
p_psMemSetPersistent(VEC->data.U8, true); \
for (int i=0;i<N;i++) { \
    VEC->data.F32[i] = 1.0; \
} \

#define PS_VECTOR_GEN_YERR_STATIC_F64(VEC, N) \
VEC = psVectorRecycle(VEC, N, PS_TYPE_F64); \
p_psMemSetPersistent(VEC, true); \
p_psMemSetPersistent(VEC->data.U8, true); \
for (int i=0;i<N;i++) { \
    VEC->data.F64[i] = 1.0; \
} \

#define PS_VECTOR_GEN_X_INDEX_STATIC_F32(VEC, N) \
VEC = psVectorRecycle(VEC, N, PS_TYPE_F32); \
p_psMemSetPersistent(VEC, true); \
p_psMemSetPersistent(VEC->data.U8, true); \
for (int i=0;i<N;i++) { \
    VEC->data.F32[i] = (float) i; \
} \

#define PS_VECTOR_GEN_X_INDEX_STATIC_F64(VEC, N) \
VEC = psVectorRecycle(VEC, N, PS_TYPE_F64); \
p_psMemSetPersistent(VEC, true); \
p_psMemSetPersistent(VEC->data.U8, true); \
for (int i=0;i<N;i++) { \
    VEC->data.F64[i] = (float) i; \
} \

#define PS_VECTOR_GEN_STATIC_RECYCLED(NAME, SIZE, TYPE) \
static psVector *NAME = NULL; \
(NAME) = psVectorRecycle((NAME), SIZE, TYPE); \
p_psMemSetPersistent((NAME), true); \
p_psMemSetPersistent((NAME)->data.U8, true); \

#define PS_VECTOR_DECLARE_ALLOC_STATIC(NAME, SIZE, TYPE) \
static psVector *(NAME) = NULL; \
if ((NAME) == NULL) { \
    (NAME) = psVectorAlloc(SIZE, TYPE); \
    p_psMemSetPersistent((NAME), true); \
} \

#define PS_VECTOR_SET_U8(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->n ; i++) { \
    (NAME)->data.U8[i] = VALUE; \
}\

#define PS_VECTOR_SET_U16(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->n ; i++) { \
    (NAME)->data.U16[i] = VALUE; \
}\

#define PS_VECTOR_SET_U32(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->n ; i++) { \
    (NAME)->data.U32[i] = VALUE; \
}\

#define PS_VECTOR_SET_U64(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->n ; i++) { \
    (NAME)->data.U64[i] = VALUE; \
}\

#define PS_VECTOR_SET_S8(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->n ; i++) { \
    (NAME)->data.S8[i] = VALUE; \
}\

#define PS_VECTOR_SET_S16(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->n ; i++) { \
    (NAME)->data.S16[i] = VALUE; \
}\

#define PS_VECTOR_SET_S32(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->n ; i++) { \
    (NAME)->data.S32[i] = VALUE; \
}\

#define PS_VECTOR_SET_S64(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->n ; i++) { \
    (NAME)->data.S64[i] = VALUE; \
}\

#define PS_VECTOR_SET_F64(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->n ; i++) { \
    (NAME)->data.F64[i] = VALUE; \
}\

#define PS_VECTOR_SET_F32(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->n ; i++) { \
    (NAME)->data.F32[i] = VALUE; \
}\


/*****************************************************************************
    PS_POLY macros:
*****************************************************************************/
#define PS_ASSERT_POLY_NON_NULL(NAME, RVAL) \
if ((NAME) == NULL || (NAME)->coeff == NULL) { \
    psError(PS_ERR_BAD_PARAMETER_NULL, true, \
            "Unallowable operation: polynomial %s or its coeffs is NULL.", \
            #NAME); \
    return(RVAL); \
} \

#define PS_ASSERT_POLY_TYPE(NAME, TYPE, RVAL) \
if ((NAME)->type != TYPE) { \
    psError(PS_ERR_BAD_PARAMETER_TYPE, true, \
            "Unallowable operation: polynomial %s has wrong type.", #NAME); \
    return(RVAL); \
} \

// The following macros declare and allocate a static polynomial of the
// specified order and type.

#define PS_POLY_1D_DECLARE_ALLOC_STATIC(NAME, ORDER, TYPE) \
static psPolynomial1D *(NAME) = NULL; \
if ((NAME) == NULL) { \
    (NAME) = psPolynomial1DAlloc(ORDER, TYPE); \
    p_psMemSetPersistent((NAME), true); \
    p_psMemSetPersistent((NAME)->coeff, true); \
    p_psMemSetPersistent((NAME)->coeffErr, true); \
    p_psMemSetPersistent((NAME)->mask, true); \
} \

#define PS_POLY_2D_DECLARE_ALLOC_STATIC(NAME, ORDER, TYPE) \
static psPolynomial2D *(NAME) = NULL; \
if ((NAME) == NULL) { \
    (NAME) = psPolynomial2DAlloc(ORDER, TYPE); \
    p_psMemSetPersistent((NAME), true); \
} \

#define PS_POLY_3D_DECLARE_ALLOC_STATIC(NAME, ORDER, TYPE) \
static psPolynomial3D *(NAME) = NULL; \
if ((NAME) == NULL) { \
    (NAME) = psPolynomial3DAlloc(ORDER, TYPE); \
    p_psMemSetPersistent((NAME), true); \
} \

#define PS_POLY_4D_DECLARE_ALLOC_STATIC(NAME, ORDER, TYPE) \
static psPolynomial4D *(NAME) = NULL; \
if ((NAME) == NULL) { \
    (NAME) = psPolynomial4DAlloc(ORDER, TYPE); \
    p_psMemSetPersistent((NAME), true); \
} \

#define PS_POLY_1D_D_DECLARE_ALLOC_STATIC(NAME, ORDER, TYPE) \
static psPolynomial1D *(NAME) = NULL; \
if ((NAME) == NULL) { \
    (NAME) = psPolynomial1DAlloc(ORDER, TYPE); \
    p_psMemSetPersistent((NAME), true); \
} \

#define PS_POLY_2D_D_DECLARE_ALLOC_STATIC(NAME, ORDER, TYPE) \
static psPolynomial2D *(NAME) = NULL; \
if ((NAME) == NULL) { \
    (NAME) = psPolynomial2DAlloc(ORDER, TYPE); \
    p_psMemSetPersistent((NAME), true); \
} \

#define PS_POLY_3D_D_DECLARE_ALLOC_STATIC(NAME, ORDER, TYPE) \
static psPolynomial3D *(NAME) = NULL; \
if ((NAME) == NULL) { \
    (NAME) = psPolynomial3DAlloc(ORDER, TYPE); \
    p_psMemSetPersistent((NAME), true); \
} \

#define PS_POLY_4D_D_DECLARE_ALLOC_STATIC(NAME, ORDER, TYPE) \
static psPolynomial4D *(NAME) = NULL; \
if ((NAME) == NULL) { \
    (NAME) = psPolynomial4DAlloc(ORDER, TYPE); \
    p_psMemSetPersistent((NAME), true); \
} \

/*****************************************************************************
    PS_IMAGE macros:
*****************************************************************************/
#define PS_ASSERT_IMAGE_NON_NULL(NAME, RVAL) PS_ASSERT_GENERAL_IMAGE_NON_NULL(NAME, return RVAL)
#define PS_ASSERT_GENERAL_IMAGE_NON_NULL(NAME, CLEANUP) \
if ((NAME) == NULL || (NAME)->data.V == NULL) { \
    psError(PS_ERR_BAD_PARAMETER_NULL, true, \
            "Unallowable operation: psImage %s or its data is NULL.", \
            #NAME); \
    CLEANUP; \
}

#define PS_ASSERT_IMAGE_NON_EMPTY(NAME, RVAL) PS_ASSERT_GENERAL_IMAGE_NON_EMPTY(NAME, return RVAL)
#define PS_ASSERT_GENERAL_IMAGE_NON_EMPTY(NAME, CLEANUP) \
if ((NAME)->numCols < 1 || (NAME)->numRows < 1) { \
    psError(PS_ERR_BAD_PARAMETER_SIZE, true, \
            "Unallowable operation: psImage %s has zero rows or columns (%dx%d).", \
            #NAME, (NAME)->numCols, (NAME)->numRows); \
    CLEANUP; \
}

#define PS_ASSERT_IMAGE_TYPE(NAME, TYPE, RVAL) \
if ((NAME)->type.type != TYPE) { \
    psError(PS_ERR_BAD_PARAMETER_TYPE, true, \
            "Unallowable operation: psImage %s has incorrect type.", \
            #NAME); \
    return(RVAL); \
}

#define PS_ASSERT_IMAGES_SIZE_EQUAL(NAME1, NAME2, RVAL) \
if (((NAME1)->numCols != (NAME2)->numCols) || \
        ((NAME1)->numRows != (NAME2)->numRows)) { \
    psError(PS_ERR_BAD_PARAMETER_SIZE, true, \
            "Unallowable operation: psImages %s and %s are not the same size.", \
            #NAME1, #NAME2); \
    return(RVAL); \
}

#define PS_ASSERT_IMAGE_SIZE(NAME1, NUM_COLS, NUM_ROWS, RVAL) \
if (((NAME1)->numCols != NUM_COLS) || \
        ((NAME1)->numRows != NUM_ROWS)) { \
    psError(PS_ERR_BAD_PARAMETER_SIZE, true, \
            "Unallowable operation: psImages %s is not the correct size.", \
            #NAME1); \
    return(RVAL); \
}

#define PS_IMAGE_PRINT_F32(NAME) \
printf("======== printing %s ========\n", #NAME); \
for (int i = 0 ; i < (NAME)->numRows ; i++) { \
    for (int j = 0 ; j < (NAME)->numCols ; j++) { \
        printf("%.2f ", (NAME)->data.F32[i][j]); \
    } \
    printf("\n"); \
}\

#define PS_IMAGE_SET_U8(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->numRows ; i++) { \
    for (int j = 0 ; j < (NAME)->numCols ; j++) { \
        (NAME)->data.U8[i][j] = (VALUE); \
    } \
}\

#define PS_IMAGE_SET_U16(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->numRows ; i++) { \
    for (int j = 0 ; j < (NAME)->numCols ; j++) { \
        (NAME)->data.U16[i][j] = (VALUE); \
    } \
}\

#define PS_IMAGE_SET_U32(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->numRows ; i++) { \
    for (int j = 0 ; j < (NAME)->numCols ; j++) { \
        (NAME)->data.U32[i][j] = (VALUE); \
    } \
}\

#define PS_IMAGE_SET_U64(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->numRows ; i++) { \
    for (int j = 0 ; j < (NAME)->numCols ; j++) { \
        (NAME)->data.U64[i][j] = (VALUE); \
    } \
}\

#define PS_IMAGE_SET_S8(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->numRows ; i++) { \
    for (int j = 0 ; j < (NAME)->numCols ; j++) { \
        (NAME)->data.S8[i][j] = (VALUE); \
    } \
}\

#define PS_IMAGE_SET_S16(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->numRows ; i++) { \
    for (int j = 0 ; j < (NAME)->numCols ; j++) { \
        (NAME)->data.S16[i][j] = (VALUE); \
    } \
}\

#define PS_IMAGE_SET_S32(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->numRows ; i++) { \
    for (int j = 0 ; j < (NAME)->numCols ; j++) { \
        (NAME)->data.S32[i][j] = (VALUE); \
    } \
}\

#define PS_IMAGE_SET_S64(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->numRows ; i++) { \
    for (int j = 0 ; j < (NAME)->numCols ; j++) { \
        (NAME)->data.S64[i][j] = (VALUE); \
    } \
}\

#define PS_IMAGE_SET_F32(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->numRows ; i++) { \
    for (int j = 0 ; j < (NAME)->numCols ; j++) { \
        (NAME)->data.F32[i][j] = (VALUE); \
    } \
}\

#define PS_IMAGE_SET_F64(NAME, VALUE) \
for (int i = 0 ; i < (NAME)->numRows ; i++) { \
    for (int j = 0 ; j < (NAME)->numCols ; j++) { \
        (NAME)->data.F64[i][j] = (VALUE); \
    } \
}\

/*****************************************************************************
    PS_READOUT macros:
*****************************************************************************/
#define PS_ASSERT_READOUT_NON_NULL(NAME, RVAL) \
if ((NAME) == NULL || (NAME)->image == NULL) { \
    psError(PS_ERR_BAD_PARAMETER_NULL, true, \
            "Unallowable operation: psReadout %s or its data is NULL.", \
            #NAME); \
    return(RVAL); \
}

#define PS_ASSERT_READOUT_NON_EMPTY(NAME, RVAL) \
if ((NAME)->image->numCols < 1 || (NAME)->image->numRows < 1) { \
    psError(PS_ERR_BAD_PARAMETER_SIZE, true, \
            "Unallowable operation: psReadout %s or its data is NULL.", #NAME); \
    return(RVAL); \
}

#define PS_ASSERT_READOUT_TYPE(NAME, TYPE, RVAL) \
if ((NAME)->image->type.type != TYPE) { \
    psError(PS_ERR_BAD_PARAMETER_TYPE, true, \
            "Unallowable operation: psImage %s has incorrect type.", #NAME); \
    return(RVAL); \
}

/*****************************************************************************
    Misc. macros:
 *****************************************************************************/
#define PS_MAX(A, B) \
(((A) > (B)) ? (A) : (B))

#define PS_MIN(A, B) \
(((A) < (B)) ? (A) : (B))

#define PS_SQR(A) \
((A) * (A))
