/*
******************************************************************************** * sinc.c -- Sinc function objects * * User routines: * Pro_SincFunc() creates a sinc function object. * Pro_SincInfo() returns parameters of sinc function. * * Mark Showalter, PDS Ring-Moon Systems Node, March 1998 *******************************************************************************/ #include <stdio.h> #include <math.h> #include "profile.h" #include "fortran.h" #define PI 3.141592653589793238462643 /******************** * Type definitions * ********************/ typedef struct ZPRO_SINC_STRUCT { RL_INT4 flag; RL_FLT8 height, step, halfwidth, shoulder, scale, threshold, shoulder_scale; } ZPRO_SINC; /********************** * Symbol definitions * **********************/ #define ZPRO_SINC_FLAG 'sinc' /************************************ * Prototypes of internal functions * ************************************/ static RL_FLT8 ZPro_SincCalc RL_PROTO((RL_FLT8 x, RL_VOID *params)); static void ZPro_PrintSinc RL_PROTO((RL_VOID *params)); /* ******************************************************************************** * EXPORTED USER ROUTINES ******************************************************************************** *$ Component_name: * Pro_SincFunc (sinc.c) *$ Abstract: * This routine generates a sinc function object. A sinc function returns * sin(pi*x)/(pi*x), scaled by an envelope function that goes smoothly to * zero outside its domain. *$ Keywords: * PROFILE, FUNCTION, PSF * C, PUBLIC, SUBROUTINE *$ Declarations: * PRO_OBJECT *Pro_SincFunc(height, step, halfwidth, shoulder) * RL_FLT8 height, step, halfwidth, shoulder; *$ Inputs: * height peak y-value of function (at origin). * step x-step between origin/zeroes of function. * halfwidth x-value above which y-value is zero. * shoulder x-width of truncation, modeled as a quarter-sinusoid. *$ Outputs: * none *$ Returns: * pointer to a new sinc function object, or NULL on non-fatal error. *$ Detailed_description: * This routine generates a sinc function object. A sinc function returns * sin(pi*x)/(pi*x), scaled by an envelope function that goes smoothly to * zero outside its domain. The x and y coordinates are initially unnamed. * * The envelope function used is as follows: * 1 |x| < halfwidth - shoulder * 0 |x| > halfwidth * (1 + cos((|x|-(halfwidth-shoulder)) * pi/shoulder))/2 * otherwise * This particular form has the advantage of being smooth throughout the * domain. *$ External_references: * Profile toolkit *$ Side_effects: * Memory is allocated. *$ Examples: * This line of code creates a sinc function with height 2 and with zeros * spaced by unity, with a full width of 10 and no shoulder transition * region. * * sinc = Pro_SincFunc(2., 1., 5., 0.); *$ Error_handling: * Profile toolkit error handling is in effect. * * Conditions raised: * RL_MEMORY_ERROR on memory allocation error. *$ Limitations: * none *$ Author_and_institution: * Mark R. Showalter * NASA/Ames Research Center *$ Version_and_date: * 1.0: March 1998 *$ Change_history: * none *******************************************************************************/ PRO_OBJECT *Pro_SincFunc(height, step, halfwidth, shoulder) RL_FLT8 height, step, halfwidth, shoulder; { ZPRO_SINC *params; PRO_OBJECT *func; RL_INT4 paramsize; paramsize = sizeof(ZPRO_SINC); params = (ZPRO_SINC *) XRL_Malloc(paramsize); if (params == NULL) return NULL; params->flag = ZPRO_SINC_FLAG; params->height = height; params->step = step; params->halfwidth = halfwidth; params->shoulder = shoulder; params->scale = PI / step; params->threshold = halfwidth - shoulder; params->shoulder_scale = 0.; if (params->shoulder != 0.) params->shoulder_scale = PI / shoulder; func = Pro_SoftFunc(ZPro_SincCalc, -halfwidth, halfwidth, (RL_VOID *) params, paramsize, TRUE, ZPro_PrintSinc); XRL_Free((RL_VOID *) params); return func; } /* ******************************************************************************** *$ Component_name: * Pro_SincInfo (sinc.c) *$ Abstract: * This routine returns information about a sinc function object. *$ Keywords: * PROFILE, FUNCTION, PSF * C, PUBLIC, SUBROUTINE *$ Declarations: * RL_BOOL Pro_SincInfo(sinc, height, step, halfwidth, shoulder) * PRO_OBJECT *sinc; * RL_FLT8 *height, *step, *halfwidth, *shoulder; *$ Inputs: * sinc pointer to a sinc function object. *$ Outputs: * *height peak y-value of function (at origin); * not returned if height==NULL. * *step x-step between origin/zeroes of function; * not returned if height==NULL. * *halfwidth x-value above which y-value is zero; * not returned if height==NULL. * *shoulder x-width of truncation, modeled as a half-sinusoid; * not returned if height==NULL. *$ Returns: * TRUE if this is a sinc function; FALSE otherwise. *$ Detailed_description: * This routine returns information about a sinc function object. *$ External_references: * Profile toolkit *$ Side_effects: * none *$ Examples: * This snippet of code of code prints the step size in a sinc function. * * PRO_OBJECT *sinc; * RL_FLT8 step; * * Pro_SincInfo(sinc, NULL, &step, NULL, NULL); * printf("Step size = %g\n", step); *$ Error_handling: * Profile toolkit error handling is in effect. * * Conditions raised: * PRO_CLASS_ERROR if this is not a software function object. *$ Limitations: * none *$ Author_and_institution: * Mark R. Showalter * NASA/Ames Research Center *$ Version_and_date: * 1.0: March 1998 *$ Change_history: * none *******************************************************************************/ RL_BOOL Pro_SincInfo(sinc, height, step, halfwidth, shoulder) PRO_OBJECT *sinc; RL_FLT8 *height, *step, *halfwidth, *shoulder; { ZPRO_SINC *params; params = (ZPRO_SINC *) XPro_SoftParams(sinc); if (params == NULL) return FALSE; if (params->flag != ZPRO_SINC_FLAG) return FALSE; if (height != NULL) *height = params->height; if (step != NULL) *step = params->step; if (halfwidth != NULL) *halfwidth = params->halfwidth; if (shoulder != NULL) *shoulder = params->shoulder; return TRUE; } /* ******************************************************************************** * INTERNAL ROUTINES ******************************************************************************** * ZPro_SincCalc(x, params) * * This internal routine calculates the value of a sinc function. * * Inputs: * x location at which to evaluate function. * params pointer to the ZPRO_SINC data structure. * * Return: value of function at x. *******************************************************************************/ static RL_FLT8 ZPro_SincCalc(x, params) RL_FLT8 x; RL_VOID *params; { ZPRO_SINC *sincpar; RL_FLT8 arg, value; sincpar = (ZPRO_SINC *) params; /* Handle mathematical singularity at origin */ if (x == 0.) return sincpar->height; /* Evaluate sinc function */ if (x < 0.) x = -x; arg = x * sincpar->scale; value = sincpar->height * sin(arg)/arg; /* If we're not in shoulder region, return */ if (x <= sincpar->threshold) return value; /* Otherwise, scale by shoulder function and return */ return value * 0.5 * (1. + cos((x - sincpar->threshold) * sincpar->shoulder_scale)); } /* ******************************************************************************** * ZPro_PrintSinc(params) * * This internal routine prints information about a sinc function. * * Inputs: * params pointer to the ZPRO_SINC data structure. *******************************************************************************/ static void ZPro_PrintSinc(params) RL_VOID *params; { ZPRO_SINC *sincpar; sincpar = (ZPRO_SINC *) params; /* Print sinc object info... */ printf("\nSinc software function parameters...\n"); printf(" height = %#g\n", sincpar->height); printf(" step = %#g\n", sincpar->step); printf("halfwidth = %#g\n", sincpar->halfwidth); printf(" shoulder = %#g\n", sincpar->shoulder); } /* ******************************************************************************** * FORTRAN INTERFACE ROUTINES ******************************************************************************** *$ Component_name: * FPro_SincFunc (sinc.c) *$ Abstract: * This routine generates a sinc function object. A sinc function returns * sin(pi*x)/(pi*x), scaled by an envelope function that goes smoothly to * zero outside its domain. *$ Keywords: * PROFILE, FUNCTION, PSF * FORTRAN, PUBLIC, SUBROUTINE *$ Declarations: * integer*4 function FPro_SincFunc(height, step, halfwidth, shoulder) * real*8 height, step, halfwidth, shoulder *$ Inputs: * height peak y-value of function (at origin). * step x-step between origin/zeroes of function. * halfwidth x-value above which y-value is zero. * shoulder x-width of truncation, modeled as a quarter-sinusoid. *$ Outputs: * none *$ Returns: * FORTRAN pointer to a new sinc function object, or 0 on non-fatal error. *$ Detailed_description: * This routine generates a sinc function object. A sinc function returns * sin(pi*x)/(pi*x), scaled by an envelope function that goes smoothly to * zero outside its domain. The x and y coordinates are initially unnamed. * * The envelope function used is as follows: * 1 |x| < halfwidth - shoulder * 0 |x| > halfwidth * (1 + cos((|x|-(halfwidth-shoulder)) * pi/shoulder))/2 * otherwise * This particular form has the advantage of being smooth throughout the * domain. *$ External_references: * Profile toolkit *$ Side_effects: * Memory is allocated. *$ Examples: * This line of code creates a sinc function with height 2 and with zeros * spaced by unity, with a full width of 10 and no shoulder transition * region. * * sinc = Pro_SincFunc(2.d0, 1.d0, 5.d0, 0.d0) *$ Error_handling: * Profile toolkit error handling is in effect. * * Conditions raised: * RL_MEMORY_ERROR on memory allocation error. *$ Limitations: * none *$ Author_and_institution: * Mark R. Showalter * NASA/Ames Research Center *$ Version_and_date: * 1.0: March 1998 *$ Change_history: * none *******************************************************************************/ RL_INT4 FORTRAN_NAME(fpro_sincfunc) (height, step, halfwidth, shoulder) RL_FLT8 *height, *step, *halfwidth, *shoulder; { RL_VOID *ptr; RL_INT4 index; ptr = (RL_VOID *) Pro_SincFunc(*height, *step, *halfwidth, *shoulder); if (ptr == NULL) return 0; index = FORT_AddPointer(ptr); if (index == 0) Pro_FreeObject((PRO_OBJECT *) ptr); return index; } /* ******************************************************************************** *$ Component_name: * FPro_SincInfo (sinc.c) *$ Abstract: * This routine returns information about a sinc function object. *$ Keywords: * PROFILE, FUNCTION, PSF * FORTRAN, PUBLIC, SUBROUTINE *$ Declarations: * logical*4 function FPro_SincInfo(sinc, height, step, halfwidth, * shoulder) * integer*4 sinc * real*8 height, step, halfwidth, shoulder *$ Inputs: * sinc FORTRAN pointer to a sinc function object. *$ Outputs: * height peak y-value of function (at origin). * step x-step between origin/zeroes of function. * halfwidth x-value above which y-value is zero. * shoulder x-width of truncation, modeled as a half-sinusoid. *$ Returns: * .TRUE. if this is a sinc function; .FALSE. otherwise. *$ Detailed_description: * This routine returns information about a sinc function object. *$ External_references: * Profile toolkit *$ Side_effects: * none *$ Examples: * This snippet of code of code prints the step size in a sinc function. * * integer*4 sinc * real*8 height, step, halfwidth, shoulder * * call FPro_SincInfo(sinc, height, step, halfwidth, shoulder) * write(*,*) 'Step size = ', step *$ Error_handling: * Profile toolkit error handling is in effect. * * Conditions raised: * PRO_CLASS_ERROR if this is not a software function object. * FORTRAN_POINTER_ERROR if sinc is not a valid FORTRAN object pointer. *$ Limitations: * none *$ Author_and_institution: * Mark R. Showalter * NASA/Ames Research Center *$ Version_and_date: * 1.0: March 1998 *$ Change_history: * none *******************************************************************************/ RL_INT4 FORTRAN_NAME(fpro_sincinfo) (sinc, height, step, halfwidth, shoulder) RL_INT4 *sinc; RL_FLT8 *height, *step, *halfwidth, *shoulder; { RL_VOID *ptr; RL_BOOL status; /* Look up sinc function pointer */ ptr = FORT_GetPointer(*sinc); if (ptr == NULL) return 0; /* Return info */ status = Pro_SincInfo(ptr, height, step, halfwidth, shoulder); return status ? FTRUE:FFALSE; } /********************************************************************************/