diff options
Diffstat (limited to 'ast/chebymap.c')
-rw-r--r-- | ast/chebymap.c | 2404 |
1 files changed, 0 insertions, 2404 deletions
diff --git a/ast/chebymap.c b/ast/chebymap.c deleted file mode 100644 index 29f544b..0000000 --- a/ast/chebymap.c +++ /dev/null @@ -1,2404 +0,0 @@ -/* To do: - - - what to do about input positions that fall outside the bounding box. - Have an attribute that can be used to select "set bad" or "extrapolate"? - - - what about overriding astRate ? - - - Providing an iterative inverse requires the Jacobian to be defined. - - for a PolyMap: if y = C.x^n then y' = n.C.x^n-1 - for a ChebyMap: if y = C.Tn(x) then y' = n.C.Un-1(x) - - where Un-1(x) is the Chebyshev polynomial of the second kind, degree - (n-1), evaluated at x. Since PolyMap.GetJacobian function uses PolyMaps - to express the Jacobian of a PolyMap, it would need to use ChebyMaps - to express the Jacobian of a ChebyMap. But this would mean that - ChebyMap needs to be able to represent Chebyshev polynomials of the - second type. In fact the set of powers associated with each coefficient - would need to indicate somehow whether to use type 1 or type 2 for - each power. Could use negative powers to indicate type 2, but PolyMap.StoreArrays - objects to negative powers. StoreArrays could just store them - without checking, and then call a virtual function to verify the powers - are OK. - - Simpler for the moment just to disable iterative inverses in - ChebyMap. - -*/ - - -/* -*class++ -* Name: -* ChebyMap - -* Purpose: -* Map coordinates using Chebyshev polynomial functions. - -* Constructor Function: -c astChebyMap -f AST_CHEBYMAP - -* Description: -* A ChebyMap is a form of Mapping which performs a Chebyshev polynomial -* transformation. Each output coordinate is a linear combination of -* Chebyshev polynomials of the first kind, of order zero up to a -* specified maximum order, evaluated at the input coordinates. The -* coefficients to be used in the linear combination are specified -* separately for each output coordinate. -* -* For a 1-dimensional ChebyMap, the forward transformation is defined -* as follows: -* -* f(x) = c0.T0(x') + c1.T1(x') + c2.T2(x') + ... -* -* where: -* - Tn(x') is the nth Chebyshev polynomial of the first kind: -* - T0(x') = 1 -* - T1(x') = x' -* - Tn+1(x') = 2.x'.Tn(x') + Tn-1(x') -* - x' is the inpux axis value, x, offset and scaled to the range -* [-1, 1] as x ranges over a specified bounding box, given when the -* ChebyMap is created. The input positions, x, supplied to the -* forward transformation must fall within the bounding box - bad -* axis values (AST__BAD) are generated for points outside the -* bounding box. -* -* For an N-dimensional ChebyMap, the forward transformation is a -* generalisation of the above form. Each output axis value is the sum -c of "ncoeff" -f of NCOEFF -* terms, where each term is the product of a single coefficient -* value and N factors of the form Tn(x'_i), where "x'_i" is the -* normalised value of the i'th input axis value. -* -* The forward and inverse transformations are defined independantly -* by separate sets of coefficients, supplied when the ChebyMap is -* created. If no coefficients are supplied to define the inverse -* transformation, the -c astPolyTran -f AST_POLYTRAN -* method of the parent PolyMap class can instead be used to create an -* inverse transformation. The inverse transformation so generated -* will be a Chebyshev polynomial with coefficients chosen to minimise -* the residuals left by a round trip (forward transformation followed -* by inverse transformation). - -* Inheritance: -* The ChebyMap class inherits from the PolyMap class. - -* Attributes: -* The ChebyMap class does not define any new attributes beyond those -* which are applicable to all PolyMaps. - -* Functions: -c In addition to those functions applicable to all PolyMap, the -c following functions may also be applied to all ChebyMaps: -f In addition to those routines applicable to all PolyMap, the -f following routines may also be applied to all ChebyMaps: -* -c - astChebyDomain: Get the bounds of the domain of the ChebyMap -f - AST_CHEBYDOMAIN: Get the bounds of the domain of the ChebyMap - -* Copyright: -* Copyright (C) 2017 East Asian Observatory. -* All Rights Reserved. - -* Licence: -* This program is free software: you can redistribute it and/or -* modify it under the terms of the GNU Lesser General Public -* License as published by the Free Software Foundation, either -* version 3 of the License, or (at your option) any later -* version. -* -* This program is distributed in the hope that it will be useful, -* but WITHOUT ANY WARRANTY; without even the implied warranty of -* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -* GNU Lesser General Public License for more details. -* -* You should have received a copy of the GNU Lesser General -* License along with this program. If not, see -* <http://www.gnu.org/licenses/>. - -* Authors: -* DSB: D.S. Berry (EAO) - -* History: -* 1-MAR-2017 (DSB): -* Original version. -* 30-MAR-2017 (DSB): -* Over-ride the astFitPoly1DInit and astFitPoly2DInit virtual -* functions inherited form the PolyMap class. -*class-- -*/ - -/* Module Macros. */ -/* ============== */ -/* Set the name of the class we are implementing. This indicates to - the header files that define class interfaces that they should make - "protected" symbols available. */ -#define astCLASS ChebyMap - -/* Include files. */ -/* ============== */ -/* Interface definitions. */ -/* ---------------------- */ - -#include "globals.h" /* Thread-safe global data access */ -#include "error.h" /* Error reporting facilities */ -#include "memory.h" /* Memory allocation facilities */ -#include "object.h" /* Base Object class */ -#include "pointset.h" /* Sets of points/coordinates */ -#include "polymap.h" /* Polynomial mappings (parent class) */ -#include "cmpmap.h" /* Compound mappings */ -#include "chebymap.h" /* Interface definition for this class */ -#include "unitmap.h" /* Unit mappings */ - -/* Error code definitions. */ -/* ----------------------- */ -#include "ast_err.h" /* AST error codes */ - -/* C header files. */ -/* --------------- */ -#include <ctype.h> -#include <math.h> -#include <stdio.h> -#include <stdlib.h> -#include <string.h> -#include <limits.h> -#include <float.h> - -/* Module Variables. */ -/* ================= */ - -/* Address of this static variable is used as a unique identifier for - member of this class. */ -static int class_check; - -/* Pointers to parent class methods which are extended by this class. */ -static int (* parent_getobjsize)( AstObject *, int * ); -static int (* parent_equal)( AstObject *, AstObject *, int * ); -static void (* parent_polypowers)( AstPolyMap *, double **, int, const int *, double **, int, int, int * ); -static AstPolyMap *(*parent_polytran)( AstPolyMap *, int, double, double, int, const double *, const double *, int * ); - - -#ifdef THREAD_SAFE -/* Define how to initialise thread-specific globals. */ -#define GLOBAL_inits \ - globals->Class_Init = 0; \ - -/* Create the function that initialises global data for this module. */ -astMAKE_INITGLOBALS(ChebyMap) - -/* Define macros for accessing each item of thread specific global data. */ -#define class_init astGLOBAL(ChebyMap,Class_Init) -#define class_vtab astGLOBAL(ChebyMap,Class_Vtab) - -#include <pthread.h> - - -#else - -/* Define the class virtual function table and its initialisation flag - as static variables. */ -static AstChebyMapVtab class_vtab; /* Virtual function table */ -static int class_init = 0; /* Virtual function table initialised? */ - -#endif - - -/* External Interface Function Prototypes. */ -/* ======================================= */ -/* The following functions have public prototypes only (i.e. no - protected prototypes), so we must provide local prototypes for use - within this module. */ -AstChebyMap *astChebyMapId_( int, int, int, const double[], int, const double[], - const double[], const double[], const double[], - const double[], const char *, ... ); - - -/* Prototypes for Private Member Functions. */ -/* ======================================== */ -static AstPolyMap *PolyTran( AstPolyMap *, int, double, double, int, const double *, const double *, int * ); -static int Equal( AstObject *, AstObject *, int * ); -static int GetIterInverse( AstPolyMap *, int * ); -static int GetObjSize( AstObject *, int * ); -static void ChebyDomain( AstChebyMap *, int, double *, double *, int * ); -static void Copy( const AstObject *, AstObject *, int * ); -static void Delete( AstObject *obj, int * ); -static void Dump( AstObject *, AstChannel *, int * ); -static void PolyPowers( AstPolyMap *, double **, int, const int *, double **, int, int, int *); -static void FitPoly1DInit( AstPolyMap *, int, double **, AstMinPackData *, double *, int *); -static void FitPoly2DInit( AstPolyMap *, int, double **, AstMinPackData *, double *, int *); - -/* Member functions. */ -/* ================= */ - -static void ChebyDomain( AstChebyMap *this, int forward, double *lbnd, - double *ubnd, int *status ){ -/* -*++ -* Name: -c astChebyDomain -f AST_CHEBYDOMAIN - -* Purpose: -* Returns the bounding box of the domain of a ChebyMap. - -* Type: -* Public virtual function. - -* Synopsis: -c #include "chebymap.h" -c void astChebyDomain( AstChebyMap *this, int forward, double *lbnd, -c double *ubnd ) -f CALL AST_CHEBYDOMAIN( THIS, FORWARD, LBND, UBND, STATUS ) - -* Class Membership: -* ChebyMap method. - -* Description: -c This function -f This routine -* returns the upper and lower limits of the box defining the domain -* of either the forward or inverse transformation of a ChebyMap. These -* are the values that were supplied when the ChebyMap was created. - -* Parameters: -c this -f THIS = INTEGER (Given) -* Pointer to the ChebyMap. -c forward -f FORWARD = LOGICAL (Given) -c A non-zero -f A .TRUE. -* value indicates that the domain of the ChebyMap's -* forward transformation is to be returned, while a zero -* value indicates that the domain of the inverse transformation -* should be returned. -c lbnd -f LBND() = DOUBLE PRECISION (Returned) -c Pointer to an -f An -* array in which to return the lower axis bounds of the ChebyMap -* domain. The number of elements should be at least equal to the -* number of ChebyMap inputs (if -c "forward" is non-zero), or outputs (if "forward" is zero). -f FORWARD is .TRUE.), or outputs (if FORWARD is .FALSE.). -c ubnd -f UBND() = DOUBLE PRECISION (Returned) -c Pointer to an -f An -* array in which to return the upper axis bounds of the ChebyMap -* domain. The number of elements should be at least equal to the -* number of ChebyMap inputs (if -c "forward" is non-zero), or outputs (if "forward" is zero). -f FORWARD is .TRUE.), or outputs (if FORWARD is .FALSE.). -f STATUS = INTEGER (Given and Returned) -f The global status. - -* Notes: -* - If the requested transformation is undefined (i.e. no -* transformation coefficients were specified when the ChebyMap was -* created), this method returns a box determined using the -c astMapBox -f AST_MAPBOX -* method on the opposite transformation, if the opposite -* transformation is defined. -* - If the above procedure fails to determine a bounding box, the supplied -* arrays are filled with AST__BAD values but no error is reported. - -*-- -*/ - -/* Local Variables: */ - double *lbnd_o; - double *offset_o; - double *offset; - double *scale_o; - double *scale; - double *ubnd_o; - int fwd_o; - int iax; - int nax; - int nax_o; - -/* Check the inherited status. */ - if( !astOK ) return; - -/* Get the scales and offsets to use, depending on the value of "forward" - and whether the ChebyMap has been inverted. */ - if( forward != astGetInvert( this ) ) { - scale = this->scale_f; - offset = this->offset_f; - nax = astGetNin( this ); - scale_o = this->scale_i; - offset_o = this->offset_i; - nax_o = astGetNout( this ); - fwd_o = 0; - } else { - scale = this->scale_i; - offset = this->offset_i; - nax = astGetNout( this ); - scale_o = this->scale_f; - offset_o = this->offset_f; - nax_o = astGetNin( this ); - fwd_o = 1; - } - -/* Check the domain is defined. */ - if( scale && offset ) { - for( iax = 0; iax < nax; iax++ ) { - if( scale[ iax ] != 0.0 ) { - lbnd[ iax ] = ( -1.0 - offset[ iax ] ) / scale[ iax ]; - ubnd[ iax ] = ( 1.0 - offset[ iax ] ) / scale[ iax ]; - } else { - lbnd[ iax ] = AST__BAD; - ubnd[ iax ] = AST__BAD; - } - } - -/* If the requested domain is not defined, see if it can be determined - by transforming the domain of the other transformation into the - requested input ot putput space. */ - } else if( scale_o && offset_o ){ - -/* Allocate memory to hold the bounding box in the other space (input or - output), and then store the bounding box values. */ - lbnd_o = astMalloc( nax_o*sizeof( *lbnd_o ) ); - ubnd_o = astMalloc( nax_o*sizeof( *ubnd_o ) ); - if( astOK ) { - for( iax = 0; iax < nax_o; iax++ ) { - if( scale_o[ iax ] != 0.0 ) { - lbnd_o[ iax ] = ( -1.0 - offset_o[ iax ] ) / scale_o[ iax ]; - ubnd_o[ iax ] = ( 1.0 - offset_o[ iax ] ) / scale_o[ iax ]; - } else { - lbnd_o[ iax ] = AST__BAD; - ubnd_o[ iax ] = AST__BAD; - } - } - -/* Loop round finding the bounds on each input axis of the requested - transformation. */ - for( iax = 0; iax < nax; iax++ ) { - astMapBox( this, lbnd_o, ubnd_o, fwd_o, iax, lbnd + iax, - ubnd + iax, NULL, NULL ); - } - -/* Free resources */ - lbnd_o = astFree( lbnd_o ); - ubnd_o = astFree( ubnd_o ); - } - - -/* If the domain of the other transformation is not defined, return bad values. */ - } else { - for( iax = 0; iax < nax; iax++ ) { - lbnd[ iax ] = AST__BAD; - ubnd[ iax ] = AST__BAD; - } - } -} - -static int Equal( AstObject *this_object, AstObject *that_object, int *status ) { -/* -* Name: -* Equal - -* Purpose: -* Test if two ChebyMaps are equivalent. - -* Type: -* Private function. - -* Synopsis: -* #include "chebymap.h" -* int Equal( AstObject *this, AstObject *that, int *status ) - -* Class Membership: -* ChebyMap member function (over-rides the astEqual protected -* method inherited from the astPolyMap class). - -* Description: -* This function returns a boolean result (0 or 1) to indicate whether -* two ChebyMaps are equivalent. - -* Parameters: -* this -* Pointer to the first Object (a ChebyMap). -* that -* Pointer to the second Object. -* status -* Pointer to the inherited status variable. - -* Returned Value: -* One if the ChebyMaps are equivalent, zero otherwise. - -* Notes: -* - A value of zero will be returned if this function is invoked -* with the global status set, or if it should fail for any reason. -*/ - -/* Local Variables: */ - AstChebyMap *that; - AstChebyMap *this; - int i; - int nin; - int nout; - int result; - -/* Initialise. */ - result = 0; - -/* Check the global error status. */ - if ( !astOK ) return result; - -/* Invoke the Equal method inherited from the parent PolyMap class. This - checks that the PolyMaps are equal. */ - result = (*parent_equal)( this_object, that_object, status ); - if( result ) { - -/* Obtain pointers to the two ChebyMap structures. */ - this = (AstChebyMap *) this_object; - that = (AstChebyMap *) that_object; - -/* Check the second object is a ChebyMap. We know the first is a - ChebyMap since we have arrived at this implementation of the virtual - function. */ - if( astIsAChebyMap( that ) ) { - -/* Get the number of axes in the input bounding box (the original input space). */ - nin = astGetInvert( this ) ? astGetNout( this ) : astGetNin( this ); - -/* Check the bounding box is the same for both ChebyMaps. */ - if( this->scale_f && that->scale_f && - this->offset_f && that->offset_f ) { - for( i = 0; i < nin && result; i++ ) { - if( !astEQUAL( this->scale_f[ i ], that->scale_f[ i ] ) || - !astEQUAL( this->offset_f[ i ], that->offset_f[ i ] )){ - result = 0; - } - } - } else if( this->scale_f || that->scale_f || - this->offset_f || that->offset_f ) { - result = 0; - } - -/* Get the number of axes in the output bounding box (the original output space). */ - nout = astGetInvert( this ) ? astGetNin( this ) : astGetNout( this ); - -/* Check the bounding box is the same for both ChebyMaps. */ - if( this->scale_i && that->scale_i && - this->offset_i && that->offset_i ) { - for( i = 0; i < nout && result; i++ ) { - if( !astEQUAL( this->scale_i[ i ], that->scale_i[ i ] ) || - !astEQUAL( this->offset_i[ i ], that->offset_i[ i ] )){ - result = 0; - } - } - } else if( this->scale_i || that->scale_i || - this->offset_i || that->offset_i ) { - result = 0; - } - } - } - -/* If an error occurred, clear the result value. */ - if ( !astOK ) result = 0; - -/* Return the result, */ - return result; -} - -static void FitPoly1DInit( AstPolyMap *this_polymap, int forward, double **table, - AstMinPackData *data, double *scales, int *status ){ -/* -* Name: -* FitPoly1DInit - -* Purpose: -* Perform initialisation needed for FitPoly1D - -* Type: -* Private function. - -* Synopsis: -* #include "chebymap.h" -* void FitPoly1DInit( AstPolyMap *this, int forward, double **table, -* AstMinPackData *data, double *scales, int *status ) - -* Class Membership: -* ChebyMap member function (over-rides the astFitPoly1DInit protected -* method inherited from the PolyMap class). - -* Description: -* This function performs initialisation needed for FitPoly1D in the -* PolyMap class. - -* Parameters: -* this -* Pointer to the PolyMap. -* forward -* Non-zero if the forward transformation of "this" is being -* replaced. Zero if the inverse transformation of "this" is being -* replaced. -* table -* Pointer to an array of 2 pointers. Each of these pointers points -* to an array of "nsamp" doubles, being the scaled and sampled values -* for x1 and y1 in that order. -* data -* Pointer to a structure holding information to pass the the -* service function invoked by the minimisation function. -* scales -* Array holding the scaling factors for the two columns of the table. -* Multiplying the table values by the scale factor produces PolyMap -* input or output axis values. The scales are modified on exit to -* take account of the scaling performed by the ChebyMap Transform -* method. -*/ - -/* Local Variables; */ - AstChebyMap *this; - double *px1; - double *pxp1; - double maxx; - double minx; - double off; - double pmax; - double pmin; - double scl; - double x1; - int k; - int w1; - -/* Check the local error status. */ - if ( !astOK ) return; - -/* Get a pointer to the ChebyMap structure. */ - this = (AstChebyMap *) this_polymap; - -/* Find the bounds of the supplied x1 values. */ - px1 = table[ 0 ]; - minx = *px1; - maxx = *px1; - px1++; - for( k = 1; k < data->nsamp; k++,px1++ ) { - if( *px1 > maxx ) { - maxx = *px1; - } else if( *px1 < minx ) { - minx = *px1; - } - } - -/* Transform the above limits from table values into PolyMap axis values. */ - pmax = maxx*scales[ 0 ]; - pmin = minx*scales[ 0 ]; - -/* Calculate the scale and offset that map the above bounds onto the range - [-1,+1], and store them in the ChebyMap. */ - if( pmax != pmin ) { - scl = 2.0/( pmax - pmin ); - off = -( pmax + pmin )/( pmax - pmin ); - } else if( astOK ){ - astError( AST__BADBX, "astPolyTran(%s): New bounding box has zero width " - "on axis 1.", status, astGetClass(this)); - } - - if( forward ) { - this->scale_f = (double *) astFree( this->scale_f ); - this->offset_f = (double *) astFree( this->offset_f ); - - this->scale_f = (double *) astMalloc( sizeof( double ) ); - this->offset_f = (double *) astMalloc( sizeof( double ) ); - if( astOK ) { - this->scale_f[ 0 ] = scl; - this->offset_f[ 0 ] = off; - } - } else { - this->scale_i = (double *) astFree( this->scale_i ); - this->offset_i = (double *) astFree( this->offset_i ); - - this->scale_i = (double *) astMalloc( sizeof( double ) ); - this->offset_i = (double *) astMalloc( sizeof( double ) ); - if( astOK ) { - this->scale_i[ 0 ] = scl; - this->offset_i[ 0 ] = off; - } - } - -/* Get pointers to the supplied x1 values. */ - px1 = table[ 0 ]; - -/* Get pointers to the location for the next "power" of x1. Here "X to - the power N" is a metaphor for Tn(x). */ - pxp1 = data->xp1; - -/* Loop round all samples. */ - for( k = 0; k < data->nsamp; k++ ) { - -/* Get the current x1 value, and scale it into the range [-1,+1]. */ - x1 = *(px1++)*scl*scales[0] + off; - -/* Find all the required "powers" of x1 and store them in the "xp1" - component of the data structure. */ - *(pxp1++) = 1.0; - *(pxp1++) = x1; - for( w1 = 2; w1 < data->order; w1++,pxp1++ ) { - pxp1[ 0 ] = 2.0*x1*pxp1[ -1 ] - pxp1[ -2 ]; - } - } - -/* The scaling representing by the scales[0] value will be performed by - the astTransform method of the ChebyMap class, so reset teh scales[0] - value to unity, to avoid the scaling being applied twice. */ - scales[ 0 ] = 1.0; - -} - -static void FitPoly2DInit( AstPolyMap *this_polymap, int forward, double **table, - AstMinPackData *data, double *scales, int *status ){ -/* -* Name: -* FitPoly2DInit - -* Purpose: -* Perform initialisation needed for FitPoly2D - -* Type: -* Private function. - -* Synopsis: -* #include "chebymap.h" -* void FitPoly2DInit( AstPolyMap *this, int forward, double **table, -* AstMinPackData *data, double *scales, int *status ) - -* Class Membership: -* ChebyMap member function (over-rides the astFitPoly2DInit protected -* method inherited from the PolyMap class). - -* Description: -* This function performs initialisation needed for FitPoly2D in the -* PolyMap class.. - -* Parameters: -* this -* Pointer to the PolyMap. -* forward -* Non-zero if the forward transformation of "this" is being -* replaced. Zero if the inverse transformation of "this" is being -* replaced. -* table -* Pointer to an array of 4 pointers. Each of these pointers points -* to an array of "nsamp" doubles, being the scaled and sampled values -* for x1, x2, y1 or y2 in that order. -* data -* Pointer to a structure holding information to pass the the -* service function invoked by the minimisation function. -* scales -* Array holding the scaling factors for the four columns of the table. -* Multiplying the table values by the scale factor produces PolyMap -* input or output axis values. -*/ - -/* Local Variables; */ - AstChebyMap *this; - double *px1; - double *px2; - double *pxp1; - double *pxp2; - double maxx; - double maxy; - double minx; - double miny; - double off[ 2 ]; - double pxmax; - double pxmin; - double pymax; - double pymin; - double scl[ 2 ]; - double x1; - double x2; - int k; - int w1; - int w2; - -/* Check the local error status. */ - if ( !astOK ) return; - -/* Get a pointer to the ChebyMap structure. */ - this = (AstChebyMap *) this_polymap; - -/* Find the bounds of the supplied x1 and x2 values. */ - px1 = table[ 0 ]; - px2 = table[ 1 ]; - minx = *px1; - maxx = *px1; - miny = *px2; - maxy = *px2; - px1++; - px2++; - for( k = 1; k < data->nsamp; k++,px1++,px2++ ) { - if( *px1 > maxx ) { - maxx = *px1; - } else if( *px1 < minx ) { - minx = *px1; - } - if( *px2 > maxy ) { - maxy = *px2; - } else if( *px2 < miny ) { - miny = *px2; - } - } - -/* Transform the above limits from table values into PolyMap axis values. */ - pxmax = maxx*scales[ 0 ]; - pxmin = minx*scales[ 0 ]; - pymax = maxy*scales[ 1 ]; - pymin = miny*scales[ 1 ]; - -/* Calculate the scale and offset that map the above bounds onto the range - [-1,+1], and store them in the ChebyMap. */ - if( pxmax != pxmin && pymax != pymin ) { - scl[ 0 ] = 2.0/( pxmax - pxmin ); - off[ 0 ] = -( pxmax + pxmin )/( pxmax - pxmin ); - scl[ 1 ] = 2.0/( pymax - pymin ); - off[ 1 ] = -( pymax + pymin )/( pymax - pymin ); - } else if( astOK ){ - astError( AST__BADBX, "astPolyTran(%s): New bounding box has zero width " - "on or both axes.", status, astGetClass(this)); - } - - if( forward ) { - this->scale_f = (double *) astFree( this->scale_f ); - this->offset_f = (double *) astFree( this->offset_f ); - - this->scale_f = (double *) astMalloc( 2*sizeof( double ) ); - this->offset_f = (double *) astMalloc( 2*sizeof( double ) ); - if( astOK ) { - this->scale_f[ 0 ] = scl[ 0 ]; - this->offset_f[ 0 ] = off[ 0 ]; - this->scale_f[ 1 ] = scl[ 1 ]; - this->offset_f[ 1 ] = off[ 1 ]; - } - } else { - this->scale_i = (double *) astFree( this->scale_i ); - this->offset_i = (double *) astFree( this->offset_i ); - - this->scale_i = (double *) astMalloc( 2*sizeof( double ) ); - this->offset_i = (double *) astMalloc( 2*sizeof( double ) ); - if( astOK ) { - this->scale_i[ 0 ] = scl[ 0 ]; - this->offset_i[ 0 ] = off[ 0 ]; - this->scale_i[ 1 ] = scl[ 1 ]; - this->offset_i[ 1 ] = off[ 1 ]; - } - } - -/* Get pointers to the supplied x1 and x2 values. */ - px1 = table[ 0 ]; - px2 = table[ 1 ]; - -/* Get pointers to the location for the next "power" of x1 anmd x2. Here "X to - the power N" is a metaphor for Tn(x). */ - pxp1 = data->xp1; - pxp2 = data->xp2; - -/* Loop round all samples. */ - for( k = 0; k < data->nsamp; k++ ) { - -/* Get the current x1 and x2 values, and scale them into the range [-1,+1]. */ - x1 = *(px1++)*scl[0]*scales[0] + off[0]; - x2 = *(px2++)*scl[1]*scales[1] + off[1]; - -/* Find all the required "powers" of x1 and store them in the "xp1" - component of the data structure. */ - *(pxp1++) = 1.0; - *(pxp1++) = x1; - for( w1 = 2; w1 < data->order; w1++,pxp1++ ) { - pxp1[ 0 ] = 2.0*x1*pxp1[ -1 ] - pxp1[ -2 ]; - } - -/* Find all the required "powers" of x2 and store them in the "xp2" - component of the data structure. */ - *(pxp2++) = 1.0; - *(pxp2++) = x2; - for( w2 = 2; w2 < data->order; w2++,pxp2++ ) { - pxp2[ 0 ] = 2.0*x2*pxp2[ -1 ] - pxp2[ -2 ]; - } - } - -/* The scaling representing by the scales[0] and scales[1] values will be - performed by the astTransform method of the ChebyMap class, so reset the - scales[0] and scales[1] values to unity, to avoid the scaling being - applied twice. */ - scales[ 0 ] = 1.0; - scales[ 1 ] = 1.0; - -} - -static int GetIterInverse( AstPolyMap *this, int *status ) { -/* -* Name: -* GetIterInverse - -* Purpose: -* Return the value of the IterInverse attribute. - -* Type: -* Private function. - -* Synopsis: -* #include "polymap.h" -* int GetIterInverse( AstObject *this, int *status ) - -* Class Membership: -* ChebyMap member function (over-rides the astGetIterInverse protected -* method inherited from the parent PolyMap class). - -* Description: -* This function returns the value of the IterInverse attribute, which -* is always zero for a ChebyMap. - -* Parameters: -* this -* Pointer to the ChebyMap. -* status -* Pointer to the inherited status variable. - -* Returned Value: -* The IterInverse value. - -* Notes: -* - A value of zero will be returned if this function is invoked -* with the global status set, or if it should fail for any reason. -*/ - - return 0; -} - -static int GetObjSize( AstObject *this_object, int *status ) { -/* -* Name: -* GetObjSize - -* Purpose: -* Return the in-memory size of an Object. - -* Type: -* Private function. - -* Synopsis: -* #include "chebymap.h" -* int GetObjSize( AstObject *this, int *status ) - -* Class Membership: -* ChebyMap member function (over-rides the astGetObjSize protected -* method inherited from the parent class). - -* Description: -* This function returns the in-memory size of the supplied ChebyMap, -* in bytes. - -* Parameters: -* this -* Pointer to the ChebyMap. -* status -* Pointer to the inherited status variable. - -* Returned Value: -* The Object size, in bytes. - -* Notes: -* - A value of zero will be returned if this function is invoked -* with the global status set, or if it should fail for any reason. -*/ - -/* Local Variables: */ - AstChebyMap *this; - int nin; - int nout; - int result; - -/* Initialise. */ - result = 0; - -/* Check the global error status. */ - if ( !astOK ) return result; - -/* Obtain a pointers to the ChebyMap structure. */ - this = (AstChebyMap *) this_object; - -/* Get the number of input and output axes. */ - nin = astGetInvert( this ) ? astGetNout( this ) : astGetNin( this ); - nout = astGetInvert( this ) ? astGetNin( this ) : astGetNout( this ); - -/* Invoke the GetObjSize method inherited from the parent class, and then - add on any components of the class structure defined by this class - which are stored in dynamically allocated memory. */ - result = (*parent_getobjsize)( this_object, status ); - - if( this->scale_f ) result += nin*sizeof( *this->scale_f ); - if( this->offset_f ) result += nin*sizeof( *this->offset_f ); - if( this->scale_i ) result += nout*sizeof( *this->scale_i ); - if( this->offset_i ) result += nout*sizeof( *this->offset_i ); - -/* If an error occurred, clear the result value. */ - if ( !astOK ) result = 0; - -/* Return the result, */ - return result; -} - -void astInitChebyMapVtab_( AstChebyMapVtab *vtab, const char *name, int *status ) { -/* -*+ -* Name: -* astInitChebyMapVtab - -* Purpose: -* Initialise a virtual function table for a ChebyMap. - -* Type: -* Protected function. - -* Synopsis: -* #include "chebymap.h" -* void astInitChebyMapVtab( AstChebyMapVtab *vtab, const char *name ) - -* Class Membership: -* ChebyMap vtab initialiser. - -* Description: -* This function initialises the component of a virtual function -* table which is used by the ChebyMap class. - -* Parameters: -* vtab -* Pointer to the virtual function table. The components used by -* all ancestral classes will be initialised if they have not already -* been initialised. -* name -* Pointer to a constant null-terminated character string which contains -* the name of the class to which the virtual function table belongs (it -* is this pointer value that will subsequently be returned by the Object -* astClass function). -*- -*/ - -/* Local Variables: */ - astDECLARE_GLOBALS /* Pointer to thread-specific global data */ - AstObjectVtab *object; /* Pointer to Object component of Vtab */ - AstPolyMapVtab *polymap; /* Pointer to PolyMap component of Vtab */ - -/* Check the local error status. */ - if ( !astOK ) return; - -/* Get a pointer to the thread specific global data structure. */ - astGET_GLOBALS(NULL); - -/* Initialize the component of the virtual function table used by the - parent class. */ - astInitPolyMapVtab( (AstPolyMapVtab *) vtab, name ); - -/* Store a unique "magic" value in the virtual function table. This - will be used (by astIsAChebyMap) to determine if an object belongs - to this class. We can conveniently use the address of the (static) - class_check variable to generate this unique value. */ - vtab->id.check = &class_check; - vtab->id.parent = &(((AstPolyMapVtab *) vtab)->id); - -/* Initialise member function pointers. */ -/* ------------------------------------ */ -/* Store pointers to the member functions (implemented here) that provide - virtual methods for this class. */ -/* none */ - -/* Save the inherited pointers to methods that will be extended, and - replace them with pointers to the new member functions. */ - object = (AstObjectVtab *) vtab; - polymap = (AstPolyMapVtab *) vtab; - - polymap->GetIterInverse = GetIterInverse; - - parent_getobjsize = object->GetObjSize; - object->GetObjSize = GetObjSize; - - parent_polypowers = polymap->PolyPowers; - polymap->PolyPowers = PolyPowers; - - parent_polytran = polymap->PolyTran; - polymap->PolyTran = PolyTran; - - parent_equal = object->Equal; - object->Equal = Equal; - - polymap->FitPoly1DInit = FitPoly1DInit; - polymap->FitPoly2DInit = FitPoly2DInit; - -/* Store pointers to the member functions (implemented here) that - provide virtual methods for this class. */ - vtab->ChebyDomain = ChebyDomain; - -/* Declare the destructor and copy constructor. */ - astSetDelete( (AstObjectVtab *) vtab, Delete ); - astSetCopy( (AstObjectVtab *) vtab, Copy ); - -/* Declare the class dump function. */ - astSetDump( vtab, Dump, "ChebyMap", "Chebyshev polynomial transformation" ); - -/* If we have just initialised the vtab for the current class, indicate - that the vtab is now initialised, and store a pointer to the class - identifier in the base "object" level of the vtab. */ - if( vtab == &class_vtab ) { - class_init = 1; - astSetVtabClassIdentifier( vtab, &(vtab->id) ); - } -} - -static void PolyPowers( AstPolyMap *this_polymap, double **work, int ncoord, - const int *mxpow, double **ptr, int point, int fwd, - int *status ){ -/* -* Name: -* PolyPowers - -* Purpose: -* Find the required powers of the input axis values. - -* Type: -* Private function. - -* Synopsis: -* #include "chebymap.h" -* void PolyPowers( AstPolyMap *this, double **work, int ncoord, -* const int *mxpow, double **ptr, int point, -* int fwd, int *status ) - -* Class Membership: -* ChebyMap member function (over-rides the astPolyPowers protected -* method inherited from the PolyMap class). - -* Description: -* This function is used by astTransform to calculate the powers of -* the axis values for a single input position. In the case of -* sub-classes, the powers may not be simply powers of the supplied -* axis values but may be more complex quantities such as a Chebyshev -* polynomial of the required degree evaluated at the input axis values. - -* Parameters: -* this -* Pointer to the PolyMap. -* work -* An array of "ncoord" pointers, each pointing to an array of -* length "max(2,mxpow)". The required values are placed in this -* array on exit. -* ncoord -* The number of axes. -* mxpow -* Pointer to an array holding the maximum power required of each -* axis value. Should have "ncoord" elements. -* ptr -* An array of "ncoord" pointers, each pointing to an array holding -* the axis values. Each of these arrays of axis values must have -* at least "point+1" elements. -* point -* The zero based index of the point within "ptr" that holds the -* axis values to be exponentiated. -* fwd -* Do the supplied coefficients define the foward transformation of -* the PolyMap? -*/ - -/* Local Variables; */ - AstChebyMap *this; - double *scales; - double *offsets; - double *pwork; - double *t; - double x; - int coord; - int ip; - -/* Check the local error status. */ - if ( !astOK ) return; - -/* Get a pointer to the ChebyMap structure. */ - this = (AstChebyMap *) this_polymap; - -/* Either transformation of a ChebyMap (forward or inverse) can be - defined either as Chebyshev polynomial or as a standard polynomial. - Chebyshev polynomials always have non-NULL scale array pointers. - If the scale array pointer is NULL, then the transformation is a - standard polynomial. If the coefficients relate to a standard - polynomial, then invoke the astPolyPowers implementation of the parent - class (PolyMap). */ - if( (fwd && !this->scale_f) || (!fwd && !this->scale_i) ) { - (*parent_polypowers)( this_polymap, work, ncoord, mxpow, ptr, point, - fwd, status ); - -/* If the coefficients relate to a Chebyshev polynomial... */ - } else { - scales = fwd ? this->scale_f : this->scale_i; - offsets = fwd ? this->offset_f : this->offset_i; - -/* This method uses a Chebyshev polynomial of the first kind of degree "i" - evaluated at "x'" instead of "x raised to the power i". Here, "x'" is - the input axis value scaled and shifted into the range [-1,+1] on each - axis. Loop over all input axes. */ - for( coord = 0; coord < ncoord; coord++ ) { - -/* Get a pointer to the array in which the powers of the current axis - value are to be returned. */ - pwork = work[ coord ]; - -/* The Chebyshev function (type 1) of degree zero is always 1.0, regardless - of the value of x. */ - pwork[ 0 ] = 1.0; - -/* Get the input axis value. If it is bad, store bad values for all - remaining powers. */ - x = ptr[ coord ][ point ]; - if( x == AST__BAD ) { - for( ip = 1; ip <= mxpow[ coord ]; ip++ ) pwork[ ip ] = AST__BAD; - -/* Otherwise, apply the required scaling to the input */ - } else { - x = x*scales[ coord ] + offsets[ coord ]; - -/* Return bad values for input positions outside the bounding box - associated with the transformation. */ - if( fabs( x ) <= 1.0 ) { - -/* The Chebyshev function of degree one is equal to x. */ - t = pwork + 1; - *t = x; - -/* Form and store the remaining Chebyshev polynomial values at the input axis value. - Use the standard recurrence relation: Tn+1(x') = 2.x'.Tn(x') + Tn-1(x'). */ - for( ip = 2; ip <= mxpow[ coord ]; ip++,t++ ) { - t[ 1 ] = 2.0*x*t[ 0 ] - t[ -1 ]; - } - } else { - for( ip = 1; ip <= mxpow[ coord ]; ip++ ) pwork[ ip ] = AST__BAD; - } - } - } - } -} - -static AstPolyMap *PolyTran( AstPolyMap *this_polymap, int forward, double acc, - double maxacc, int maxorder, const double *lbnd, - const double *ubnd, int *status ){ -/* -* Name: -* PolyTran - -* Purpose: -* Fit a PolyMap inverse or forward transformation. - -* Type: -* Private function. - -* Synopsis: -* #include "polymap.h" -* AstPolyMap *PolyTran( AstPolyMap *this, int forward, double acc, -* double maxacc, int maxorder, const double *lbnd, -* const double *ubnd ) - -* Class Membership: -* ChebyMap member function (over-rides the astPolyTran method inherited -* from the PolyMap class). - -* Description: -* This function creates a new PolyMap which is a copy of the supplied -* PolyMap, in which a specified transformation (forward or inverse) -* has been replaced by a new polynomial transformation. The -* coefficients of the new transformation are estimated by sampling -* the other transformation and performing a least squares polynomial -* fit in the opposite direction to the sampled positions and values. -* -* This method can only be used on (1-input,1-output) or (2-input,2-output) -* PolyMaps. -* -* The transformation to create is specified by the "forward" parameter. -* In what follows "X" refers to the inputs of the PolyMap, and "Y" to -* the outputs of the PolyMap. The forward transformation transforms -* input values (X) into output values (Y), and the inverse transformation -* transforms output values (Y) into input values (X). Within a PolyMap, -* each transformation is represented by an independent set of -* polynomials, P_f or P_i: Y=P_f(X) for the forward transformation and -* X=P_i(Y) for the inverse transformation. -* -* The "forward" parameter specifies the transformation to be replaced. If -* it is non-zero, a new forward transformation is created by first finding -* the input values (X) using the inverse transformation -* (which must be available) at a regular grid of points (Y) covering a -* rectangular region of the PolyMap's output space. The coefficients of -* the required forward polynomial, Y=P_f(X), are chosen in order to -* minimise the sum of the squared residuals between the sampled values -* of Y and P_f(X). -* -* If "forward" is zero (probably the most likely case), -* a new inverse transformation is created by -* first finding the output values (Y) using the forward transformation -* (which must be available) at a regular grid of points (X) covering a -* rectangular region of the PolyMap's input space. The coefficients of -* the required inverse polynomial, X=P_i(Y), are chosen in order to -* minimise the sum of the squared residuals between the sampled values -* of X and P_i(Y). -* -* This fitting process is performed repeatedly with increasing -* polynomial orders (starting with linear) until the target -* accuracy is achieved, or a specified maximum order is reached. If -* the target accuracy cannot be achieved even with this maximum-order -* polynomial, the best fitting maximum-order polynomial is returned so -* long as its accuracy is better than "maxacc". -* If it is not, an error is reported. - -* Parameters: -* this -* Pointer to the original Mapping. -* forward -* If non-zero, the forward PolyMap transformation is replaced. -* Otherwise the inverse transformation is replaced. -* acc -* The target accuracy, expressed as a geodesic distance within -* the PolyMap's input space (if "forward" is zero) or output -* space (if "forward" is non-zero). -* maxacc -* The maximum allowed accuracy for an acceptable polynomial, -* expressed as a geodesic distance within the PolyMap's input -* space (if "forward" is zero) or output space (if "forward" is -* non-zero). -* maxorder -* The maximum allowed polynomial order. This is one more than the -* maximum power of either input axis. So for instance, a value of -* 3 refers to a quadratic polynomial. Note, cross terms with total -* powers greater than or equal to maxorder are not inlcuded in the -* fit. So the maximum number of terms in each of the fitted -* polynomials is maxorder*(maxorder+1)/2. -* lbnd -* Pointer to an array holding the lower bounds of a rectangular -* region within the PolyMap's input space (if "forward" is zero) -* or output space (if "forward" is non-zero). The new polynomial -* will be evaluated over this rectangle. The length of this array -* should equal the value of the PolyMap's Nin or Nout attribute, -* depending on "forward". If a NULL pointer is supplied, the lower -* bounds of the box supplied when the ChebyMap was constructed is -* used. -* ubnd -* Pointer to an -* array holding the upper bounds of a rectangular region within -* the PolyMap's input space (if "forward" is zero) or output space -* (if "forward" is non-zero). The new polynomial will be evaluated -* over this rectangle. The length of this array should equal the -* value of the PolyMap's Nin or Nout attribute, depending on "forward". -* If a NULL pointer is supplied, the upper bounds of the box supplied -* when the ChebyMap was constructed is used. - -* Returned Value: -* astPolyTran() -* A pointer to the new PolyMap. A NULL pointer will be returned if -* the fit fails to achieve the accuracy specified by "maxacc", but -* no error will be reported. - -* Notes: -* - This function can only be used on 1D or 2D PolyMaps which have -* the same number of inputs and outputs. -* - A null Object pointer (AST__NULL) will be returned if this -* function is invoked with the AST error status set, or if it -* should fail for any reason. -*/ - -/* Local Variables: */ - AstChebyMap *this; - AstPolyMap *result; - const char *word; - double *offset; - double *scale; - double this_lbnd[ 2 ]; - double this_ubnd[ 2 ]; - int inverted; - int nax; - -/* Initialise. */ - result = NULL; - -/* Check the inherited status. */ - if ( !astOK ) return result; - -/* Get a pointer to the CHebyMap structure. */ - this = (AstChebyMap *) this_polymap; - -/* Select the ChebyMap scales and offsets to be used. */ - inverted = astGetInvert( this ); - if( ( inverted && !forward ) || ( !inverted && forward ) ) { - word = "inverse"; - scale = this->scale_i; - offset = this->offset_i; - nax = ((AstMapping *)this)->nout; - } else { - word = "forward"; - scale = this->scale_f; - offset = this->offset_f; - nax = ((AstMapping *)this)->nin; - } - -/* The scaled box for a Chebyshev polynomial spans [-1,+1] on each axis. - Create the corresponding unscaled box. If the user supplies any bounds, - use them in preference to the bounds in the ChebyMap. */ - if( lbnd ) { - this_lbnd[ 0 ] = lbnd[ 0 ]; - if( nax > 1 ) this_lbnd[ 1 ] = lbnd[ 1 ]; - - } else if( scale && offset ) { - this_lbnd[ 0 ] = ( -1.0 - offset[ 0 ] )/scale[ 0 ]; - if( nax > 1 ) this_lbnd[ 1 ] = ( -1.0 - offset[ 1 ] )/scale[ 1 ]; - - } else if( astOK ) { - astError( AST__NOBOX, "astPolyTran(%s): The %s transformation is " - "not a Chebyshev polynomial and therefore requires a " - "user-supplied bounding box. But no lower bounds were " - "supplied. ", status, astGetClass( this ), word ); - } - - - if( ubnd ) { - this_ubnd[ 0 ] = ubnd[ 0 ]; - if( nax > 1 ) this_ubnd[ 1 ] = ubnd[ 1 ]; - } else if( scale && offset ) { - this_ubnd[ 0 ] = ( 1.0 - offset[ 0 ] )/scale[ 0 ]; - if( nax > 1 ) this_ubnd[ 1 ] = ( 1.0 - offset[ 1 ] )/scale[ 1 ]; - } else if( astOK ) { - astError( AST__NOBOX, "astPolyTran(%s): The %s transformation is " - "not a Chebyshev polynomial and therefore requires a " - "user-supplied bounding box. But no upper bounds were " - "supplied. ", status, astGetClass( this ), word ); - } - -/* Invoke the parent astPolyMap method, using the bounding box selected - above. */ - result = (*parent_polytran)( this_polymap, forward, acc, maxacc, maxorder, - this_lbnd, this_ubnd, status ); - -/* Return the new ChebyMap. */ - return result; -} - - -/* Functions which access class attributes. */ -/* ---------------------------------------- */ -/* Implement member functions to access the attributes associated with - this class using the macros defined for this purpose in the - "object.h" file. For a description of each attribute, see the class - interface (in the associated .h file). */ - -/* Copy constructor. */ -/* ----------------- */ -static void Copy( const AstObject *objin, AstObject *objout, int *status ) { -/* -* Name: -* Copy - -* Purpose: -* Copy constructor for ChebyMap objects. - -* Type: -* Private function. - -* Synopsis: -* void Copy( const AstObject *objin, AstObject *objout, int *status ) - -* Description: -* This function implements the copy constructor for ChebyMap objects. - -* Parameters: -* objin -* Pointer to the object to be copied. -* objout -* Pointer to the object being constructed. -* status -* Pointer to the inherited status variable. - -* Returned Value: -* void - -* Notes: -* - This constructor makes a deep copy, including a copy of the -* coefficients associated with the input ChebyMap. -*/ - - -/* Local Variables: */ - AstChebyMap *in; /* Pointer to input ChebyMap */ - AstChebyMap *out; /* Pointer to output ChebyMap */ - int nin; /* No. of input coordinates */ - int nout; /* No. of output coordinates */ - -/* Check the global error status. */ - if ( !astOK ) return; - -/* Obtain pointers to the input and output ChebyMaps. */ - in = (AstChebyMap *) objin; - out = (AstChebyMap *) objout; - -/* Nullify the pointers stored in the output object since these will - currently be pointing at the input data (since the output is a simple - byte-for-byte copy of the input). Otherwise, the input data could be - freed by accidient if the output object is deleted due to an error - occuring in this function. */ - out->scale_f = NULL; - out->offset_f = NULL; - out->scale_i = NULL; - out->offset_i = NULL; - -/* Get the number of inputs and outputs of the uninverted Mapping. */ - nin = ( (AstMapping *) in )->nin; - nout = ( (AstMapping *) in )->nout; - -/* Copy the bounding box arrays. */ - if( in->scale_f ) out->scale_f = (double *) astStore( NULL, - (void *) in->scale_f, - sizeof( double )*nin ); - if( in->offset_f ) out->offset_f = (double *) astStore( NULL, - (void *) in->offset_f, - sizeof( double )*nin ); - if( in->scale_i ) out->scale_i = (double *) astStore( NULL, - (void *) in->scale_i, - sizeof( double )*nout ); - if( in->offset_i ) out->offset_i = (double *) astStore( NULL, - (void *) in->offset_i, - sizeof( double )*nout ); -} - -/* Destructor. */ -/* ----------- */ -static void Delete( AstObject *obj, int *status ) { -/* -* Name: -* Delete - -* Purpose: -* Destructor for ChebyMap objects. - -* Type: -* Private function. - -* Synopsis: -* void Delete( AstObject *obj, int *status ) - -* Description: -* This function implements the destructor for ChebyMap objects. - -* Parameters: -* obj -* Pointer to the object to be deleted. -* status -* Pointer to the inherited status variable. - -* Returned Value: -* void - -* Notes: -* This function attempts to execute even if the global error status is -* set. -*/ - -/* Local Variables: */ - AstChebyMap *this; - -/* Obtain a pointer to the ChebyMap structure. */ - this = (AstChebyMap *) obj; - -/* Free the boundib box arrays. */ - this->scale_f = astFree( this->scale_f ); - this->offset_f = astFree( this->offset_f ); - this->scale_i = astFree( this->scale_i ); - this->offset_i = astFree( this->offset_i ); -} - -/* Dump function. */ -/* -------------- */ -static void Dump( AstObject *this_object, AstChannel *channel, int *status ) { -/* -* Name: -* Dump - -* Purpose: -* Dump function for ChebyMap objects. - -* Type: -* Private function. - -* Synopsis: -* void Dump( AstObject *this, AstChannel *channel, int *status ) - -* Description: -* This function implements the Dump function which writes out data -* for the ChebyMap class to an output Channel. - -* Parameters: -* this -* Pointer to the ChebyMap whose data are being written. -* channel -* Pointer to the Channel to which the data are being written. -* status -* Pointer to the inherited status variable. -*/ - -#define KEY_LEN 50 /* Maximum length of a keyword */ - -/* Local Variables: */ - AstChebyMap *this; /* Pointer to the ChebyMap structure */ - char buff[ KEY_LEN + 1 ]; /* Buffer for keyword string */ - char comm[ 100 ]; /* Buffer for comment string */ - int i; /* Loop index */ - int nin; /* No. of input coords */ - int nout; /* No. of output coords */ - -/* Check the global error status. */ - if ( !astOK ) return; - -/* Obtain a pointer to the ChebyMap structure. */ - this = (AstChebyMap *) this_object; - -/* Find the number of inputs and outputs of the uninverted Mapping. */ - nin = ( (AstMapping *) this )->nin; - nout = ( (AstMapping *) this )->nout; - -/* Write out values representing the instance variables for the - ChebyMap class. */ - -/* The input axis scale factors. */ - if( this->scale_f ){ - for( i = 0; i < nin; i++ ){ - (void) sprintf( buff, "FSCL%d", i + 1 ); - (void) sprintf( comm, "Scale factor on input %d", i + 1 ); - astWriteDouble( channel, buff, 1, 1, (this->scale_f)[ i ], comm ); - } - } - -/* The input axis offsets. */ - if( this->offset_f ){ - for( i = 0; i < nin; i++ ){ - (void) sprintf( buff, "FOFF%d", i + 1 ); - (void) sprintf( comm, "Offset on input %d", i + 1 ); - astWriteDouble( channel, buff, 1, 1, (this->offset_f)[ i ], comm ); - } - } - -/* The output axis scale factors. */ - if( this->scale_i ){ - for( i = 0; i < nout; i++ ){ - (void) sprintf( buff, "ISCL%d", i + 1 ); - (void) sprintf( comm, "Scale factor on output %d", i + 1 ); - astWriteDouble( channel, buff, 1, 1, (this->scale_i)[ i ], comm ); - } - } - -/* The output axis offsets. */ - if( this->offset_i ){ - for( i = 0; i < nout; i++ ){ - (void) sprintf( buff, "IOFF%d", i + 1 ); - (void) sprintf( comm, "Offset on output %d", i + 1 ); - astWriteDouble( channel, buff, 1, 1, (this->offset_i)[ i ], comm ); - } - } - -/* Undefine macros local to this function. */ -#undef KEY_LEN -} - -/* Standard class functions. */ -/* ========================= */ -/* Implement the astIsAChebyMap and astCheckChebyMap functions using the macros - defined for this purpose in the "object.h" header file. */ -astMAKE_ISA(ChebyMap,Mapping) -astMAKE_CHECK(ChebyMap) - -AstChebyMap *astChebyMap_( int nin, int nout, int ncoeff_f, const double coeff_f[], - int ncoeff_i, const double coeff_i[], - const double lbnd_f[], const double ubnd_f[], - const double lbnd_i[], const double ubnd_i[], - const char *options, int *status, ...){ -/* -*++ -* Name: -c astChebyMap -f AST_CHEBYMAP - -* Purpose: -* Create a ChebyMap. - -* Type: -* Public function. - -* Synopsis: -c #include "chebymap.h" -c AstChebyMap *astChebyMap( int nin, int nout, int ncoeff_f, const double coeff_f[], -c int ncoeff_i, const double coeff_i[], -c const double lbnd_f[], const double ubnd_f[], -c const double lbnd_i[], const double ubnd_i[], -c const char *options, ... ) -f RESULT = AST_CHEBYMAP( NIN, NOUT, NCOEFF_F, COEFF_F, NCOEFF_I, COEFF_I, -f LBND_F, UBND_F, LBND_I, UBND_I, OPTIONS, STATUS ) - -* Class Membership: -* ChebyMap constructor. - -* Description: -* This function creates a new ChebyMap and optionally initialises -* its attributes. -* -* A ChebyMap is a form of Mapping which performs a Chebyshev polynomial -* transformation. Each output coordinate is a linear combination of -* Chebyshev polynomials of the first kind, of order zero up to a -* specified maximum order, evaluated at the input coordinates. The -* coefficients to be used in the linear combination are specified -* separately for each output coordinate. -* -* For a 1-dimensional ChebyMap, the forward transformation is defined -* as follows: -* -* f(x) = c0.T0(x') + c1.T1(x') + c2.T2(x') + ... -* -* where: -* - Tn(x') is the nth Chebyshev polynomial of the first kind: -* - T0(x') = 1 -* - T1(x') = x' -* - Tn+1(x') = 2.x'.Tn(x') + Tn-1(x') -* - x' is the inpux axis value, x, offset and scaled to the range -* [-1, 1] as x ranges over a specified bounding box, given when the -* ChebyMap is created. The input positions, x, supplied to the -* forward transformation must fall within the bounding box - bad -* axis values (AST__BAD) are generated for points outside the -* bounding box. -* -* For an N-dimensional ChebyMap, the forward transformation is a -* generalisation of the above form. Each output axis value is the sum -c of "ncoeff" -f of NCOEFF -* terms, where each term is the product of a single coefficient -* value and N factors of the form Tn(x'_i), where "x'_i" is the -* normalised value of the i'th input axis value. -* -* The forward and inverse transformations are defined independantly -* by separate sets of coefficients, supplied when the ChebyMap is -* created. If no coefficients are supplied to define the inverse -* transformation, the -c astPolyTran -f AST_POLYTRAN -* method of the parent PolyMap class can instead be used to create an -* inverse transformation. The inverse transformation so generated -* will be a Chebyshev polynomial with coefficients chosen to minimise -* the residuals left by a round trip (forward transformation followed -* by inverse transformation). - -* Parameters: -c nin -f NIN = INTEGER (Given) -* The number of input coordinates. -c nout -f NOUT = INTEGER (Given) -* The number of output coordinates. -c ncoeff_f -f NCOEFF_F = INTEGER (Given) -* The number of non-zero coefficients necessary to define the -* forward transformation of the ChebyMap. If zero is supplied, the -* forward transformation will be undefined. -c coeff_f -f COEFF_F( * ) = DOUBLE PRECISION (Given) -* An array containing -c "ncoeff_f*( 2 + nin )" elements. Each group of "2 + nin" -f "NCOEFF_F*( 2 + NIN )" elements. Each group of "2 + NIN" -* adjacent elements describe a single coefficient of the forward -* transformation. Within each such group, the first element is the -* coefficient value; the next element is the integer index of the -* ChebyMap output which uses the coefficient within its defining -* expression (the first output has index 1); the remaining elements -* of the group give the integer powers to use with each input -* coordinate value (powers must not be negative, and floating -* point values are rounded to the nearest integer). -c If "ncoeff_f" is zero, a NULL pointer may be supplied for "coeff_f". -* -* For instance, if the ChebyMap has 3 inputs and 2 outputs, each group -* consisting of 5 elements, A groups such as "(1.2, 2.0, 1.0, 3.0, 0.0)" -* describes a coefficient with value 1.2 which is used within the -* definition of output 2. The output value is incremented by the -* product of the coefficient value, the value of the Chebyshev -* polynomial of power 1 evaluated at input coordinate 1, and the -* value of the Chebyshev polynomial of power 3 evaluated at input -* coordinate 2. Input coordinate 3 is not used since its power is -* specified as zero. As another example, the group "(-1.0, 1.0, -* 0.0, 0.0, 0.0 )" adds a constant value -1.0 onto output 1 (it is -* a constant value since the power for every input axis is given as -* zero). -* -c Each final output coordinate value is the sum of the "ncoeff_f" terms -c described by the "ncoeff_f" groups within the supplied array. -f Each final output coordinate value is the sum of the "NCOEFF_F" terms -f described by the "NCOEFF_F" groups within the supplied array. -c ncoeff_i -f NCOEFF_I = INTEGER (Given) -* The number of non-zero coefficients necessary to define the -* inverse transformation of the ChebyMap. If zero is supplied, the -* inverse transformation will be undefined. -c coeff_i -f COEFF_I( * ) = DOUBLE PRECISION (Given) -* An array containing -c "ncoeff_i*( 2 + nout )" elements. Each group of "2 + nout" -f "NCOEFF_I*( 2 + NOUT )" elements. Each group of "2 + NOUT" -* adjacent elements describe a single coefficient of the inverse -c transformation, using the same schame as "coeff_f", -f transformation, using the same schame as "COEFF_F", -* except that "inputs" and "outputs" are transposed. -c If "ncoeff_i" is zero, a NULL pointer may be supplied for "coeff_i". -c lbnd_f -f LBND_F( * ) = DOUBLE PRECISION (Given) -* An array containing the lower bounds of the input bounding box within -* which the ChebyMap is defined. This argument is not used or -* accessed if -c ncoeff_f is zero, and so a NULL pointer may be supplied. -f NCOEFF_F is zero. -* If supplied, the array should contain -c "nin" elements. -f "NIN" elements. -c ubnd_f -f UBND_F( * ) = DOUBLE PRECISION (Given) -* An array containing the upper bounds of the input bounding box within -* which the ChebyMap is defined. This argument is not used or -* accessed if -c ncoeff_f is zero, and so a NULL pointer may be supplied. -f NCOEFF_F is zero. -* If supplied, the array should contain -c "nin" elements. -f "NIN" elements. -c lbnd_i -f LBND_I( * ) = DOUBLE PRECISION (Given) -* An array containing the lower bounds of the output bounding box within -* which the ChebyMap is defined. This argument is not used or -* accessed if -c ncoeff_i is zero, and so a NULL pointer may be supplied. -f NCOEFF_I is zero. -* If supplied, the array should contain -c "nout" elements. -f "NOUT" elements. -c ubnd_i -f UBND_I( * ) = DOUBLE PRECISION (Given) -* An array containing the upper bounds of the output bounding box within -* which the ChebyMap is defined. This argument is not used or -* accessed if -c ncoeff_i is zero, and so a NULL pointer may be supplied. -f NCOEFF_I is zero. -* If supplied, the array should contain -c "nout" elements. -f "NOUT" elements. -c options -f OPTIONS = CHARACTER * ( * ) (Given) -c Pointer to a null-terminated string containing an optional -c comma-separated list of attribute assignments to be used for -c initialising the new ChebyMap. The syntax used is identical to -c that for the astSet function and may include "printf" format -c specifiers identified by "%" symbols in the normal way. -f A character string containing an optional comma-separated -f list of attribute assignments to be used for initialising the -f new ChebyMap. The syntax used is identical to that for the -f AST_SET routine. -c ... -c If the "options" string contains "%" format specifiers, then -c an optional list of additional arguments may follow it in -c order to supply values to be substituted for these -c specifiers. The rules for supplying these are identical to -c those for the astSet function (and for the C "printf" -c function). -f STATUS = INTEGER (Given and Returned) -f The global status. - -* Returned Value: -c astChebyMap() -f AST_CHEBYMAP = INTEGER -* A pointer to the new ChebyMap. - -* Notes: -* - A null Object pointer (AST__NULL) will be returned if this -c function is invoked with the AST error status set, or if it -f function is invoked with STATUS set to an error value, or if it -* should fail for any reason. -*-- -*/ - -/* Local Variables: */ - astDECLARE_GLOBALS /* Pointer to thread-specific global data */ - AstChebyMap *new; /* Pointer to new ChebyMap */ - va_list args; /* Variable argument list */ - -/* Check the global status. */ - if ( !astOK ) return NULL; - -/* Get a pointer to the thread specific global data structure. */ - astGET_GLOBALS(NULL); - -/* Initialise the ChebyMap, allocating memory and initialising the - virtual function table as well if necessary. */ - new = astInitChebyMap( NULL, sizeof( AstChebyMap ), !class_init, - &class_vtab, "ChebyMap", nin, nout, - ncoeff_f, coeff_f, ncoeff_i, coeff_i, - lbnd_f, ubnd_f, lbnd_i, ubnd_i ); - -/* If successful, note that the virtual function table has been - initialised. */ - if ( astOK ) { - class_init = 1; - -/* Obtain the variable argument list and pass it along with the options string - to the astVSet method to initialise the new ChebyMap's attributes. */ - va_start( args, status ); - astVSet( new, options, NULL, args ); - va_end( args ); - -/* If an error occurred, clean up by deleting the new object. */ - if ( !astOK ) new = astDelete( new ); - } - -/* Return a pointer to the new ChebyMap. */ - return new; -} - -AstChebyMap *astChebyMapId_( int nin, int nout, int ncoeff_f, const double coeff_f[], - int ncoeff_i, const double coeff_i[], - const double lbnd_f[], const double ubnd_f[], - const double lbnd_i[], const double ubnd_i[], - const char *options, ... ){ -/* -* Name: -* astChebyMapId_ - -* Purpose: -* Create a ChebyMap. - -* Type: -* Private function. - -* Synopsis: -* #include "chebymap.h" -* AstChebyMap *astChebyMap( int nin, int nout, int ncoeff_f, const double coeff_f[], -* int ncoeff_i, const double coeff_i[], const -* double lbnd_f[], const double ubnd_f[], -* double lbnd_i[], const double ubnd_i[], -* const char *options, ... ) - -* Class Membership: -* ChebyMap constructor. - -* Description: -* This function implements the external (public) interface to the -* astChebyMap constructor function. It returns an ID value (instead -* of a true C pointer) to external users, and must be provided -* because astChebyMap_ has a variable argument list which cannot be -* encapsulated in a macro (where this conversion would otherwise -* occur). -* -* The variable argument list also prevents this function from -* invoking astChebyMap_ directly, so it must be a re-implementation -* of it in all respects, except for the final conversion of the -* result to an ID value. - -* Parameters: -* As for astChebyMap_. - -* Returned Value: -* The ID value associated with the new ChebyMap. -*/ - -/* Local Variables: */ - astDECLARE_GLOBALS /* Pointer to thread-specific global data */ - AstChebyMap *new; /* Pointer to new ChebyMap */ - va_list args; /* Variable argument list */ - int *status; /* Pointer to inherited status value */ - -/* Get a pointer to the inherited status value. */ - status = astGetStatusPtr; - -/* Get a pointer to the thread specific global data structure. */ - astGET_GLOBALS(NULL); - -/* Check the global status. */ - if ( !astOK ) return NULL; - -/* Initialise the ChebyMap, allocating memory and initialising the - virtual function table as well if necessary. */ - new = astInitChebyMap( NULL, sizeof( AstChebyMap ), !class_init, - &class_vtab, "ChebyMap", nin, nout, - ncoeff_f, coeff_f, ncoeff_i, coeff_i, - lbnd_f, ubnd_f, lbnd_i, ubnd_i ); - -/* If successful, note that the virtual function table has been - initialised. */ - if ( astOK ) { - class_init = 1; - -/* Obtain the variable argument list and pass it along with the options string - to the astVSet method to initialise the new ChebyMap's attributes. */ - va_start( args, options ); - astVSet( new, options, NULL, args ); - va_end( args ); - -/* If an error occurred, clean up by deleting the new object. */ - if ( !astOK ) new = astDelete( new ); - } - -/* Return an ID value for the new ChebyMap. */ - return astMakeId( new ); -} - -AstChebyMap *astInitChebyMap_( void *mem, size_t size, int init, - AstChebyMapVtab *vtab, const char *name, - int nin, int nout, int ncoeff_f, const double coeff_f[], - int ncoeff_i, const double coeff_i[], - const double lbnd_f[], const double ubnd_f[], - const double lbnd_i[], const double ubnd_i[], - int *status ){ -/* -*+ -* Name: -* astInitChebyMap - -* Purpose: -* Initialise a ChebyMap. - -* Type: -* Protected function. - -* Synopsis: -* #include "chebymap.h" -* AstChebyMap *astInitChebyMap( void *mem, size_t size, int init, -* AstChebyMapVtab *vtab, const char *name, -* int nin, int nout, int ncoeff_f, -* const double coeff_f[], int ncoeff_i, -* const double coeff_i[] -* const double lbnd_f[], const double ubnd_f[], -* const double lbnd_i[], const double ubnd_i[] ) - -* Class Membership: -* ChebyMap initialiser. - -* Description: -* This function is provided for use by class implementations to initialise -* a new ChebyMap object. It allocates memory (if necessary) to accommodate -* the ChebyMap plus any additional data associated with the derived class. -* It then initialises a ChebyMap structure at the start of this memory. If -* the "init" flag is set, it also initialises the contents of a virtual -* function table for a ChebyMap at the start of the memory passed via the -* "vtab" parameter. - -* Parameters: -* mem -* A pointer to the memory in which the ChebyMap is to be initialised. -* This must be of sufficient size to accommodate the ChebyMap data -* (sizeof(ChebyMap)) plus any data used by the derived class. If a value -* of NULL is given, this function will allocate the memory itself using -* the "size" parameter to determine its size. -* size -* The amount of memory used by the ChebyMap (plus derived class data). -* This will be used to allocate memory if a value of NULL is given for -* the "mem" parameter. This value is also stored in the ChebyMap -* structure, so a valid value must be supplied even if not required for -* allocating memory. -* init -* A logical flag indicating if the ChebyMap's virtual function table is -* to be initialised. If this value is non-zero, the virtual function -* table will be initialised by this function. -* vtab -* Pointer to the start of the virtual function table to be associated -* with the new ChebyMap. -* name -* Pointer to a constant null-terminated character string which contains -* the name of the class to which the new object belongs (it is this -* pointer value that will subsequently be returned by the astGetClass -* method). -* nin -* The number of input coordinate values per point. This is the -* same as the number of columns in the matrix. -* nout -* The number of output coordinate values per point. This is the -* same as the number of rows in the matrix. -* ncoeff_f -* The number of non-zero coefficients necessary to define the -* forward transformation of the ChebyMap. If zero is supplied, the -* forward transformation will be undefined. -* coeff_f -* An array containing "ncoeff_f*( 2 + nin )" elements. Each group -* of "2 + nin" adjacent elements describe a single coefficient of -* the forward transformation. Within each such group, the first -* element is the coefficient value; the next element is the -* integer index of the ChebyMap output which uses the coefficient -* within its defining polynomial (the first output has index 1); -* the remaining elements of the group give the integer powers to -* use with each input coordinate value (powers must not be -* negative) -* -* For instance, if the ChebyMap has 3 inputs and 2 outputs, each group -* consisting of 5 elements, A groups such as "(1.2, 2.0, 1.0, 3.0, 0.0)" -* describes a coefficient with value 1.2 which is used within the -* definition of output 2. The output value is incremented by the -* product of the coefficient value, the value of input coordinate -* 1 raised to the power 1, and the value of input coordinate 2 raised -* to the power 3. Input coordinate 3 is not used since its power is -* specified as zero. As another example, the group "(-1.0, 1.0, -* 0.0, 0.0, 0.0 )" describes adds a constant value -1.0 onto -* output 1 (it is a constant value since the power for every input -* axis is given as zero). -* -* Each final output coordinate value is the sum of the "ncoeff_f" terms -* described by the "ncoeff_f" groups within the supplied array. -* ncoeff_i -* The number of non-zero coefficients necessary to define the -* inverse transformation of the ChebyMap. If zero is supplied, the -* inverse transformation will be undefined. -* coeff_i -* An array containing -* "ncoeff_i*( 2 + nout )" elements. Each group of "2 + nout" -* adjacent elements describe a single coefficient of the inverse -* transformation, using the same schame as "coeff_f", except that -* "inputs" and "outputs" are transposed. -* lbnd_f -* An array containing the lower bounds of the input bounding box within -* which the ChebyMap is defined. The array should contain "nin" elements. -* ubnd_f -* An array containing the upper bounds of the input bounding box within -* which the ChebyMap is defined. The array should contain "nin" elements. -* lbnd_i -* An array containing the lower bounds of the output bounding box within -* which the ChebyMap is defined. The array should contain "nout" elements. -* ubnd_i -* An array containing the upper bounds of the output bounding box within -* which the ChebyMap is defined. The array should contain "nout" elements. - -* Returned Value: -* A pointer to the new ChebyMap. - -* Notes: -* - A null pointer will be returned if this function is invoked with the -* global error status set, or if it should fail for any reason. -*- -*/ - -/* Local Variables: */ - AstChebyMap *new; - int i; - -/* Check the global status. */ - if ( !astOK ) return NULL; - -/* If necessary, initialise the virtual function table. */ - if ( init ) astInitChebyMapVtab( vtab, name ); - -/* Initialise a PolyMap structure (the parent class) as the first component - within the ChebyMap structure, allocating memory if necessary. */ - new = (AstChebyMap *) astInitPolyMap( mem, size, 0, - (AstPolyMapVtab *) vtab, name, - nin, nout, ncoeff_f, coeff_f, - ncoeff_i, coeff_i ); - if ( astOK ) { - -/* Initialise the ChebyMap data. */ -/* ---------------------------- */ - -/* First initialise the pointers in case of errors. */ - new->scale_f = NULL; - new->offset_f = NULL; - new->scale_i = NULL; - new->offset_i = NULL; - -/* Calculate the scales and offsets that map the supplied input bounding box - onto the range [-1,+1] on each input axis, and store them. */ - if( ncoeff_f > 0 ) { - new->scale_f = (double *) astMalloc( sizeof( double )*nin ); - new->offset_f = (double *) astMalloc( sizeof( double )*nin ); - if( astOK ) { - for( i = 0; i < nin; i++ ) { - if( ubnd_f[ i ] != lbnd_f[ i ] ) { - new->scale_f[ i ] = 2.0/( ubnd_f[ i ] - lbnd_f[ i ] ); - new->offset_f[ i ] = -( ubnd_f[ i ] + lbnd_f[ i ] )/( ubnd_f[ i ] - lbnd_f[ i ] ); - } else if( astOK ){ - astError( AST__BADBX, "astInitChebyMap(%s): Input bounding box " - "has zero width on input axis %d.", status, name, i + 1 ); - break; - } - } - } - } - -/* Calculate the scales and offsets that map the supplied output bounding box - onto the range [-1,+1] on each output axis, and store them. */ - if( ncoeff_i > 0 ) { - new->scale_i = (double *) astMalloc( sizeof( double )*nout ); - new->offset_i = (double *) astMalloc( sizeof( double )*nout ); - if( astOK ) { - for( i = 0; i < nout; i++ ) { - if( ubnd_i[ i ] != lbnd_i[ i ] ) { - new->scale_i[ i ] = 2.0/( ubnd_i[ i ] - lbnd_i[ i ] ); - new->offset_i[ i ] = -( ubnd_i[ i ] + lbnd_i[ i ] )/( ubnd_i[ i ] - lbnd_i[ i ] ); - } else if( astOK ){ - astError( AST__BADBX, "astInitChebyMap(%s): Output bounding box " - "has zero width on output axis %d.", status, name, i + 1 ); - break; - } - } - } - } - -/* If an error occurred, clean up by deleting the new ChebyMap. */ - if ( !astOK ) new = astDelete( new ); - } - -/* Return a pointer to the new ChebyMap. */ - return new; -} - -AstChebyMap *astLoadChebyMap_( void *mem, size_t size, - AstChebyMapVtab *vtab, const char *name, - AstChannel *channel, int *status ) { -/* -*+ -* Name: -* astLoadChebyMap - -* Purpose: -* Load a ChebyMap. - -* Type: -* Protected function. - -* Synopsis: -* #include "chebymap.h" -* AstChebyMap *astLoadChebyMap( void *mem, size_t size, -* AstChebyMapVtab *vtab, const char *name, -* AstChannel *channel ) - -* Class Membership: -* ChebyMap loader. - -* Description: -* This function is provided to load a new ChebyMap using data read -* from a Channel. It first loads the data used by the parent class -* (which allocates memory if necessary) and then initialises a -* ChebyMap structure in this memory, using data read from the input -* Channel. -* -* If the "init" flag is set, it also initialises the contents of a -* virtual function table for a ChebyMap at the start of the memory -* passed via the "vtab" parameter. - - -* Parameters: -* mem -* A pointer to the memory into which the ChebyMap is to be -* loaded. This must be of sufficient size to accommodate the -* ChebyMap data (sizeof(ChebyMap)) plus any data used by derived -* classes. If a value of NULL is given, this function will -* allocate the memory itself using the "size" parameter to -* determine its size. -* size -* The amount of memory used by the ChebyMap (plus derived class -* data). This will be used to allocate memory if a value of -* NULL is given for the "mem" parameter. This value is also -* stored in the ChebyMap structure, so a valid value must be -* supplied even if not required for allocating memory. -* -* If the "vtab" parameter is NULL, the "size" value is ignored -* and sizeof(AstChebyMap) is used instead. -* vtab -* Pointer to the start of the virtual function table to be -* associated with the new ChebyMap. If this is NULL, a pointer -* to the (static) virtual function table for the ChebyMap class -* is used instead. -* name -* Pointer to a constant null-terminated character string which -* contains the name of the class to which the new object -* belongs (it is this pointer value that will subsequently be -* returned by the astGetClass method). -* -* If the "vtab" parameter is NULL, the "name" value is ignored -* and a pointer to the string "ChebyMap" is used instead. - -* Returned Value: -* A pointer to the new ChebyMap. - -* Notes: -* - A null pointer will be returned if this function is invoked -* with the global error status set, or if it should fail for any -* reason. -*- -*/ - -#define KEY_LEN 50 /* Maximum length of a keyword */ - - astDECLARE_GLOBALS /* Pointer to thread-specific global data */ -/* Local Variables: */ - AstChebyMap *new; /* Pointer to the new ChebyMap */ - char buff[ KEY_LEN + 1 ]; /* Buffer for keyword string */ - int i; /* Loop index */ - int ngood; /* No. of non-bad values */ - int nin; /* No. of input coords */ - int nout; /* No. of output coords */ - -/* Get a pointer to the thread specific global data structure. */ - astGET_GLOBALS(channel); - -/* Initialise. */ - new = NULL; - -/* Check the global error status. */ - if ( !astOK ) return new; - -/* If a NULL virtual function table has been supplied, then this is - the first loader to be invoked for this ChebyMap. In this case the - ChebyMap belongs to this class, so supply appropriate values to be - passed to the parent class loader (and its parent, etc.). */ - if ( !vtab ) { - size = sizeof( AstChebyMap ); - vtab = &class_vtab; - name = "ChebyMap"; - -/* If required, initialise the virtual function table for this class. */ - if ( !class_init ) { - astInitChebyMapVtab( vtab, name ); - class_init = 1; - } - } - -/* Invoke the parent class loader to load data for all the ancestral - classes of the current one, returning a pointer to the resulting - partly-built ChebyMap. */ - new = astLoadPolyMap( mem, size, (AstPolyMapVtab *) vtab, name, - channel ); - - if ( astOK ) { - -/* Get the number of inputs and outputs for the uninverted Mapping. */ - nin = ( (AstMapping *) new )->nin; - nout = ( (AstMapping *) new )->nout; - -/* Initialise values */ - new->scale_f = NULL; - new->offset_f = NULL; - new->scale_i = NULL; - new->offset_i = NULL; - -/* Read input data. */ -/* ================ */ - -/* Request the input Channel to read all the input data appropriate to - this class into the internal "values list". */ - astReadClassData( channel, "ChebyMap" ); - -/* Is the forward transformation defined? */ - if( ((AstPolyMap *) new)->ncoeff_f ) { - -/* Allocate memory to hold the scales and offsets for the input axes. */ - new->scale_f = astMalloc( sizeof( double )*(size_t) nin ); - new->offset_f = astMalloc( sizeof( double )*(size_t) nin ); - if( astOK ) { - -/* Get the scale factors. */ - ngood = 0; - for( i = 0; i < nin; i++ ){ - (void) sprintf( buff, "fscl%d", i + 1 ); - (new->scale_f)[ i ] = astReadDouble( channel, buff, AST__BAD ); - if( (new->scale_f)[ i ] != AST__BAD ) ngood++; - } - -/* Get the offsets of the bounding box. */ - for( i = 0; i < nin; i++ ){ - (void) sprintf( buff, "foff%d", i + 1 ); - (new->offset_f)[ i ] = astReadDouble( channel, buff, AST__BAD ); - if( (new->offset_f)[ i ] != AST__BAD ) ngood++; - } - -/* The scale and offset values should all be AST__BAD if the transformation - is a standard polynomial. Anull the scale and offset arrays to - indicate this. */ - if( ngood == 0 ) { - new->scale_f = astFree( new->scale_f ); - new->offset_f = astFree( new->offset_f ); - -/* Otherwise, there should be no bad values. */ - } else if( ngood != 2*nin && astOK ) { - astError( AST__OBJIN, "astLoadChebyMap: insufficient scale " - "and offset values for the forward transformation " - "in loaded ChebyMap.", status ); - } - } - } - -/* Is the inverse transformation defined? */ - if( ((AstPolyMap *) new)->ncoeff_i ) { - -/* Allocate memory to hold the scales and offsets for the output axes. */ - new->scale_i = astMalloc( sizeof( double )*(size_t) nout ); - new->offset_i = astMalloc( sizeof( double )*(size_t) nout ); - if( astOK ) { - -/* Get the scale factors. */ - ngood = 0; - for( i = 0; i < nout; i++ ){ - (void) sprintf( buff, "iscl%d", i + 1 ); - (new->scale_i)[ i ] = astReadDouble( channel, buff, AST__BAD ); - if( (new->scale_i)[ i ] != AST__BAD ) ngood++; - } - -/* Get the offsets of the bounding box. */ - for( i = 0; i < nout; i++ ){ - (void) sprintf( buff, "ioff%d", i + 1 ); - (new->offset_i)[ i ] = astReadDouble( channel, buff, AST__BAD ); - if( (new->offset_i)[ i ] != AST__BAD ) ngood++; - } - -/* The scale and offset values should all be AST__BAD if the transformation - is a standard polynomial. Anull the scale and offset arrays to - indicate this. */ - if( ngood == 0 ) { - new->scale_i = astFree( new->scale_i ); - new->offset_i = astFree( new->offset_i ); - -/* Otherwise, there should be no bad values. */ - } else if( ngood != 2*nout && astOK ) { - astError( AST__OBJIN, "astLoadChebyMap: insufficient scale " - "and offset values for the inverse transformation " - "in loaded ChebyMap.", status ); - } - } - } - -/* If an error occurred, clean up by deleting the new ChebyMap. */ - if ( !astOK ) new = astDelete( new ); - } - -/* Return the new ChebyMap pointer. */ - return new; - -/* Undefine macros local to this function. */ -#undef KEY_LEN -} - -/* Virtual function interfaces. */ -/* ============================ */ -/* These provide the external interface to the virtual functions defined by - this class. Each simply checks the global error status and then locates and - executes the appropriate member function, using the function pointer stored - in the object's virtual function table (this pointer is located using the - astMEMBER macro defined in "object.h"). - - Note that the member function may not be the one defined here, as it may - have been over-ridden by a derived class. However, it should still have the - same interface. */ - - -void astChebyDomain_( AstChebyMap *this, int forward, double *lbnd, double *ubnd, int *status ){ - if ( !astOK ) return; - (**astMEMBER(this,ChebyMap,ChebyDomain))( this, forward, lbnd, ubnd, status ); -} - - |