WCSLIB 8.3
Loading...
Searching...
No Matches
wcs.h
Go to the documentation of this file.
1/*============================================================================
2 WCSLIB 8.3 - an implementation of the FITS WCS standard.
3 Copyright (C) 1995-2024, Mark Calabretta
4
5 This file is part of WCSLIB.
6
7 WCSLIB is free software: you can redistribute it and/or modify it under the
8 terms of the GNU Lesser General Public License as published by the Free
9 Software Foundation, either version 3 of the License, or (at your option)
10 any later version.
11
12 WCSLIB is distributed in the hope that it will be useful, but WITHOUT ANY
13 WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14 FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
15 more details.
16
17 You should have received a copy of the GNU Lesser General Public License
18 along with WCSLIB. If not, see http://www.gnu.org/licenses.
19
20 Author: Mark Calabretta, Australia Telescope National Facility, CSIRO.
21 http://www.atnf.csiro.au/people/Mark.Calabretta
22 $Id: wcs.h,v 8.3 2024/05/13 16:33:00 mcalabre Exp $
23*=============================================================================
24*
25* WCSLIB 8.3 - C routines that implement the FITS World Coordinate System
26* (WCS) standard. Refer to the README file provided with WCSLIB for an
27* overview of the library.
28*
29*
30* Summary of the wcs routines
31* ---------------------------
32* Routines in this suite implement the FITS World Coordinate System (WCS)
33* standard which defines methods to be used for computing world coordinates
34* from image pixel coordinates, and vice versa. The standard, and proposed
35* extensions for handling distortions, are described in
36*
37= "Representations of world coordinates in FITS",
38= Greisen, E.W., & Calabretta, M.R. 2002, A&A, 395, 1061 (WCS Paper I)
39=
40= "Representations of celestial coordinates in FITS",
41= Calabretta, M.R., & Greisen, E.W. 2002, A&A, 395, 1077 (WCS Paper II)
42=
43= "Representations of spectral coordinates in FITS",
44= Greisen, E.W., Calabretta, M.R., Valdes, F.G., & Allen, S.L.
45= 2006, A&A, 446, 747 (WCS Paper III)
46=
47= "Representations of distortions in FITS world coordinate systems",
48= Calabretta, M.R. et al. (WCS Paper IV, draft dated 2004/04/22),
49= available from http://www.atnf.csiro.au/people/Mark.Calabretta
50=
51= "Mapping on the HEALPix grid",
52= Calabretta, M.R., & Roukema, B.F. 2007, MNRAS, 381, 865 (WCS Paper V)
53=
54= "Representing the 'Butterfly' Projection in FITS -- Projection Code XPH",
55= Calabretta, M.R., & Lowe, S.R. 2013, PASA, 30, e050 (WCS Paper VI)
56=
57= "Representations of time coordinates in FITS -
58= Time and relative dimension in space",
59= Rots, A.H., Bunclark, P.S., Calabretta, M.R., Allen, S.L.,
60= Manchester, R.N., & Thompson, W.T. 2015, A&A, 574, A36 (WCS Paper VII)
61*
62* These routines are based on the wcsprm struct which contains all information
63* needed for the computations. The struct contains some members that must be
64* set by the user, and others that are maintained by these routines, somewhat
65* like a C++ class but with no encapsulation.
66*
67* wcsnpv(), wcsnps(), wcsini(), wcsinit(), wcssub(), wcsfree(), and wcstrim(),
68* are provided to manage the wcsprm struct, wcssize() computes its total size
69* including allocated memory, wcsenq() returns information about the state of
70* the struct, and wcsprt() prints its contents. wcscopy(), which does a deep
71* copy of one wcsprm struct to another, is defined as a preprocessor macro
72* function that invokes wcssub().
73*
74* wcsperr() prints the error message(s) (if any) stored in a wcsprm struct,
75* and the linprm, celprm, prjprm, spcprm, and tabprm structs that it contains.
76*
77* A setup routine, wcsset(), computes intermediate values in the wcsprm struct
78* from parameters in it that were supplied by the user. The struct always
79* needs to be set up by wcsset() but this need not be called explicitly -
80* refer to the explanation of wcsprm::flag.
81*
82* wcsp2s() and wcss2p() implement the WCS world coordinate transformations.
83* In fact, they are high level driver routines for the WCS linear,
84* logarithmic, celestial, spectral and tabular transformation routines
85* described in lin.h, log.h, cel.h, spc.h and tab.h.
86*
87* Given either the celestial longitude or latitude plus an element of the
88* pixel coordinate a hybrid routine, wcsmix(), iteratively solves for the
89* unknown elements.
90*
91* wcsccs() changes the celestial coordinate system of a wcsprm struct, for
92* example, from equatorial to galactic, and wcssptr() translates the spectral
93* axis. For example, a 'FREQ' axis may be translated into 'ZOPT-F2W' and vice
94* versa.
95*
96* wcslib_version() returns the WCSLIB version number.
97*
98* Quadcube projections:
99* ---------------------
100* The quadcube projections (TSC, CSC, QSC) may be represented in FITS in
101* either of two ways:
102*
103* a: The six faces may be laid out in one plane and numbered as follows:
104*
105= 0
106=
107= 4 3 2 1 4 3 2
108=
109= 5
110*
111* Faces 2, 3 and 4 may appear on one side or the other (or both). The
112* world-to-pixel routines map faces 2, 3 and 4 to the left but the
113* pixel-to-world routines accept them on either side.
114*
115* b: The "COBE" convention in which the six faces are stored in a
116* three-dimensional structure using a CUBEFACE axis indexed from
117* 0 to 5 as above.
118*
119* These routines support both methods; wcsset() determines which is being
120* used by the presence or absence of a CUBEFACE axis in ctype[]. wcsp2s()
121* and wcss2p() translate the CUBEFACE axis representation to the single
122* plane representation understood by the lower-level WCSLIB projection
123* routines.
124*
125*
126* wcsnpv() - Memory allocation for PVi_ma
127* ---------------------------------------
128* wcsnpv() sets or gets the value of NPVMAX (default 64). This global
129* variable controls the number of pvcard structs, for holding PVi_ma
130* keyvalues, that wcsini() should allocate space for. It is also used by
131* wcsinit() as the default value of npvmax.
132*
133* PLEASE NOTE: This function is not thread-safe.
134*
135* Given:
136* n int Value of NPVMAX; ignored if < 0. Use a value less
137* than zero to get the current value.
138*
139* Function return value:
140* int Current value of NPVMAX.
141*
142*
143* wcsnps() - Memory allocation for PSi_ma
144* ---------------------------------------
145* wcsnps() sets or gets the value of NPSMAX (default 8). This global variable
146* controls the number of pscard structs, for holding PSi_ma keyvalues, that
147* wcsini() should allocate space for. It is also used by wcsinit() as the
148* default value of npsmax.
149*
150* PLEASE NOTE: This function is not thread-safe.
151*
152* Given:
153* n int Value of NPSMAX; ignored if < 0. Use a value less
154* than zero to get the current value.
155*
156* Function return value:
157* int Current value of NPSMAX.
158*
159*
160* wcsini() - Default constructor for the wcsprm struct
161* ----------------------------------------------------
162* wcsini() is a thin wrapper on wcsinit(). It invokes it with npvmax,
163* npsmax, and ndpmax set to -1 which causes it to use the values of the
164* global variables NDPMAX, NPSMAX, and NDPMAX. It is thereby potentially
165* thread-unsafe if these variables are altered dynamically via wcsnpv(),
166* wcsnps(), and disndp(). Use wcsinit() for a thread-safe alternative in
167* this case.
168*
169*
170* wcsinit() - Default constructor for the wcsprm struct
171* -----------------------------------------------------
172* wcsinit() optionally allocates memory for arrays in a wcsprm struct and sets
173* all members of the struct to default values.
174*
175* PLEASE NOTE: every wcsprm struct should be initialized by wcsinit(),
176* possibly repeatedly. On the first invokation, and only the first
177* invokation, wcsprm::flag must be set to -1 to initialize memory management,
178* regardless of whether wcsinit() will actually be used to allocate memory.
179*
180* Given:
181* alloc int If true, allocate memory unconditionally for the
182* crpix, etc. arrays. Please note that memory is never
183* allocated by wcsinit() for the auxprm, tabprm, nor
184* wtbarr structs.
185*
186* If false, it is assumed that pointers to these arrays
187* have been set by the user except if they are null
188* pointers in which case memory will be allocated for
189* them regardless. (In other words, setting alloc true
190* saves having to initalize these pointers to zero.)
191*
192* naxis int The number of world coordinate axes. This is used to
193* determine the length of the various wcsprm vectors and
194* matrices and therefore the amount of memory to
195* allocate for them.
196*
197* Given and returned:
198* wcs struct wcsprm*
199* Coordinate transformation parameters.
200*
201* Note that, in order to initialize memory management,
202* wcsprm::flag should be set to -1 when wcs is
203* initialized for the first time (memory leaks may
204* result if it had already been initialized).
205*
206* Given:
207* npvmax int The number of PVi_ma keywords to allocate space for.
208* If set to -1, the value of the global variable NPVMAX
209* will be used. This is potentially thread-unsafe if
210* wcsnpv() is being used dynamically to alter its value.
211*
212* npsmax int The number of PSi_ma keywords to allocate space for.
213* If set to -1, the value of the global variable NPSMAX
214* will be used. This is potentially thread-unsafe if
215* wcsnps() is being used dynamically to alter its value.
216*
217* ndpmax int The number of DPja or DQia keywords to allocate space
218* for. If set to -1, the value of the global variable
219* NDPMAX will be used. This is potentially
220* thread-unsafe if disndp() is being used dynamically to
221* alter its value.
222*
223* Function return value:
224* int Status return value:
225* 0: Success.
226* 1: Null wcsprm pointer passed.
227* 2: Memory allocation failed.
228*
229* For returns > 1, a detailed error message is set in
230* wcsprm::err if enabled, see wcserr_enable().
231*
232*
233* wcsauxi() - Default constructor for the auxprm struct
234* -----------------------------------------------------
235* wcsauxi() optionally allocates memory for an auxprm struct, attaches it to
236* wcsprm, and sets all members of the struct to default values.
237*
238* Given:
239* alloc int If true, allocate memory unconditionally for the
240* auxprm struct.
241*
242* If false, it is assumed that wcsprm::aux has already
243* been set to point to an auxprm struct, in which case
244* the user is responsible for managing that memory.
245* However, if wcsprm::aux is a null pointer, memory will
246* be allocated regardless. (In other words, setting
247* alloc true saves having to initalize the pointer to
248* zero.)
249*
250* Given and returned:
251* wcs struct wcsprm*
252* Coordinate transformation parameters.
253*
254* Function return value:
255* int Status return value:
256* 0: Success.
257* 1: Null wcsprm pointer passed.
258* 2: Memory allocation failed.
259*
260*
261* wcssub() - Subimage extraction routine for the wcsprm struct
262* ------------------------------------------------------------
263* wcssub() extracts the coordinate description for a subimage from a wcsprm
264* struct. It does a deep copy, using wcsinit() to allocate memory for its
265* arrays if required. Only the "information to be provided" part of the
266* struct is extracted. Consequently, wcsset() need not have been, and won't
267* be invoked on the struct from which the subimage is extracted. A call to
268* wcsset() is required to set up the subimage struct.
269*
270* The world coordinate system of the subimage must be separable in the sense
271* that the world coordinates at any point in the subimage must depend only on
272* the pixel coordinates of the axes extracted. In practice, this means that
273* the linear transformation matrix of the original image must not contain
274* non-zero off-diagonal terms that associate any of the subimage axes with any
275* of the non-subimage axes. Likewise, if any distortions are associated with
276* the subimage axes, they must not depend on any of the axes that are not
277* being extracted.
278*
279* Note that while the required elements of the tabprm array are extracted, the
280* wtbarr array is not. (Thus it is not appropriate to call wcssub() after
281* wcstab() but before filling the tabprm structs - refer to wcshdr.h.)
282*
283* wcssub() can also add axes to a wcsprm struct. The new axes will be created
284* using the defaults set by wcsinit() which produce a simple, unnamed, linear
285* axis with world coordinate equal to the pixel coordinate. These default
286* values can be changed afterwards, before invoking wcsset().
287*
288* Given:
289* alloc int If true, allocate memory for the crpix, etc. arrays in
290* the destination. Otherwise, it is assumed that
291* pointers to these arrays have been set by the user
292* except if they are null pointers in which case memory
293* will be allocated for them regardless.
294*
295* wcssrc const struct wcsprm*
296* Struct to extract from.
297*
298* Given and returned:
299* nsub int*
300* axes int[] Vector of length *nsub containing the image axis
301* numbers (1-relative) to extract. Order is
302* significant; axes[0] is the axis number of the input
303* image that corresponds to the first axis in the
304* subimage, etc.
305*
306* Use an axis number of 0 to create a new axis using
307* the defaults set by wcsinit(). They can be changed
308* later.
309*
310* nsub (the pointer) may be set to zero, and so also may
311* *nsub, which is interpreted to mean all axes in the
312* input image; the number of axes will be returned if
313* nsub != 0x0. axes itself (the pointer) may be set to
314* zero to indicate the first *nsub axes in their
315* original order.
316*
317* Set both nsub (or *nsub) and axes to zero to do a deep
318* copy of one wcsprm struct to another.
319*
320* Subimage extraction by coordinate axis type may be
321* done by setting the elements of axes[] to the
322* following special preprocessor macro values:
323*
324* WCSSUB_LONGITUDE: Celestial longitude.
325* WCSSUB_LATITUDE: Celestial latitude.
326* WCSSUB_CUBEFACE: Quadcube CUBEFACE axis.
327* WCSSUB_SPECTRAL: Spectral axis.
328* WCSSUB_STOKES: Stokes axis.
329* WCSSUB_TIME: Time axis.
330*
331* Refer to the notes (below) for further usage examples.
332*
333* On return, *nsub will be set to the number of axes in
334* the subimage; this may be zero if there were no axes
335* of the required type(s) (in which case no memory will
336* be allocated). axes[] will contain the axis numbers
337* that were extracted, or 0 for newly created axes. The
338* vector length must be sufficient to contain all axis
339* numbers. No checks are performed to verify that the
340* coordinate axes are consistent, this is done by
341* wcsset().
342*
343* wcsdst struct wcsprm*
344* Struct describing the subimage. wcsprm::flag should
345* be set to -1 if wcsdst was not previously initialized
346* (memory leaks may result if it was previously
347* initialized).
348*
349* Function return value:
350* int Status return value:
351* 0: Success.
352* 1: Null wcsprm pointer passed.
353* 2: Memory allocation failed.
354* 12: Invalid subimage specification.
355* 13: Non-separable subimage coordinate system.
356*
357* For returns > 1, a detailed error message is set in
358* wcsprm::err if enabled, see wcserr_enable().
359*
360* Notes:
361* 1: Combinations of subimage axes of particular types may be extracted in
362* the same order as they occur in the input image by combining
363* preprocessor codes, for example
364*
365= *nsub = 1;
366= axes[0] = WCSSUB_LONGITUDE | WCSSUB_LATITUDE | WCSSUB_SPECTRAL;
367*
368* would extract the longitude, latitude, and spectral axes in the same
369* order as the input image. If one of each were present, *nsub = 3 would
370* be returned.
371*
372* For convenience, WCSSUB_CELESTIAL is defined as the combination
373* WCSSUB_LONGITUDE | WCSSUB_LATITUDE | WCSSUB_CUBEFACE.
374*
375* The codes may also be negated to extract all but the types specified,
376* for example
377*
378= *nsub = 4;
379= axes[0] = WCSSUB_LONGITUDE;
380= axes[1] = WCSSUB_LATITUDE;
381= axes[2] = WCSSUB_CUBEFACE;
382= axes[3] = -(WCSSUB_SPECTRAL | WCSSUB_STOKES);
383*
384* The last of these specifies all axis types other than spectral or
385* Stokes. Extraction is done in the order specified by axes[] a
386* longitude axis (if present) would be extracted first (via axes[0]) and
387* not subsequently (via axes[3]). Likewise for the latitude and cubeface
388* axes in this example.
389*
390* From the foregoing, it is apparent that the value of *nsub returned may
391* be less than or greater than that given. However, it will never exceed
392* the number of axes in the input image (plus the number of newly-created
393* axes if any were specified on input).
394*
395*
396* wcscompare() - Compare two wcsprm structs for equality
397* ------------------------------------------------------
398* wcscompare() compares two wcsprm structs for equality.
399*
400* Given:
401* cmp int A bit field controlling the strictness of the
402* comparison. When 0, all fields must be identical.
403*
404* The following constants may be or'ed together to
405* relax the comparison:
406* WCSCOMPARE_ANCILLARY: Ignore ancillary keywords
407* that don't change the WCS transformation, such
408* as DATE-OBS or EQUINOX.
409* WCSCOMPARE_TILING: Ignore integral differences in
410* CRPIXja. This is the 'tiling' condition, where
411* two WCSes cover different regions of the same
412* map projection and align on the same map grid.
413* WCSCOMPARE_CRPIX: Ignore any differences at all in
414* CRPIXja. The two WCSes cover different regions
415* of the same map projection but may not align on
416* the same map grid. Overrides WCSCOMPARE_TILING.
417*
418* tol double Tolerance for comparison of floating-point values.
419* For example, for tol == 1e-6, all floating-point
420* values in the structs must be equal to the first 6
421* decimal places. A value of 0 implies exact equality.
422*
423* wcs1 const struct wcsprm*
424* The first wcsprm struct to compare.
425*
426* wcs2 const struct wcsprm*
427* The second wcsprm struct to compare.
428*
429* Returned:
430* equal int* Non-zero when the given structs are equal.
431*
432* Function return value:
433* int Status return value:
434* 0: Success.
435* 1: Null pointer passed.
436*
437*
438* wcscopy() macro - Copy routine for the wcsprm struct
439* ----------------------------------------------------
440* wcscopy() does a deep copy of one wcsprm struct to another. As of
441* WCSLIB 3.6, it is implemented as a preprocessor macro that invokes
442* wcssub() with the nsub and axes pointers both set to zero.
443*
444*
445* wcsfree() - Destructor for the wcsprm struct
446* --------------------------------------------
447* wcsfree() frees memory allocated for the wcsprm arrays by wcsinit() and/or
448* wcsset(). wcsinit() records the memory it allocates and wcsfree() will only
449* attempt to free this.
450*
451* PLEASE NOTE: wcsfree() must not be invoked on a wcsprm struct that was not
452* initialized by wcsinit().
453*
454* Given and returned:
455* wcs struct wcsprm*
456* Coordinate transformation parameters.
457*
458* Function return value:
459* int Status return value:
460* 0: Success.
461* 1: Null wcsprm pointer passed.
462*
463*
464* wcstrim() - Free unused arrays in the wcsprm struct
465* ---------------------------------------------------
466* wcstrim() frees memory allocated by wcsinit() for arrays in the wcsprm
467* struct that remains unused after it has been set up by wcsset().
468*
469* The free'd array members are associated with FITS WCS keyrecords that are
470* rarely used and usually just bloat the struct: wcsprm::crota, wcsprm::colax,
471* wcsprm::cname, wcsprm::crder, wcsprm::csyer, wcsprm::czphs, and
472* wcsprm::cperi. If unused, wcsprm::pv, wcsprm::ps, and wcsprm::cd are also
473* freed.
474*
475* Once these arrays have been freed, a test such as
476=
477= if (!undefined(wcs->cname[i])) {...}
478=
479* must be protected as follows
480=
481= if (wcs->cname && !undefined(wcs->cname[i])) {...}
482=
483* In addition, if wcsprm::npv is non-zero but less than wcsprm::npvmax, then
484* the unused space in wcsprm::pv will be recovered (using realloc()).
485* Likewise for wcsprm::ps.
486*
487* Given and returned:
488* wcs struct wcsprm*
489* Coordinate transformation parameters.
490*
491* Function return value:
492* int Status return value:
493* 0: Success.
494* 1: Null wcsprm pointer passed.
495* 14: wcsprm struct is unset.
496*
497*
498* wcssize() - Compute the size of a wcsprm struct
499* -----------------------------------------------
500* wcssize() computes the full size of a wcsprm struct, including allocated
501* memory.
502*
503* Given:
504* wcs const struct wcsprm*
505* Coordinate transformation parameters.
506*
507* If NULL, the base size of the struct and the allocated
508* size are both set to zero.
509*
510* Returned:
511* sizes int[2] The first element is the base size of the struct as
512* returned by sizeof(struct wcsprm). The second element
513* is the total allocated size, in bytes, assuming that
514* the allocation was done by wcsini(). This figure
515* includes memory allocated for members of constituent
516* structs, such as wcsprm::lin.
517*
518* It is not an error for the struct not to have been set
519* up via wcsset(), which normally results in additional
520* memory allocation.
521*
522* Function return value:
523* int Status return value:
524* 0: Success.
525*
526*
527* auxsize() - Compute the size of a auxprm struct
528* -----------------------------------------------
529* auxsize() computes the full size of an auxprm struct, including allocated
530* memory.
531*
532* Given:
533* aux const struct auxprm*
534* Auxiliary coordinate information.
535*
536* If NULL, the base size of the struct and the allocated
537* size are both set to zero.
538*
539* Returned:
540* sizes int[2] The first element is the base size of the struct as
541* returned by sizeof(struct auxprm). The second element
542* is the total allocated size, in bytes, currently zero.
543*
544* Function return value:
545* int Status return value:
546* 0: Success.
547*
548*
549* wcsenq() - enquire about the state of a wcsprm struct
550* -----------------------------------------------------
551* wcsenq() may be used to obtain information about the state of a wcsprm
552* struct. The function returns a true/false answer for the enquiry asked.
553*
554* Given:
555* wcs const struct wcsprm*
556* Coordinate transformation parameters.
557*
558* enquiry int Enquiry according to the following parameters:
559* WCSENQ_MEM: memory in the struct is being managed by
560* WCSLIB (see wcsini()).
561* WCSENQ_SET: the struct has been set up by wcsset().
562* WCSENQ_BYP: the struct is in bypass mode (see
563* wcsset()).
564* WCSENQ_CHK: the struct is self-consistent in that
565* no changes have been made to any of the
566* "parameters to be given" since the last
567* call to wcsset().
568* These may be combined by logical OR, e.g.
569* WCSENQ_MEM | WCSENQ_SET. The enquiry result will be
570* the logical AND of the individual results.
571*
572* Function return value:
573* int Enquiry result:
574* 0: False.
575* 1: True.
576*
577*
578* wcsprt() - Print routine for the wcsprm struct
579* ----------------------------------------------
580* wcsprt() prints the contents of a wcsprm struct using wcsprintf(). Mainly
581* intended for diagnostic purposes.
582*
583* Given:
584* wcs const struct wcsprm*
585* Coordinate transformation parameters.
586*
587* Function return value:
588* int Status return value:
589* 0: Success.
590* 1: Null wcsprm pointer passed.
591*
592*
593* wcsperr() - Print error messages from a wcsprm struct
594* -----------------------------------------------------
595* wcsperr() prints the error message(s), if any, stored in a wcsprm struct,
596* and the linprm, celprm, prjprm, spcprm, and tabprm structs that it contains.
597* If there are no errors then nothing is printed. It uses wcserr_prt(), q.v.
598*
599* Given:
600* wcs const struct wcsprm*
601* Coordinate transformation parameters.
602*
603* prefix const char *
604* If non-NULL, each output line will be prefixed with
605* this string.
606*
607* Function return value:
608* int Status return value:
609* 0: Success.
610* 1: Null wcsprm pointer passed.
611*
612*
613* wcsbchk() - Enable/disable bounds checking
614* ------------------------------------------
615* wcsbchk() is used to control bounds checking in the projection routines.
616* Note that wcsset() always enables bounds checking. wcsbchk() will invoke
617* wcsset() on the wcsprm struct beforehand if necessary.
618*
619* Given and returned:
620* wcs struct wcsprm*
621* Coordinate transformation parameters.
622*
623* Given:
624* bounds int If bounds&1 then enable strict bounds checking for the
625* spherical-to-Cartesian (s2x) transformation for the
626* AZP, SZP, TAN, SIN, ZPN, and COP projections.
627*
628* If bounds&2 then enable strict bounds checking for the
629* Cartesian-to-spherical (x2s) transformation for the
630* HPX and XPH projections.
631*
632* If bounds&4 then enable bounds checking on the native
633* coordinates returned by the Cartesian-to-spherical
634* (x2s) transformations using prjchk().
635*
636* Zero it to disable all checking.
637*
638* Function return value:
639* int Status return value:
640* 0: Success.
641* 1: Null wcsprm pointer passed.
642*
643*
644* wcsset() - Setup routine for the wcsprm struct
645* ----------------------------------------------
646* wcsset() sets up a wcsprm struct according to information supplied within
647* it (refer to the description of the wcsprm struct).
648*
649* wcsset() recognizes the NCP projection and converts it to the equivalent SIN
650* projection and likewise translates GLS into SFL. It also translates the
651* AIPS spectral types ('FREQ-LSR', 'FELO-HEL', etc.), possibly changing the
652* input header keywords wcsprm::ctype and/or wcsprm::specsys if necessary.
653*
654* Note that this routine need not be called directly; it will be invoked by
655* wcsp2s() and wcss2p() if the wcsprm::flag is anything other than a
656* predefined magic value.
657*
658* wcsset() normally operates regardless of the value of wcsprm::flag; i.e.
659* even if a struct was previously set up it will be reset unconditionally.
660* However, a wcsprm struct may be put into "bypass" mode by invoking wcsset()
661* initially with wcsprm::flag == 1 (rather than 0). wcsset() will return
662* immediately if invoked on a struct in that state. To take a struct out of
663* bypass mode, simply reset wcsprm::flag to zero. See also wcsenq().
664*
665* Given and returned:
666* wcs struct wcsprm*
667* Coordinate transformation parameters.
668*
669* Function return value:
670* int Status return value:
671* 0: Success.
672* 1: Null wcsprm pointer passed.
673* 2: Memory allocation failed.
674* 3: Linear transformation matrix is singular.
675* 4: Inconsistent or unrecognized coordinate axis
676* types.
677* 5: Invalid parameter value.
678* 6: Invalid coordinate transformation parameters.
679* 7: Ill-conditioned coordinate transformation
680* parameters.
681*
682* For returns > 1, a detailed error message is set in
683* wcsprm::err if enabled, see wcserr_enable().
684*
685* Notes:
686* 1: wcsset() always enables strict bounds checking in the projection
687* routines (via a call to prjini()). Use wcsbchk() to modify
688* bounds-checking after wcsset() is invoked.
689*
690*
691* wcsp2s() - Pixel-to-world transformation
692* ----------------------------------------
693* wcsp2s() transforms pixel coordinates to world coordinates.
694*
695* Given and returned:
696* wcs struct wcsprm*
697* Coordinate transformation parameters.
698*
699* Given:
700* ncoord,
701* nelem int The number of coordinates, each of vector length
702* nelem but containing wcs.naxis coordinate elements.
703* Thus nelem must equal or exceed the value of the
704* NAXIS keyword unless ncoord == 1, in which case nelem
705* is not used.
706*
707* pixcrd const double[ncoord][nelem]
708* Array of pixel coordinates.
709*
710* Returned:
711* imgcrd double[ncoord][nelem]
712* Array of intermediate world coordinates. For
713* celestial axes, imgcrd[][wcs.lng] and
714* imgcrd[][wcs.lat] are the projected x-, and
715* y-coordinates in pseudo "degrees". For spectral
716* axes, imgcrd[][wcs.spec] is the intermediate spectral
717* coordinate, in SI units. For time axes,
718* imgcrd[][wcs.time] is the intermediate time
719* coordinate.
720*
721* phi,theta double[ncoord]
722* Longitude and latitude in the native coordinate system
723* of the projection [deg].
724*
725* world double[ncoord][nelem]
726* Array of world coordinates. For celestial axes,
727* world[][wcs.lng] and world[][wcs.lat] are the
728* celestial longitude and latitude [deg]. For spectral
729* axes, world[][wcs.spec] is the spectral coordinate, in
730* SI units. For time axes, world[][wcs.time] is the
731* time coordinate.
732*
733* stat int[ncoord]
734* Status return value for each coordinate:
735* 0: Success.
736* 1+: A bit mask indicating invalid pixel coordinate
737* element(s).
738*
739* Function return value:
740* int Status return value:
741* 0: Success.
742* 1: Null wcsprm pointer passed.
743* 2: Memory allocation failed.
744* 3: Linear transformation matrix is singular.
745* 4: Inconsistent or unrecognized coordinate axis
746* types.
747* 5: Invalid parameter value.
748* 6: Invalid coordinate transformation parameters.
749* 7: Ill-conditioned coordinate transformation
750* parameters.
751* 8: One or more of the pixel coordinates were
752* invalid, as indicated by the stat vector.
753*
754* For returns > 1, a detailed error message is set in
755* wcsprm::err if enabled, see wcserr_enable().
756*
757*
758* wcss2p() - World-to-pixel transformation
759* ----------------------------------------
760* wcss2p() transforms world coordinates to pixel coordinates.
761*
762* Given and returned:
763* wcs struct wcsprm*
764* Coordinate transformation parameters.
765*
766* Given:
767* ncoord,
768* nelem int The number of coordinates, each of vector length nelem
769* but containing wcs.naxis coordinate elements. Thus
770* nelem must equal or exceed the value of the NAXIS
771* keyword unless ncoord == 1, in which case nelem is not
772* used.
773*
774* world const double[ncoord][nelem]
775* Array of world coordinates. For celestial axes,
776* world[][wcs.lng] and world[][wcs.lat] are the
777* celestial longitude and latitude [deg]. For spectral
778* axes, world[][wcs.spec] is the spectral coordinate, in
779* SI units. For time axes, world[][wcs.time] is the
780* time coordinate.
781*
782* Returned:
783* phi,theta double[ncoord]
784* Longitude and latitude in the native coordinate
785* system of the projection [deg].
786*
787* imgcrd double[ncoord][nelem]
788* Array of intermediate world coordinates. For
789* celestial axes, imgcrd[][wcs.lng] and
790* imgcrd[][wcs.lat] are the projected x-, and
791* y-coordinates in pseudo "degrees". For quadcube
792* projections with a CUBEFACE axis the face number is
793* also returned in imgcrd[][wcs.cubeface]. For
794* spectral axes, imgcrd[][wcs.spec] is the intermediate
795* spectral coordinate, in SI units. For time axes,
796* imgcrd[][wcs.time] is the intermediate time
797* coordinate.
798*
799* pixcrd double[ncoord][nelem]
800* Array of pixel coordinates.
801*
802* stat int[ncoord]
803* Status return value for each coordinate:
804* 0: Success.
805* 1+: A bit mask indicating invalid world coordinate
806* element(s).
807*
808* Function return value:
809* int Status return value:
810* 0: Success.
811* 1: Null wcsprm pointer passed.
812* 2: Memory allocation failed.
813* 3: Linear transformation matrix is singular.
814* 4: Inconsistent or unrecognized coordinate axis
815* types.
816* 5: Invalid parameter value.
817* 6: Invalid coordinate transformation parameters.
818* 7: Ill-conditioned coordinate transformation
819* parameters.
820* 9: One or more of the world coordinates were
821* invalid, as indicated by the stat vector.
822*
823* For returns > 1, a detailed error message is set in
824* wcsprm::err if enabled, see wcserr_enable().
825*
826*
827* wcsmix() - Hybrid coordinate transformation
828* -------------------------------------------
829* wcsmix(), given either the celestial longitude or latitude plus an element
830* of the pixel coordinate, solves for the remaining elements by iterating on
831* the unknown celestial coordinate element using wcss2p(). Refer also to the
832* notes below.
833*
834* Given and returned:
835* wcs struct wcsprm*
836* Indices for the celestial coordinates obtained
837* by parsing the wcsprm::ctype[].
838*
839* Given:
840* mixpix int Which element of the pixel coordinate is given.
841*
842* mixcel int Which element of the celestial coordinate is given:
843* 1: Celestial longitude is given in
844* world[wcs.lng], latitude returned in
845* world[wcs.lat].
846* 2: Celestial latitude is given in
847* world[wcs.lat], longitude returned in
848* world[wcs.lng].
849*
850* vspan const double[2]
851* Solution interval for the celestial coordinate [deg].
852* The ordering of the two limits is irrelevant.
853* Longitude ranges may be specified with any convenient
854* normalization, for example [-120,+120] is the same as
855* [240,480], except that the solution will be returned
856* with the same normalization, i.e. lie within the
857* interval specified.
858*
859* vstep const double
860* Step size for solution search [deg]. If zero, a
861* sensible, although perhaps non-optimal default will be
862* used.
863*
864* viter int If a solution is not found then the step size will be
865* halved and the search recommenced. viter controls how
866* many times the step size is halved. The allowed range
867* is 5 - 10.
868*
869* Given and returned:
870* world double[naxis]
871* World coordinate elements. world[wcs.lng] and
872* world[wcs.lat] are the celestial longitude and
873* latitude [deg]. Which is given and which returned
874* depends on the value of mixcel. All other elements
875* are given.
876*
877* Returned:
878* phi,theta double[naxis]
879* Longitude and latitude in the native coordinate
880* system of the projection [deg].
881*
882* imgcrd double[naxis]
883* Image coordinate elements. imgcrd[wcs.lng] and
884* imgcrd[wcs.lat] are the projected x-, and
885* y-coordinates in pseudo "degrees".
886*
887* Given and returned:
888* pixcrd double[naxis]
889* Pixel coordinate. The element indicated by mixpix is
890* given and the remaining elements are returned.
891*
892* Function return value:
893* int Status return value:
894* 0: Success.
895* 1: Null wcsprm pointer passed.
896* 2: Memory allocation failed.
897* 3: Linear transformation matrix is singular.
898* 4: Inconsistent or unrecognized coordinate axis
899* types.
900* 5: Invalid parameter value.
901* 6: Invalid coordinate transformation parameters.
902* 7: Ill-conditioned coordinate transformation
903* parameters.
904* 10: Invalid world coordinate.
905* 11: No solution found in the specified interval.
906*
907* For returns > 1, a detailed error message is set in
908* wcsprm::err if enabled, see wcserr_enable().
909*
910* Notes:
911* 1: Initially the specified solution interval is checked to see if it's a
912* "crossing" interval. If it isn't, a search is made for a crossing
913* solution by iterating on the unknown celestial coordinate starting at
914* the upper limit of the solution interval and decrementing by the
915* specified step size. A crossing is indicated if the trial value of the
916* pixel coordinate steps through the value specified. If a crossing
917* interval is found then the solution is determined by a modified form of
918* "regula falsi" division of the crossing interval. If no crossing
919* interval was found within the specified solution interval then a search
920* is made for a "non-crossing" solution as may arise from a point of
921* tangency. The process is complicated by having to make allowance for
922* the discontinuities that occur in all map projections.
923*
924* Once one solution has been determined others may be found by subsequent
925* invokations of wcsmix() with suitably restricted solution intervals.
926*
927* Note the circumstance that arises when the solution point lies at a
928* native pole of a projection in which the pole is represented as a
929* finite curve, for example the zenithals and conics. In such cases two
930* or more valid solutions may exist but wcsmix() only ever returns one.
931*
932* Because of its generality wcsmix() is very compute-intensive. For
933* compute-limited applications more efficient special-case solvers could
934* be written for simple projections, for example non-oblique cylindrical
935* projections.
936*
937*
938* wcsccs() - Change celestial coordinate system
939* ---------------------------------------------
940* wcsccs() changes the celestial coordinate system of a wcsprm struct. For
941* example, from equatorial to galactic coordinates.
942*
943* Parameters that define the spherical coordinate transformation, essentially
944* being three Euler angles, must be provided. Thereby wcsccs() does not need
945* prior knowledge of specific celestial coordinate systems. It also has the
946* advantage of making it completely general.
947*
948* Auxiliary members of the wcsprm struct relating to equatorial celestial
949* coordinate systems may also be changed.
950*
951* Only orthodox spherical coordinate systems are supported. That is, they
952* must be right-handed, with latitude increasing from zero at the equator to
953* +90 degrees at the pole. This precludes systems such as aziumuth and zenith
954* distance, which, however, could be handled as negative azimuth and
955* elevation.
956*
957* PLEASE NOTE: Information in the wcsprm struct relating to the original
958* coordinate system will be overwritten and therefore lost. If this is
959* undesirable, invoke wcsccs() on a copy of the struct made with wcssub().
960* The wcsprm struct is reset on return with an explicit call to wcsset().
961*
962* Given and returned:
963* wcs struct wcsprm*
964* Coordinate transformation parameters. Particular
965* "values to be given" elements of the wcsprm struct
966* are modified.
967*
968* Given:
969* lng2p1,
970* lat2p1 double Longitude and latitude in the new celestial coordinate
971* system of the pole (i.e. latitude +90) of the original
972* system [deg]. See notes 1 and 2 below.
973*
974* lng1p2 double Longitude in the original celestial coordinate system
975* of the pole (i.e. latitude +90) of the new system
976* [deg]. See note 1 below.
977*
978* clng,clat const char*
979* Longitude and latitude identifiers of the new CTYPEia
980* celestial axis codes, without trailing dashes. For
981* example, "RA" and "DEC" or "GLON" and "GLAT". Up to
982* four characters are used, longer strings need not be
983* null-terminated.
984*
985* radesys const char*
986* Used when transforming to equatorial coordinates,
987* identified by clng == "RA" and clat = "DEC". May be
988* set to the null pointer to preserve the current value.
989* Up to 71 characters are used, longer strings need not
990* be null-terminated.
991*
992* If the new coordinate system is anything other than
993* equatorial, then wcsprm::radesys will be cleared.
994*
995* equinox double Used when transforming to equatorial coordinates. May
996* be set to zero to preserve the current value.
997*
998* If the new coordinate system is not equatorial, then
999* wcsprm::equinox will be marked as undefined.
1000*
1001* alt const char*
1002* Character code for alternate coordinate descriptions
1003* (i.e. the 'a' in keyword names such as CTYPEia). This
1004* is blank for the primary coordinate description, or
1005* one of the 26 upper-case letters, A-Z. May be set to
1006* the null pointer, or null string if no change is
1007* required.
1008*
1009* Function return value:
1010* int Status return value:
1011* 0: Success.
1012* 1: Null wcsprm pointer passed.
1013* 12: Invalid subimage specification (no celestial
1014* axes).
1015*
1016* Notes:
1017* 1: Follows the prescription given in WCS Paper II, Sect. 2.7 for changing
1018* celestial coordinates.
1019*
1020* The implementation takes account of indeterminacies that arise in that
1021* prescription in the particular cases where one of the poles of the new
1022* system is at the fiducial point, or one of them is at the native pole.
1023*
1024* 2: If lat2p1 == +90, i.e. where the poles of the two coordinate systems
1025* coincide, then the spherical coordinate transformation becomes a simple
1026* change in origin of longitude given by
1027* lng2 = lng1 + (lng2p1 - lng1p2 - 180), and lat2 = lat1, where
1028* (lng2,lat2) are coordinates in the new system, and (lng1,lat1) are
1029* coordinates in the original system.
1030*
1031* Likewise, if lat2p1 == -90, then lng2 = -lng1 + (lng2p1 + lng1p2), and
1032* lat2 = -lat1.
1033*
1034* 3: For example, if the original coordinate system is B1950 equatorial and
1035* the desired new coordinate system is galactic, then
1036*
1037* - (lng2p1,lat2p1) are the galactic coordinates of the B1950 celestial
1038* pole, defined by the IAU to be (123.0,+27.4), and lng1p2 is the B1950
1039* right ascension of the galactic pole, defined as 192.25. Clearly
1040* these coordinates are fixed for a particular coordinate
1041* transformation.
1042*
1043* - (clng,clat) would be 'GLON' and 'GLAT', these being the FITS standard
1044* identifiers for galactic coordinates.
1045*
1046* - Since the new coordinate system is not equatorial, wcsprm::radesys
1047* and wcsprm::equinox will be cleared.
1048*
1049* 4. The coordinates required for some common transformations (obtained from
1050* https://ned.ipac.caltech.edu/coordinate_calculator) are as follows:
1051*
1052= (123.0000,+27.4000) galactic coordinates of B1950 celestial pole,
1053= (192.2500,+27.4000) B1950 equatorial coordinates of galactic pole.
1054*
1055= (122.9319,+27.1283) galactic coordinates of J2000 celestial pole,
1056= (192.8595,+27.1283) J2000 equatorial coordinates of galactic pole.
1057*
1058= (359.6774,+89.7217) B1950 equatorial coordinates of J2000 pole,
1059= (180.3162,+89.7217) J2000 equatorial coordinates of B1950 pole.
1060*
1061= (270.0000,+66.5542) B1950 equatorial coordinates of B1950 ecliptic pole,
1062= ( 90.0000,+66.5542) B1950 ecliptic coordinates of B1950 celestial pole.
1063*
1064= (270.0000,+66.5607) J2000 equatorial coordinates of J2000 ecliptic pole,
1065= ( 90.0000,+66.5607) J2000 ecliptic coordinates of J2000 celestial pole.
1066*
1067= ( 26.7315,+15.6441) supergalactic coordinates of B1950 celestial pole,
1068= (283.1894,+15.6441) B1950 equatorial coordinates of supergalactic pole.
1069*
1070= ( 26.4505,+15.7089) supergalactic coordinates of J2000 celestial pole,
1071= (283.7542,+15.7089) J2000 equatorial coordinates of supergalactic pole.
1072*
1073*
1074* wcssptr() - Spectral axis translation
1075* -------------------------------------
1076* wcssptr() translates the spectral axis in a wcsprm struct. For example, a
1077* 'FREQ' axis may be translated into 'ZOPT-F2W' and vice versa.
1078*
1079* PLEASE NOTE: Information in the wcsprm struct relating to the original
1080* coordinate system will be overwritten and therefore lost. If this is
1081* undesirable, invoke wcssptr() on a copy of the struct made with wcssub().
1082* The wcsprm struct is reset on return with an explicit call to wcsset().
1083*
1084* Given and returned:
1085* wcs struct wcsprm*
1086* Coordinate transformation parameters.
1087*
1088* i int* Index of the spectral axis (0-relative). If given < 0
1089* it will be set to the first spectral axis identified
1090* from the ctype[] keyvalues in the wcsprm struct.
1091*
1092* ctype char[9] Desired spectral CTYPEia. Wildcarding may be used as
1093* for the ctypeS2 argument to spctrn() as described in
1094* the prologue of spc.h, i.e. if the final three
1095* characters are specified as "???", or if just the
1096* eighth character is specified as '?', the correct
1097* algorithm code will be substituted and returned.
1098*
1099* Function return value:
1100* int Status return value:
1101* 0: Success.
1102* 1: Null wcsprm pointer passed.
1103* 2: Memory allocation failed.
1104* 3: Linear transformation matrix is singular.
1105* 4: Inconsistent or unrecognized coordinate axis
1106* types.
1107* 5: Invalid parameter value.
1108* 6: Invalid coordinate transformation parameters.
1109* 7: Ill-conditioned coordinate transformation
1110* parameters.
1111* 12: Invalid subimage specification (no spectral
1112* axis).
1113*
1114* For returns > 1, a detailed error message is set in
1115* wcsprm::err if enabled, see wcserr_enable().
1116*
1117*
1118* wcslib_version() - WCSLIB version number
1119* ----------------------------------------
1120* wcslib_version() returns the WCSLIB version number.
1121*
1122* The major version number changes when the ABI changes or when the license
1123* conditions change. ABI changes typically result from a change to the
1124* contents of one of the structs. The major version number is used to
1125* distinguish between incompatible versions of the sharable library.
1126*
1127* The minor version number changes with new functionality or bug fixes that do
1128* not involve a change in the ABI.
1129*
1130* The auxiliary version number (which is often absent) signals changes to the
1131* documentation, test suite, build procedures, or any other change that does
1132* not affect the compiled library.
1133*
1134* Returned:
1135* vers[3] int[3] The broken-down version number:
1136* 0: Major version number.
1137* 1: Minor version number.
1138* 2: Auxiliary version number (zero if absent).
1139* May be given as a null pointer if not required.
1140*
1141* Function return value:
1142* char* A null-terminated, statically allocated string
1143* containing the version number in the usual form, i.e.
1144* "<major>.<minor>.<auxiliary>".
1145*
1146*
1147* wcsprm struct - Coordinate transformation parameters
1148* ----------------------------------------------------
1149* The wcsprm struct contains information required to transform world
1150* coordinates. It consists of certain members that must be set by the user
1151* ("given") and others that are set by the WCSLIB routines ("returned").
1152* While the addresses of the arrays themselves may be set by wcsinit() if it
1153* (optionally) allocates memory, their contents must be set by the user.
1154*
1155* Some parameters that are given are not actually required for transforming
1156* coordinates. These are described as "auxiliary"; the struct simply provides
1157* a place to store them, though they may be used by wcshdo() in constructing a
1158* FITS header from a wcsprm struct. Some of the returned values are supplied
1159* for informational purposes and others are for internal use only as
1160* indicated.
1161*
1162* In practice, it is expected that a WCS parser would scan the FITS header to
1163* determine the number of coordinate axes. It would then use wcsinit() to
1164* allocate memory for arrays in the wcsprm struct and set default values.
1165* Then as it reread the header and identified each WCS keyrecord it would load
1166* the value into the relevant wcsprm array element. This is essentially what
1167* wcspih() does - refer to the prologue of wcshdr.h. As the final step,
1168* wcsset() is invoked, either directly or indirectly, to set the derived
1169* members of the wcsprm struct. wcsset() strips off trailing blanks in all
1170* string members and null-fills the character array.
1171*
1172* int flag
1173* (Given and returned) This flag must be set to zero (or 1, see wcsset())
1174* whenever any of the following wcsprm members are set or changed:
1175*
1176* - wcsprm::naxis (q.v., not normally set by the user),
1177* - wcsprm::crpix,
1178* - wcsprm::pc,
1179* - wcsprm::cdelt,
1180* - wcsprm::crval,
1181* - wcsprm::cunit,
1182* - wcsprm::ctype,
1183* - wcsprm::lonpole,
1184* - wcsprm::latpole,
1185* - wcsprm::restfrq,
1186* - wcsprm::restwav,
1187* - wcsprm::npv,
1188* - wcsprm::pv,
1189* - wcsprm::nps,
1190* - wcsprm::ps,
1191* - wcsprm::cd,
1192* - wcsprm::crota,
1193* - wcsprm::altlin,
1194* - wcsprm::ntab,
1195* - wcsprm::nwtb,
1196* - wcsprm::tab,
1197* - wcsprm::wtb.
1198*
1199* This signals the initialization routine, wcsset(), to recompute the
1200* returned members of the linprm, celprm, spcprm, and tabprm structs.
1201* wcsset() will reset flag to indicate that this has been done.
1202*
1203* PLEASE NOTE: flag should be set to -1 when wcsinit() is called for the
1204* first time for a particular wcsprm struct in order to initialize memory
1205* management. It must ONLY be used on the first initialization otherwise
1206* memory leaks may result.
1207*
1208* int naxis
1209* (Given or returned) Number of pixel and world coordinate elements.
1210*
1211* If wcsinit() is used to initialize the linprm struct (as would normally
1212* be the case) then it will set naxis from the value passed to it as a
1213* function argument. The user should not subsequently modify it.
1214*
1215* double *crpix
1216* (Given) Address of the first element of an array of double containing
1217* the coordinate reference pixel, CRPIXja.
1218*
1219* double *pc
1220* (Given) Address of the first element of the PCi_ja (pixel coordinate)
1221* transformation matrix. The expected order is
1222*
1223= struct wcsprm wcs;
1224= wcs.pc = {PC1_1, PC1_2, PC2_1, PC2_2};
1225*
1226* This may be constructed conveniently from a 2-D array via
1227*
1228= double m[2][2] = {{PC1_1, PC1_2},
1229= {PC2_1, PC2_2}};
1230*
1231* which is equivalent to
1232*
1233= double m[2][2];
1234= m[0][0] = PC1_1;
1235= m[0][1] = PC1_2;
1236= m[1][0] = PC2_1;
1237= m[1][1] = PC2_2;
1238*
1239* The storage order for this 2-D array is the same as for the 1-D array,
1240* whence
1241*
1242= wcs.pc = *m;
1243*
1244* would be legitimate.
1245*
1246* double *cdelt
1247* (Given) Address of the first element of an array of double containing
1248* the coordinate increments, CDELTia.
1249*
1250* double *crval
1251* (Given) Address of the first element of an array of double containing
1252* the coordinate reference values, CRVALia.
1253*
1254* char (*cunit)[72]
1255* (Given) Address of the first element of an array of char[72] containing
1256* the CUNITia keyvalues which define the units of measurement of the
1257* CRVALia, CDELTia, and CDi_ja keywords.
1258*
1259* As CUNITia is an optional header keyword, cunit[][72] may be left blank
1260* but otherwise is expected to contain a standard units specification as
1261* defined by WCS Paper I. Utility function wcsutrn(), described in
1262* wcsunits.h, is available to translate commonly used non-standard units
1263* specifications but this must be done as a separate step before invoking
1264* wcsset().
1265*
1266* For celestial axes, if cunit[][72] is not blank, wcsset() uses
1267* wcsunits() to parse it and scale cdelt[], crval[], and cd[][*] to
1268* degrees. It then resets cunit[][72] to "deg".
1269*
1270* For spectral axes, if cunit[][72] is not blank, wcsset() uses wcsunits()
1271* to parse it and scale cdelt[], crval[], and cd[][*] to SI units. It
1272* then resets cunit[][72] accordingly.
1273*
1274* wcsset() ignores cunit[][72] for other coordinate types; cunit[][72] may
1275* be used to label coordinate values.
1276*
1277* These variables accomodate the longest allowed string-valued FITS
1278* keyword, being limited to 68 characters, plus the null-terminating
1279* character.
1280*
1281* char (*ctype)[72]
1282* (Given) Address of the first element of an array of char[72] containing
1283* the coordinate axis types, CTYPEia.
1284*
1285* The ctype[][72] keyword values must be in upper case and there must be
1286* zero or one pair of matched celestial axis types, and zero or one
1287* spectral axis. The ctype[][72] strings should be padded with blanks on
1288* the right and null-terminated so that they are at least eight characters
1289* in length.
1290*
1291* These variables accomodate the longest allowed string-valued FITS
1292* keyword, being limited to 68 characters, plus the null-terminating
1293* character.
1294*
1295* double lonpole
1296* (Given and returned) The native longitude of the celestial pole, phi_p,
1297* given by LONPOLEa [deg] or by PVi_2a [deg] attached to the longitude
1298* axis which takes precedence if defined, and ...
1299* double latpole
1300* (Given and returned) ... the native latitude of the celestial pole,
1301* theta_p, given by LATPOLEa [deg] or by PVi_3a [deg] attached to the
1302* longitude axis which takes precedence if defined.
1303*
1304* lonpole and latpole may be left to default to values set by wcsinit()
1305* (see celprm::ref), but in any case they will be reset by wcsset() to
1306* the values actually used. Note therefore that if the wcsprm struct is
1307* reused without resetting them, whether directly or via wcsinit(), they
1308* will no longer have their default values.
1309*
1310* double restfrq
1311* (Given) The rest frequency [Hz], and/or ...
1312* double restwav
1313* (Given) ... the rest wavelength in vacuo [m], only one of which need be
1314* given, the other should be set to zero.
1315*
1316* int npv
1317* (Given) The number of entries in the wcsprm::pv[] array.
1318*
1319* int npvmax
1320* (Given or returned) The length of the wcsprm::pv[] array.
1321*
1322* npvmax will be set by wcsinit() if it allocates memory for wcsprm::pv[],
1323* otherwise it must be set by the user. See also wcsnpv().
1324*
1325* struct pvcard *pv
1326* (Given) Address of the first element of an array of length npvmax of
1327* pvcard structs.
1328*
1329* As a FITS header parser encounters each PVi_ma keyword it should load it
1330* into a pvcard struct in the array and increment npv. wcsset()
1331* interprets these as required.
1332*
1333* Note that, if they were not given, wcsset() resets the entries for
1334* PVi_1a, PVi_2a, PVi_3a, and PVi_4a for longitude axis i to match
1335* phi_0 and theta_0 (the native longitude and latitude of the reference
1336* point), LONPOLEa and LATPOLEa respectively.
1337*
1338* int nps
1339* (Given) The number of entries in the wcsprm::ps[] array.
1340*
1341* int npsmax
1342* (Given or returned) The length of the wcsprm::ps[] array.
1343*
1344* npsmax will be set by wcsinit() if it allocates memory for wcsprm::ps[],
1345* otherwise it must be set by the user. See also wcsnps().
1346*
1347* struct pscard *ps
1348* (Given) Address of the first element of an array of length npsmax of
1349* pscard structs.
1350*
1351* As a FITS header parser encounters each PSi_ma keyword it should load it
1352* into a pscard struct in the array and increment nps. wcsset()
1353* interprets these as required (currently no PSi_ma keyvalues are
1354* recognized).
1355*
1356* double *cd
1357* (Given) For historical compatibility, the wcsprm struct supports two
1358* alternate specifications of the linear transformation matrix, those
1359* associated with the CDi_ja keywords, and ...
1360* double *crota
1361* (Given) ... those associated with the CROTAi keywords. Although these
1362* may not formally co-exist with PCi_ja, the approach taken here is simply
1363* to ignore them if given in conjunction with PCi_ja.
1364*
1365* int altlin
1366* (Given) altlin is a bit flag that denotes which of the PCi_ja, CDi_ja
1367* and CROTAi keywords are present in the header:
1368*
1369* - Bit 0: PCi_ja is present.
1370*
1371* - Bit 1: CDi_ja is present.
1372*
1373* Matrix elements in the IRAF convention are equivalent to the product
1374* CDi_ja = CDELTia * PCi_ja, but the defaults differ from that of the
1375* PCi_ja matrix. If one or more CDi_ja keywords are present then all
1376* unspecified CDi_ja default to zero. If no CDi_ja (or CROTAi) keywords
1377* are present, then the header is assumed to be in PCi_ja form whether
1378* or not any PCi_ja keywords are present since this results in an
1379* interpretation of CDELTia consistent with the original FITS
1380* specification.
1381*
1382* While CDi_ja may not formally co-exist with PCi_ja, it may co-exist
1383* with CDELTia and CROTAi which are to be ignored.
1384*
1385* - Bit 2: CROTAi is present.
1386*
1387* In the AIPS convention, CROTAi may only be associated with the
1388* latitude axis of a celestial axis pair. It specifies a rotation in
1389* the image plane that is applied AFTER the CDELTia; any other CROTAi
1390* keywords are ignored.
1391*
1392* CROTAi may not formally co-exist with PCi_ja.
1393*
1394* CROTAi and CDELTia may formally co-exist with CDi_ja but if so are to
1395* be ignored.
1396*
1397* - Bit 3: PCi_ja + CDELTia was derived from CDi_ja by wcspcx().
1398*
1399* This bit is set by wcspcx() when it derives PCi_ja and CDELTia from
1400* CDi_ja via an orthonormal decomposition. In particular, it signals
1401* wcsset() not to replace PCi_ja by a copy of CDi_ja with CDELTia set
1402* to unity.
1403*
1404* CDi_ja and CROTAi keywords, if found, are to be stored in the wcsprm::cd
1405* and wcsprm::crota arrays which are dimensioned similarly to wcsprm::pc
1406* and wcsprm::cdelt. FITS header parsers should use the following
1407* procedure:
1408*
1409* - Whenever a PCi_ja keyword is encountered: altlin |= 1;
1410*
1411* - Whenever a CDi_ja keyword is encountered: altlin |= 2;
1412*
1413* - Whenever a CROTAi keyword is encountered: altlin |= 4;
1414*
1415* If none of these bits are set the PCi_ja representation results, i.e.
1416* wcsprm::pc and wcsprm::cdelt will be used as given.
1417*
1418* These alternate specifications of the linear transformation matrix are
1419* translated immediately to PCi_ja by wcsset() and are invisible to the
1420* lower-level WCSLIB routines. In particular, unless bit 3 is also set,
1421* wcsset() resets wcsprm::cdelt to unity if CDi_ja is present (and no
1422* PCi_ja).
1423*
1424* If CROTAi are present but none is associated with the latitude axis
1425* (and no PCi_ja or CDi_ja), then wcsset() reverts to a unity PCi_ja
1426* matrix.
1427*
1428* int velref
1429* (Given) AIPS velocity code VELREF, refer to spcaips().
1430*
1431* It is not necessary to reset the wcsprm struct (via wcsset()) when
1432* wcsprm::velref is changed.
1433*
1434* char alt[4]
1435* (Given, auxiliary) Character code for alternate coordinate descriptions
1436* (i.e. the 'a' in keyword names such as CTYPEia). This is blank for the
1437* primary coordinate description, or one of the 26 upper-case letters,
1438* A-Z.
1439*
1440* An array of four characters is provided for alignment purposes, only the
1441* first is used.
1442*
1443* It is not necessary to reset the wcsprm struct (via wcsset()) when
1444* wcsprm::alt is changed.
1445*
1446* int colnum
1447* (Given, auxiliary) Where the coordinate representation is associated
1448* with an image-array column in a FITS binary table, this variable may be
1449* used to record the relevant column number.
1450*
1451* It should be set to zero for an image header or pixel list.
1452*
1453* It is not necessary to reset the wcsprm struct (via wcsset()) when
1454* wcsprm::colnum is changed.
1455*
1456* int *colax
1457* (Given, auxiliary) Address of the first element of an array of int
1458* recording the column numbers for each axis in a pixel list.
1459*
1460* The array elements should be set to zero for an image header or image
1461* array in a binary table.
1462*
1463* It is not necessary to reset the wcsprm struct (via wcsset()) when
1464* wcsprm::colax is changed.
1465*
1466* char (*cname)[72]
1467* (Given, auxiliary) The address of the first element of an array of
1468* char[72] containing the coordinate axis names, CNAMEia.
1469*
1470* These variables accomodate the longest allowed string-valued FITS
1471* keyword, being limited to 68 characters, plus the null-terminating
1472* character.
1473*
1474* It is not necessary to reset the wcsprm struct (via wcsset()) when
1475* wcsprm::cname is changed.
1476*
1477* double *crder
1478* (Given, auxiliary) Address of the first element of an array of double
1479* recording the random error in the coordinate value, CRDERia.
1480*
1481* It is not necessary to reset the wcsprm struct (via wcsset()) when
1482* wcsprm::crder is changed.
1483*
1484* double *csyer
1485* (Given, auxiliary) Address of the first element of an array of double
1486* recording the systematic error in the coordinate value, CSYERia.
1487*
1488* It is not necessary to reset the wcsprm struct (via wcsset()) when
1489* wcsprm::csyer is changed.
1490*
1491* double *czphs
1492* (Given, auxiliary) Address of the first element of an array of double
1493* recording the time at the zero point of a phase axis, CZPHSia.
1494*
1495* It is not necessary to reset the wcsprm struct (via wcsset()) when
1496* wcsprm::czphs is changed.
1497*
1498* double *cperi
1499* (Given, auxiliary) Address of the first element of an array of double
1500* recording the period of a phase axis, CPERIia.
1501*
1502* It is not necessary to reset the wcsprm struct (via wcsset()) when
1503* wcsprm::cperi is changed.
1504*
1505* char wcsname[72]
1506* (Given, auxiliary) The name given to the coordinate representation,
1507* WCSNAMEa. This variable accomodates the longest allowed string-valued
1508* FITS keyword, being limited to 68 characters, plus the null-terminating
1509* character.
1510*
1511* It is not necessary to reset the wcsprm struct (via wcsset()) when
1512* wcsprm::wcsname is changed.
1513*
1514* char timesys[72]
1515* (Given, auxiliary) TIMESYS keyvalue, being the time scale (UTC, TAI,
1516* etc.) in which all other time-related auxiliary header values are
1517* recorded. Also defines the time scale for an image axis with CTYPEia
1518* set to 'TIME'.
1519*
1520* It is not necessary to reset the wcsprm struct (via wcsset()) when
1521* wcsprm::timesys is changed.
1522*
1523* char trefpos[72]
1524* (Given, auxiliary) TREFPOS keyvalue, being the location in space where
1525* the recorded time is valid.
1526*
1527* It is not necessary to reset the wcsprm struct (via wcsset()) when
1528* wcsprm::trefpos is changed.
1529*
1530* char trefdir[72]
1531* (Given, auxiliary) TREFDIR keyvalue, being the reference direction used
1532* in calculating a pathlength delay.
1533*
1534* It is not necessary to reset the wcsprm struct (via wcsset()) when
1535* wcsprm::trefdir is changed.
1536*
1537* char plephem[72]
1538* (Given, auxiliary) PLEPHEM keyvalue, being the Solar System ephemeris
1539* used for calculating a pathlength delay.
1540*
1541* It is not necessary to reset the wcsprm struct (via wcsset()) when
1542* wcsprm::plephem is changed.
1543*
1544* char timeunit[72]
1545* (Given, auxiliary) TIMEUNIT keyvalue, being the time units in which
1546* the following header values are expressed: TSTART, TSTOP, TIMEOFFS,
1547* TIMSYER, TIMRDER, TIMEDEL. It also provides the default value for
1548* CUNITia for time axes.
1549*
1550* It is not necessary to reset the wcsprm struct (via wcsset()) when
1551* wcsprm::timeunit is changed.
1552*
1553* char dateref[72]
1554* (Given, auxiliary) DATEREF keyvalue, being the date of a reference epoch
1555* relative to which other time measurements refer.
1556*
1557* It is not necessary to reset the wcsprm struct (via wcsset()) when
1558* wcsprm::dateref is changed.
1559*
1560* double mjdref[2]
1561* (Given, auxiliary) MJDREF keyvalue, equivalent to DATEREF expressed as
1562* a Modified Julian Date (MJD = JD - 2400000.5). The value is given as
1563* the sum of the two-element vector, allowing increased precision.
1564*
1565* It is not necessary to reset the wcsprm struct (via wcsset()) when
1566* wcsprm::mjdref is changed.
1567*
1568* double timeoffs
1569* (Given, auxiliary) TIMEOFFS keyvalue, being a time offset, which may be
1570* used, for example, to provide a uniform clock correction for times
1571* referenced to DATEREF.
1572*
1573* It is not necessary to reset the wcsprm struct (via wcsset()) when
1574* wcsprm::timeoffs is changed.
1575*
1576* char dateobs[72]
1577* (Given, auxiliary) DATE-OBS keyvalue, being the date at the start of the
1578* observation unless otherwise explained in the DATE-OBS keycomment, in
1579* ISO format, yyyy-mm-ddThh:mm:ss.
1580*
1581* It is not necessary to reset the wcsprm struct (via wcsset()) when
1582* wcsprm::dateobs is changed.
1583*
1584* char datebeg[72]
1585* (Given, auxiliary) DATE-BEG keyvalue, being the date at the start of the
1586* observation in ISO format, yyyy-mm-ddThh:mm:ss.
1587*
1588* It is not necessary to reset the wcsprm struct (via wcsset()) when
1589* wcsprm::datebeg is changed.
1590*
1591* char dateavg[72]
1592* (Given, auxiliary) DATE-AVG keyvalue, being the date at a representative
1593* mid-point of the observation in ISO format, yyyy-mm-ddThh:mm:ss.
1594*
1595* It is not necessary to reset the wcsprm struct (via wcsset()) when
1596* wcsprm::dateavg is changed.
1597*
1598* char dateend[72]
1599* (Given, auxiliary) DATE-END keyvalue, baing the date at the end of the
1600* observation in ISO format, yyyy-mm-ddThh:mm:ss.
1601*
1602* It is not necessary to reset the wcsprm struct (via wcsset()) when
1603* wcsprm::dateend is changed.
1604*
1605* double mjdobs
1606* (Given, auxiliary) MJD-OBS keyvalue, equivalent to DATE-OBS expressed
1607* as a Modified Julian Date (MJD = JD - 2400000.5).
1608*
1609* It is not necessary to reset the wcsprm struct (via wcsset()) when
1610* wcsprm::mjdobs is changed.
1611*
1612* double mjdbeg
1613* (Given, auxiliary) MJD-BEG keyvalue, equivalent to DATE-BEG expressed
1614* as a Modified Julian Date (MJD = JD - 2400000.5).
1615*
1616* It is not necessary to reset the wcsprm struct (via wcsset()) when
1617* wcsprm::mjdbeg is changed.
1618*
1619* double mjdavg
1620* (Given, auxiliary) MJD-AVG keyvalue, equivalent to DATE-AVG expressed
1621* as a Modified Julian Date (MJD = JD - 2400000.5).
1622*
1623* It is not necessary to reset the wcsprm struct (via wcsset()) when
1624* wcsprm::mjdavg is changed.
1625*
1626* double mjdend
1627* (Given, auxiliary) MJD-END keyvalue, equivalent to DATE-END expressed
1628* as a Modified Julian Date (MJD = JD - 2400000.5).
1629*
1630* It is not necessary to reset the wcsprm struct (via wcsset()) when
1631* wcsprm::mjdend is changed.
1632*
1633* double jepoch
1634* (Given, auxiliary) JEPOCH keyvalue, equivalent to DATE-OBS expressed
1635* as a Julian epoch.
1636*
1637* It is not necessary to reset the wcsprm struct (via wcsset()) when
1638* wcsprm::jepoch is changed.
1639*
1640* double bepoch
1641* (Given, auxiliary) BEPOCH keyvalue, equivalent to DATE-OBS expressed
1642* as a Besselian epoch
1643*
1644* It is not necessary to reset the wcsprm struct (via wcsset()) when
1645* wcsprm::bepoch is changed.
1646*
1647* double tstart
1648* (Given, auxiliary) TSTART keyvalue, equivalent to DATE-BEG expressed
1649* as a time in units of TIMEUNIT relative to DATEREF+TIMEOFFS.
1650*
1651* It is not necessary to reset the wcsprm struct (via wcsset()) when
1652* wcsprm::tstart is changed.
1653*
1654* double tstop
1655* (Given, auxiliary) TSTOP keyvalue, equivalent to DATE-END expressed
1656* as a time in units of TIMEUNIT relative to DATEREF+TIMEOFFS.
1657*
1658* It is not necessary to reset the wcsprm struct (via wcsset()) when
1659* wcsprm::tstop is changed.
1660*
1661* double xposure
1662* (Given, auxiliary) XPOSURE keyvalue, being the effective exposure time
1663* in units of TIMEUNIT.
1664*
1665* It is not necessary to reset the wcsprm struct (via wcsset()) when
1666* wcsprm::xposure is changed.
1667*
1668* double telapse
1669* (Given, auxiliary) TELAPSE keyvalue, equivalent to the elapsed time
1670* between DATE-BEG and DATE-END, in units of TIMEUNIT.
1671*
1672* It is not necessary to reset the wcsprm struct (via wcsset()) when
1673* wcsprm::telapse is changed.
1674*
1675* double timsyer
1676* (Given, auxiliary) TIMSYER keyvalue, being the absolute error of the
1677* time values, in units of TIMEUNIT.
1678*
1679* It is not necessary to reset the wcsprm struct (via wcsset()) when
1680* wcsprm::timsyer is changed.
1681*
1682* double timrder
1683* (Given, auxiliary) TIMRDER keyvalue, being the accuracy of time stamps
1684* relative to each other, in units of TIMEUNIT.
1685*
1686* It is not necessary to reset the wcsprm struct (via wcsset()) when
1687* wcsprm::timrder is changed.
1688*
1689* double timedel
1690* (Given, auxiliary) TIMEDEL keyvalue, being the resolution of the time
1691* stamps.
1692*
1693* It is not necessary to reset the wcsprm struct (via wcsset()) when
1694* wcsprm::timedel is changed.
1695*
1696* double timepixr
1697* (Given, auxiliary) TIMEPIXR keyvalue, being the relative position of the
1698* time stamps in binned time intervals, a value between 0.0 and 1.0.
1699*
1700* It is not necessary to reset the wcsprm struct (via wcsset()) when
1701* wcsprm::timepixr is changed.
1702*
1703* double obsgeo[6]
1704* (Given, auxiliary) Location of the observer in a standard terrestrial
1705* reference frame. The first three give ITRS Cartesian coordinates
1706* OBSGEO-X [m], OBSGEO-Y [m], OBSGEO-Z [m], and the second three give
1707* OBSGEO-L [deg], OBSGEO-B [deg], OBSGEO-H [m], which are related through
1708* a standard transformation.
1709*
1710* It is not necessary to reset the wcsprm struct (via wcsset()) when
1711* wcsprm::obsgeo is changed.
1712*
1713* char obsorbit[72]
1714* (Given, auxiliary) OBSORBIT keyvalue, being the URI, URL, or name of an
1715* orbit ephemeris file giving spacecraft coordinates relating to TREFPOS.
1716*
1717* It is not necessary to reset the wcsprm struct (via wcsset()) when
1718* wcsprm::obsorbit is changed.
1719*
1720* char radesys[72]
1721* (Given, auxiliary) The equatorial or ecliptic coordinate system type,
1722* RADESYSa.
1723*
1724* It is not necessary to reset the wcsprm struct (via wcsset()) when
1725* wcsprm::radesys is changed.
1726*
1727* double equinox
1728* (Given, auxiliary) The equinox associated with dynamical equatorial or
1729* ecliptic coordinate systems, EQUINOXa (or EPOCH in older headers). Not
1730* applicable to ICRS equatorial or ecliptic coordinates.
1731*
1732* It is not necessary to reset the wcsprm struct (via wcsset()) when
1733* wcsprm::equinox is changed.
1734*
1735* char specsys[72]
1736* (Given, auxiliary) Spectral reference frame (standard of rest),
1737* SPECSYSa.
1738*
1739* It is not necessary to reset the wcsprm struct (via wcsset()) when
1740* wcsprm::specsys is changed.
1741*
1742* char ssysobs[72]
1743* (Given, auxiliary) The spectral reference frame in which there is no
1744* differential variation in the spectral coordinate across the
1745* field-of-view, SSYSOBSa.
1746*
1747* It is not necessary to reset the wcsprm struct (via wcsset()) when
1748* wcsprm::ssysobs is changed.
1749*
1750* double velosys
1751* (Given, auxiliary) The relative radial velocity [m/s] between the
1752* observer and the selected standard of rest in the direction of the
1753* celestial reference coordinate, VELOSYSa.
1754*
1755* It is not necessary to reset the wcsprm struct (via wcsset()) when
1756* wcsprm::velosys is changed.
1757*
1758* double zsource
1759* (Given, auxiliary) The redshift, ZSOURCEa, of the source.
1760*
1761* It is not necessary to reset the wcsprm struct (via wcsset()) when
1762* wcsprm::zsource is changed.
1763*
1764* char ssyssrc[72]
1765* (Given, auxiliary) The spectral reference frame (standard of rest),
1766* SSYSSRCa, in which wcsprm::zsource was measured.
1767*
1768* It is not necessary to reset the wcsprm struct (via wcsset()) when
1769* wcsprm::ssyssrc is changed.
1770*
1771* double velangl
1772* (Given, auxiliary) The angle [deg] that should be used to decompose an
1773* observed velocity into radial and transverse components.
1774*
1775* It is not necessary to reset the wcsprm struct (via wcsset()) when
1776* wcsprm::velangl is changed.
1777*
1778* struct auxprm *aux
1779* (Given, auxiliary) This struct holds auxiliary coordinate system
1780* information of a specialist nature. While these parameters may be
1781* widely recognized within particular fields of astronomy, they differ
1782* from the above auxiliary parameters in not being defined by any of the
1783* FITS WCS standards. Collecting them together in a separate struct that
1784* is allocated only when required helps to control bloat in the size of
1785* the wcsprm struct.
1786*
1787* int ntab
1788* (Given) See wcsprm::tab.
1789*
1790* int nwtb
1791* (Given) See wcsprm::wtb.
1792*
1793* struct tabprm *tab
1794* (Given) Address of the first element of an array of ntab tabprm structs
1795* for which memory has been allocated. These are used to store tabular
1796* transformation parameters.
1797*
1798* Although technically wcsprm::ntab and tab are "given", they will
1799* normally be set by invoking wcstab(), whether directly or indirectly.
1800*
1801* The tabprm structs contain some members that must be supplied and others
1802* that are derived. The information to be supplied comes primarily from
1803* arrays stored in one or more FITS binary table extensions. These
1804* arrays, referred to here as "wcstab arrays", are themselves located by
1805* parameters stored in the FITS image header.
1806*
1807* struct wtbarr *wtb
1808* (Given) Address of the first element of an array of nwtb wtbarr structs
1809* for which memory has been allocated. These are used in extracting
1810* wcstab arrays from a FITS binary table.
1811*
1812* Although technically wcsprm::nwtb and wtb are "given", they will
1813* normally be set by invoking wcstab(), whether directly or indirectly.
1814*
1815* char lngtyp[8]
1816* (Returned) Four-character WCS celestial longitude and ...
1817* char lattyp[8]
1818* (Returned) ... latitude axis types. e.g. "RA", "DEC", "GLON", "GLAT",
1819* etc. extracted from 'RA--', 'DEC-', 'GLON', 'GLAT', etc. in the first
1820* four characters of CTYPEia but with trailing dashes removed. (Declared
1821* as char[8] for alignment reasons.)
1822*
1823* int lng
1824* (Returned) Index for the longitude coordinate, and ...
1825* int lat
1826* (Returned) ... index for the latitude coordinate, and ...
1827* int spec
1828* (Returned) ... index for the spectral coordinate, and ...
1829* int time
1830* (Returned) ... index for the time coordinate in the imgcrd[][] and
1831* world[][] arrays in the API of wcsp2s(), wcss2p() and wcsmix().
1832*
1833* These may also serve as indices into the pixcrd[][] array provided that
1834* the PCi_ja matrix does not transpose axes.
1835*
1836* int cubeface
1837* (Returned) Index into the pixcrd[][] array for the CUBEFACE axis. This
1838* is used for quadcube projections where the cube faces are stored on a
1839* separate axis (see wcs.h).
1840*
1841* int chksum
1842* (Returned) Checksum of keyvalues provided (see wcsprm::flag). Used by
1843* wcsenq() to validate the self-consistency of the struct. Note that
1844* the checksum incorporates addresses and is therefore highly specific to
1845* the instance of the wcsprm struct.
1846*
1847* int *types
1848* (Returned) Address of the first element of an array of int containing a
1849* four-digit type code for each axis.
1850*
1851* - First digit (i.e. 1000s):
1852* - 0: Non-specific coordinate type.
1853* - 1: Stokes coordinate.
1854* - 2: Celestial coordinate (including CUBEFACE).
1855* - 3: Spectral coordinate.
1856* - 4: Time coordinate.
1857*
1858* - Second digit (i.e. 100s):
1859* - 0: Linear axis.
1860* - 1: Quantized axis (STOKES, CUBEFACE).
1861* - 2: Non-linear celestial axis.
1862* - 3: Non-linear spectral axis.
1863* - 4: Logarithmic axis.
1864* - 5: Tabular axis.
1865*
1866* - Third digit (i.e. 10s):
1867* - 0: Group number, e.g. lookup table number, being an index into the
1868* tabprm array (see above).
1869*
1870* - The fourth digit is used as a qualifier depending on the axis type.
1871*
1872* - For celestial axes:
1873* - 0: Longitude coordinate.
1874* - 1: Latitude coordinate.
1875* - 2: CUBEFACE number.
1876*
1877* - For lookup tables: the axis number in a multidimensional table.
1878*
1879* CTYPEia in "4-3" form with unrecognized algorithm code will have its
1880* type set to -1 and generate an error.
1881*
1882* struct linprm lin
1883* (Returned) Linear transformation parameters (usage is described in the
1884* prologue to lin.h).
1885*
1886* struct celprm cel
1887* (Returned) Celestial transformation parameters (usage is described in
1888* the prologue to cel.h).
1889*
1890* struct spcprm spc
1891* (Returned) Spectral transformation parameters (usage is described in the
1892* prologue to spc.h).
1893*
1894* struct wcserr *err
1895* (Returned) If enabled, when an error status is returned, this struct
1896* contains detailed information about the error, see wcserr_enable().
1897*
1898* int m_flag
1899* (For internal use only.)
1900* int m_naxis
1901* (For internal use only.)
1902* double *m_crpix
1903* (For internal use only.)
1904* double *m_pc
1905* (For internal use only.)
1906* double *m_cdelt
1907* (For internal use only.)
1908* double *m_crval
1909* (For internal use only.)
1910* char (*m_cunit)[72]
1911* (For internal use only.)
1912* char (*m_ctype)[72]
1913* (For internal use only.)
1914* struct pvcard *m_pv
1915* (For internal use only.)
1916* struct pscard *m_ps
1917* (For internal use only.)
1918* double *m_cd
1919* (For internal use only.)
1920* double *m_crota
1921* (For internal use only.)
1922* int *m_colax
1923* (For internal use only.)
1924* char (*m_cname)[72]
1925* (For internal use only.)
1926* double *m_crder
1927* (For internal use only.)
1928* double *m_csyer
1929* (For internal use only.)
1930* double *m_czphs
1931* (For internal use only.)
1932* double *m_cperi
1933* (For internal use only.)
1934* struct tabprm *m_tab
1935* (For internal use only.)
1936* struct wtbarr *m_wtb
1937* (For internal use only.)
1938*
1939*
1940* pvcard struct - Store for PVi_ma keyrecords
1941* -------------------------------------------
1942* The pvcard struct is used to pass the parsed contents of PVi_ma keyrecords
1943* to wcsset() via the wcsprm struct.
1944*
1945* All members of this struct are to be set by the user.
1946*
1947* int i
1948* (Given) Axis number (1-relative), as in the FITS PVi_ma keyword. If
1949* i == 0, wcsset() will replace it with the latitude axis number.
1950*
1951* int m
1952* (Given) Parameter number (non-negative), as in the FITS PVi_ma keyword.
1953*
1954* double value
1955* (Given) Parameter value.
1956*
1957*
1958* pscard struct - Store for PSi_ma keyrecords
1959* -------------------------------------------
1960* The pscard struct is used to pass the parsed contents of PSi_ma keyrecords
1961* to wcsset() via the wcsprm struct.
1962*
1963* All members of this struct are to be set by the user.
1964*
1965* int i
1966* (Given) Axis number (1-relative), as in the FITS PSi_ma keyword.
1967*
1968* int m
1969* (Given) Parameter number (non-negative), as in the FITS PSi_ma keyword.
1970*
1971* char value[72]
1972* (Given) Parameter value.
1973*
1974*
1975* auxprm struct - Additional auxiliary parameters
1976* -----------------------------------------------
1977* The auxprm struct holds auxiliary coordinate system information of a
1978* specialist nature. It is anticipated that this struct will expand in future
1979* to accomodate additional parameters.
1980*
1981* All members of this struct are to be set by the user.
1982*
1983* double rsun_ref
1984* (Given, auxiliary) Reference radius of the Sun used in coordinate
1985* calculations (m).
1986*
1987* double dsun_obs
1988* (Given, auxiliary) Distance between the centre of the Sun and the
1989* observer (m).
1990*
1991* double crln_obs
1992* (Given, auxiliary) Carrington heliographic longitude of the observer
1993* (deg).
1994*
1995* double hgln_obs
1996* (Given, auxiliary) Stonyhurst heliographic longitude of the observer
1997* (deg).
1998*
1999* double hglt_obs
2000* (Given, auxiliary) Heliographic latitude (Carrington or Stonyhurst) of
2001* the observer (deg).
2002*
2003* double a_radius
2004* Length of the semi-major axis of a triaxial ellipsoid approximating the
2005* shape of a body (e.g. planet) in the solar system (m).
2006*
2007* double b_radius
2008* Length of the intermediate axis, normal to the semi-major and semi-minor
2009* axes, of a triaxial ellipsoid approximating the shape of a body (m).
2010*
2011* double c_radius
2012* Length of the semi-minor axis, normal to the semi-major axis, of a
2013* triaxial ellipsoid approximating the shape of a body (m).
2014*
2015* double blon_obs
2016* Bodycentric longitude of the observer in the coordinate system fixed to
2017* the planet or other solar system body (deg, in range 0 to 360).
2018*
2019* double blat_obs
2020* Bodycentric latitude of the observer in the coordinate system fixed to
2021* the planet or other solar system body (deg).
2022*
2023* double bdis_obs
2024* Bodycentric distance of the observer (m).
2025*
2026*
2027* Global variable: const char *wcs_errmsg[] - Status return messages
2028* ------------------------------------------------------------------
2029* Error messages to match the status value returned from each function.
2030*
2031*===========================================================================*/
2032
2033#ifndef WCSLIB_WCS
2034#define WCSLIB_WCS
2035
2036#include "lin.h"
2037#include "cel.h"
2038#include "spc.h"
2039
2040#ifdef __cplusplus
2041extern "C" {
2042#define wtbarr wtbarr_s // See prologue of wtbarr.h.
2043#endif
2044
2046 WCSENQ_MEM = 1, // wcsprm struct memory is managed by WCSLIB.
2047 WCSENQ_SET = 2, // wcsprm struct has been set up.
2048 WCSENQ_BYP = 4, // wcsprm struct is in bypass mode.
2049 WCSENQ_CHK = 8, // wcsprm struct is self-consistent.
2050};
2051
2052#define WCSSUB_LONGITUDE 0x1001
2053#define WCSSUB_LATITUDE 0x1002
2054#define WCSSUB_CUBEFACE 0x1004
2055#define WCSSUB_CELESTIAL 0x1007
2056#define WCSSUB_SPECTRAL 0x1008
2057#define WCSSUB_STOKES 0x1010
2058#define WCSSUB_TIME 0x1020
2059
2060
2061#define WCSCOMPARE_ANCILLARY 0x0001
2062#define WCSCOMPARE_TILING 0x0002
2063#define WCSCOMPARE_CRPIX 0x0004
2064
2065
2066extern const char *wcs_errmsg[];
2067
2069 WCSERR_SUCCESS = 0, // Success.
2070 WCSERR_NULL_POINTER = 1, // Null wcsprm pointer passed.
2071 WCSERR_MEMORY = 2, // Memory allocation failed.
2072 WCSERR_SINGULAR_MTX = 3, // Linear transformation matrix is singular.
2073 WCSERR_BAD_CTYPE = 4, // Inconsistent or unrecognized coordinate
2074 // axis type.
2075 WCSERR_BAD_PARAM = 5, // Invalid parameter value.
2076 WCSERR_BAD_COORD_TRANS = 6, // Unrecognized coordinate transformation
2077 // parameter.
2078 WCSERR_ILL_COORD_TRANS = 7, // Ill-conditioned coordinate transformation
2079 // parameter.
2080 WCSERR_BAD_PIX = 8, // One or more of the pixel coordinates were
2081 // invalid.
2082 WCSERR_BAD_WORLD = 9, // One or more of the world coordinates were
2083 // invalid.
2084 WCSERR_BAD_WORLD_COORD = 10, // Invalid world coordinate.
2085 WCSERR_NO_SOLUTION = 11, // No solution found in the specified
2086 // interval.
2087 WCSERR_BAD_SUBIMAGE = 12, // Invalid subimage specification.
2088 WCSERR_NON_SEPARABLE = 13, // Non-separable subimage coordinate system.
2089 WCSERR_UNSET = 14 // wcsprm struct is unset.
2091
2092
2093// Struct used for storing PVi_ma keywords.
2094struct pvcard {
2095 int i; // Axis number, as in PVi_ma (1-relative).
2096 int m; // Parameter number, ditto (0-relative).
2097 double value; // Parameter value.
2098};
2099
2100// Size of the pvcard struct in int units, used by the Fortran wrappers.
2101#define PVLEN (sizeof(struct pvcard)/sizeof(int))
2102
2103// Struct used for storing PSi_ma keywords.
2104struct pscard {
2105 int i; // Axis number, as in PSi_ma (1-relative).
2106 int m; // Parameter number, ditto (0-relative).
2107 char value[72]; // Parameter value.
2108};
2109
2110// Size of the pscard struct in int units, used by the Fortran wrappers.
2111#define PSLEN (sizeof(struct pscard)/sizeof(int))
2112
2113// Struct used to hold additional auxiliary parameters.
2114struct auxprm {
2115 double rsun_ref; // Solar radius.
2116 double dsun_obs; // Distance from Sun centre to observer.
2117 double crln_obs; // Carrington heliographic lng of observer.
2118 double hgln_obs; // Stonyhurst heliographic lng of observer.
2119 double hglt_obs; // Heliographic latitude of observer.
2120
2121 double a_radius; // Semi-major axis of solar system body.
2122 double b_radius; // Semi-intermediate axis of solar system body.
2123 double c_radius; // Semi-minor axis of solar system body.
2124 double blon_obs; // Bodycentric longitude of observer.
2125 double blat_obs; // Bodycentric latitude of observer.
2126 double bdis_obs; // Bodycentric distance of observer.
2127 double dummy[2]; // Reserved for future use.
2128};
2129
2130// Size of the auxprm struct in int units, used by the Fortran wrappers.
2131#define AUXLEN (sizeof(struct auxprm)/sizeof(int))
2132
2133
2134struct wcsprm {
2135 // Initialization flag (see the prologue above).
2136 //--------------------------------------------------------------------------
2137 int flag; // Set to zero to force initialization.
2138
2139 // FITS header keyvalues to be provided (see the prologue above).
2140 //--------------------------------------------------------------------------
2141 int naxis; // Number of axes (pixel and coordinate).
2142 double *crpix; // CRPIXja keyvalues for each pixel axis.
2143 double *pc; // PCi_ja linear transformation matrix.
2144 double *cdelt; // CDELTia keyvalues for each coord axis.
2145 double *crval; // CRVALia keyvalues for each coord axis.
2146
2147 char (*cunit)[72]; // CUNITia keyvalues for each coord axis.
2148 char (*ctype)[72]; // CTYPEia keyvalues for each coord axis.
2149
2150 double lonpole; // LONPOLEa keyvalue.
2151 double latpole; // LATPOLEa keyvalue.
2152
2153 double restfrq; // RESTFRQa keyvalue.
2154 double restwav; // RESTWAVa keyvalue.
2155
2156 int npv; // Number of PVi_ma keywords, and the
2157 int npvmax; // number for which space was allocated.
2158 struct pvcard *pv; // PVi_ma keywords for each i and m.
2159
2160 int nps; // Number of PSi_ma keywords, and the
2161 int npsmax; // number for which space was allocated.
2162 struct pscard *ps; // PSi_ma keywords for each i and m.
2163
2164 // Alternative header keyvalues (see the prologue above).
2165 //--------------------------------------------------------------------------
2166 double *cd; // CDi_ja linear transformation matrix.
2167 double *crota; // CROTAi keyvalues for each coord axis.
2168 int altlin; // Alternative representations
2169 // Bit 0: PCi_ja is present,
2170 // Bit 1: CDi_ja is present,
2171 // Bit 2: CROTAi is present.
2172 int velref; // AIPS velocity code, VELREF.
2173
2174 // Auxiliary coordinate system information of a general nature. Not
2175 // used by WCSLIB. Refer to the prologue comments above for a brief
2176 // explanation of these values.
2177 char alt[4];
2179 int *colax;
2180 // Auxiliary coordinate axis information.
2181 char (*cname)[72];
2182 double *crder;
2183 double *csyer;
2184 double *czphs;
2185 double *cperi;
2186
2187 char wcsname[72];
2188 // Time reference system and measurement.
2189 char timesys[72], trefpos[72], trefdir[72], plephem[72];
2190 char timeunit[72];
2191 char dateref[72];
2192 double mjdref[2];
2193 double timeoffs;
2194 // Data timestamps and durations.
2195 char dateobs[72], datebeg[72], dateavg[72], dateend[72];
2198 double tstart, tstop;
2200 // Timing accuracy.
2203 // Spatial & celestial reference frame.
2204 double obsgeo[6];
2205 char obsorbit[72];
2206 char radesys[72];
2207 double equinox;
2208 char specsys[72];
2209 char ssysobs[72];
2210 double velosys;
2211 double zsource;
2212 char ssyssrc[72];
2213 double velangl;
2214
2215 // Additional auxiliary coordinate system information of a specialist
2216 // nature. Not used by WCSLIB. Refer to the prologue comments above.
2217 struct auxprm *aux;
2218
2219 // Coordinate lookup tables (see the prologue above).
2220 //--------------------------------------------------------------------------
2221 int ntab; // Number of separate tables.
2222 int nwtb; // Number of wtbarr structs.
2223 struct tabprm *tab; // Tabular transformation parameters.
2224 struct wtbarr *wtb; // Array of wtbarr structs.
2225
2226 //--------------------------------------------------------------------------
2227 // Information derived from the FITS header keyvalues by wcsset().
2228 //--------------------------------------------------------------------------
2229 char lngtyp[8], lattyp[8]; // Celestial axis types, e.g. RA, DEC.
2230 int lng, lat, spec, time; // Longitude, latitude, spectral, and time
2231 // axis indices (0-relative).
2232 int cubeface; // True if there is a CUBEFACE axis.
2233 int chksum; // Checksum of keyvalues provided.
2234 int *types; // Coordinate type codes for each axis.
2235
2236 struct linprm lin; // Linear transformation parameters.
2237 struct celprm cel; // Celestial transformation parameters.
2238 struct spcprm spc; // Spectral transformation parameters.
2239
2240 //--------------------------------------------------------------------------
2241 // THE REMAINDER OF THE WCSPRM STRUCT IS PRIVATE.
2242 //--------------------------------------------------------------------------
2243
2244 // Error handling, if enabled.
2245 //--------------------------------------------------------------------------
2246 struct wcserr *err;
2247
2248 // Memory management.
2249 //--------------------------------------------------------------------------
2252 char (*m_cunit)[72], (*m_ctype)[72];
2253 struct pvcard *m_pv;
2254 struct pscard *m_ps;
2255 double *m_cd, *m_crota;
2257 char (*m_cname)[72];
2259 struct auxprm *m_aux;
2260 struct tabprm *m_tab;
2261 struct wtbarr *m_wtb;
2262};
2263
2264// Size of the wcsprm struct in int units, used by the Fortran wrappers.
2265#define WCSLEN (sizeof(struct wcsprm)/sizeof(int))
2266
2267
2268int wcsnpv(int n);
2269
2270int wcsnps(int n);
2271
2272int wcsini(int alloc, int naxis, struct wcsprm *wcs);
2273
2274int wcsinit(int alloc, int naxis, struct wcsprm *wcs, int npvmax, int npsmax,
2275 int ndpmax);
2276
2277int wcsauxi(int alloc, struct wcsprm *wcs);
2278
2279int wcssub(int alloc, const struct wcsprm *wcssrc, int *nsub, int axes[],
2280 struct wcsprm *wcsdst);
2281
2282int wcscompare(int cmp, double tol, const struct wcsprm *wcs1,
2283 const struct wcsprm *wcs2, int *equal);
2284
2285int wcsfree(struct wcsprm *wcs);
2286
2287int wcstrim(struct wcsprm *wcs);
2288
2289int wcssize(const struct wcsprm *wcs, int sizes[2]);
2290
2291int auxsize(const struct auxprm *aux, int sizes[2]);
2292
2293int wcsenq(const struct wcsprm *wcs, int enquiry);
2294
2295int wcsprt(const struct wcsprm *wcs);
2296
2297int wcsperr(const struct wcsprm *wcs, const char *prefix);
2298
2299int wcsbchk(struct wcsprm *wcs, int bounds);
2300
2301int wcsset(struct wcsprm *wcs);
2302
2303int wcsp2s(struct wcsprm *wcs, int ncoord, int nelem, const double pixcrd[],
2304 double imgcrd[], double phi[], double theta[], double world[],
2305 int stat[]);
2306
2307int wcss2p(struct wcsprm *wcs, int ncoord, int nelem, const double world[],
2308 double phi[], double theta[], double imgcrd[], double pixcrd[],
2309 int stat[]);
2310
2311int wcsmix(struct wcsprm *wcs, int mixpix, int mixcel, const double vspan[2],
2312 double vstep, int viter, double world[], double phi[],
2313 double theta[], double imgcrd[], double pixcrd[]);
2314
2315int wcsccs(struct wcsprm *wcs, double lng2p1, double lat2p1, double lng1p2,
2316 const char *clng, const char *clat, const char *radesys,
2317 double equinox, const char *alt);
2318
2319int wcssptr(struct wcsprm *wcs, int *i, char ctype[9]);
2320
2321const char* wcslib_version(int vers[3]);
2322
2323// Defined mainly for backwards compatibility, use wcssub() instead.
2324#define wcscopy(alloc, wcssrc, wcsdst) wcssub(alloc, wcssrc, 0x0, 0x0, wcsdst)
2325
2326
2327// Deprecated.
2328#define wcsini_errmsg wcs_errmsg
2329#define wcssub_errmsg wcs_errmsg
2330#define wcscopy_errmsg wcs_errmsg
2331#define wcsfree_errmsg wcs_errmsg
2332#define wcsprt_errmsg wcs_errmsg
2333#define wcsset_errmsg wcs_errmsg
2334#define wcsp2s_errmsg wcs_errmsg
2335#define wcss2p_errmsg wcs_errmsg
2336#define wcsmix_errmsg wcs_errmsg
2337
2338#ifdef __cplusplus
2339#undef wtbarr
2340}
2341#endif
2342
2343#endif // WCSLIB_WCS
Additional auxiliary parameters.
Definition wcs.h:2114
double dsun_obs
Definition wcs.h:2116
double hglt_obs
Definition wcs.h:2119
double c_radius
Definition wcs.h:2123
double hgln_obs
Definition wcs.h:2118
double dummy[2]
Definition wcs.h:2127
double crln_obs
Definition wcs.h:2117
double b_radius
Definition wcs.h:2122
double blon_obs
Definition wcs.h:2124
double a_radius
Definition wcs.h:2121
double rsun_ref
Definition wcs.h:2115
double bdis_obs
Definition wcs.h:2126
double blat_obs
Definition wcs.h:2125
Celestial transformation parameters.
Definition cel.h:456
Linear transformation parameters.
Definition lin.h:718
Store for PSi_ma keyrecords.
Definition wcs.h:2104
int i
Definition wcs.h:2105
int m
Definition wcs.h:2106
char value[72]
Definition wcs.h:2107
Store for PVi_ma keyrecords.
Definition wcs.h:2094
double value
Definition wcs.h:2097
int i
Definition wcs.h:2095
int m
Definition wcs.h:2096
Spectral transformation parameters.
Definition spc.h:869
Tabular transformation parameters.
Definition tab.h:611
Error message handling.
Definition wcserr.h:243
Coordinate transformation parameters.
Definition wcs.h:2134
char timeunit[72]
Definition wcs.h:2190
struct pscard * m_ps
Definition wcs.h:2254
char timesys[72]
Definition wcs.h:2189
struct pvcard * pv
Definition wcs.h:2158
char dateref[72]
Definition wcs.h:2191
double mjdavg
Definition wcs.h:2196
int lng
Definition wcs.h:2230
double * czphs
Definition wcs.h:2184
char(* m_cname)[72]
Definition wcs.h:2257
double zsource
Definition wcs.h:2211
double * m_crder
Definition wcs.h:2258
int npv
Definition wcs.h:2156
double timrder
Definition wcs.h:2201
char trefpos[72]
Definition wcs.h:2189
double telapse
Definition wcs.h:2199
double * m_csyer
Definition wcs.h:2258
double tstart
Definition wcs.h:2198
double * csyer
Definition wcs.h:2183
double * m_crpix
Definition wcs.h:2251
double * cperi
Definition wcs.h:2185
double mjdend
Definition wcs.h:2196
char wcsname[72]
Definition wcs.h:2187
struct tabprm * tab
Definition wcs.h:2223
struct linprm lin
Definition wcs.h:2236
double * pc
Definition wcs.h:2143
int flag
Definition wcs.h:2137
struct auxprm * aux
Definition wcs.h:2217
int npsmax
Definition wcs.h:2161
double * m_cdelt
Definition wcs.h:2251
double timeoffs
Definition wcs.h:2193
double * crder
Definition wcs.h:2182
int nps
Definition wcs.h:2160
int * m_colax
Definition wcs.h:2256
double * m_crval
Definition wcs.h:2251
double timepixr
Definition wcs.h:2202
double * m_crota
Definition wcs.h:2255
double * m_cperi
Definition wcs.h:2258
int m_flag
Definition wcs.h:2250
char lngtyp[8]
Definition wcs.h:2229
double restwav
Definition wcs.h:2154
double latpole
Definition wcs.h:2151
int m_naxis
Definition wcs.h:2250
int chksum
Definition wcs.h:2233
char radesys[72]
Definition wcs.h:2206
double * m_pc
Definition wcs.h:2251
struct pvcard * m_pv
Definition wcs.h:2253
int naxis
Definition wcs.h:2141
double * m_czphs
Definition wcs.h:2258
int * colax
Definition wcs.h:2179
double * crval
Definition wcs.h:2145
double * m_cd
Definition wcs.h:2255
double tstop
Definition wcs.h:2198
double obsgeo[6]
Definition wcs.h:2204
double timsyer
Definition wcs.h:2201
int nwtb
Definition wcs.h:2222
char ssyssrc[72]
Definition wcs.h:2212
double equinox
Definition wcs.h:2207
int altlin
Definition wcs.h:2168
struct wtbarr * wtb
Definition wcs.h:2224
int npvmax
Definition wcs.h:2157
char(* cname)[72]
Definition wcs.h:2181
double jepoch
Definition wcs.h:2197
int ntab
Definition wcs.h:2221
char ssysobs[72]
Definition wcs.h:2209
struct pscard * ps
Definition wcs.h:2162
int colnum
Definition wcs.h:2178
char datebeg[72]
Definition wcs.h:2195
double velangl
Definition wcs.h:2213
double mjdbeg
Definition wcs.h:2196
char(* cunit)[72]
Definition wcs.h:2147
double bepoch
Definition wcs.h:2197
double mjdref[2]
Definition wcs.h:2192
int time
Definition wcs.h:2230
char trefdir[72]
Definition wcs.h:2189
char plephem[72]
Definition wcs.h:2189
char dateobs[72]
Definition wcs.h:2195
double * crpix
Definition wcs.h:2142
int * types
Definition wcs.h:2234
int lat
Definition wcs.h:2230
int spec
Definition wcs.h:2230
char specsys[72]
Definition wcs.h:2208
double mjdobs
Definition wcs.h:2196
double xposure
Definition wcs.h:2199
int velref
Definition wcs.h:2172
struct celprm cel
Definition wcs.h:2237
struct wtbarr * m_wtb
Definition wcs.h:2261
char dateend[72]
Definition wcs.h:2195
struct auxprm * m_aux
Definition wcs.h:2259
double restfrq
Definition wcs.h:2153
double * cdelt
Definition wcs.h:2144
int cubeface
Definition wcs.h:2232
struct tabprm * m_tab
Definition wcs.h:2260
char obsorbit[72]
Definition wcs.h:2205
char(* ctype)[72]
Definition wcs.h:2148
char lattyp[8]
Definition wcs.h:2229
char dateavg[72]
Definition wcs.h:2195
char alt[4]
Definition wcs.h:2177
struct spcprm spc
Definition wcs.h:2238
double timedel
Definition wcs.h:2202
double * crota
Definition wcs.h:2167
char(* m_cunit)[72]
Definition wcs.h:2252
double velosys
Definition wcs.h:2210
struct wcserr * err
Definition wcs.h:2246
double lonpole
Definition wcs.h:2150
double * cd
Definition wcs.h:2166
Extraction of coordinate lookup tables from BINTABLE.
Definition getwcstab.h:167
int i
Definition getwcstab.h:168
wcs_errmsg_enum
Definition wcs.h:2068
@ WCSERR_BAD_WORLD
Definition wcs.h:2082
@ WCSERR_BAD_PIX
Definition wcs.h:2080
@ WCSERR_SINGULAR_MTX
Definition wcs.h:2072
@ WCSERR_NON_SEPARABLE
Definition wcs.h:2088
@ WCSERR_BAD_CTYPE
Definition wcs.h:2073
@ WCSERR_MEMORY
Definition wcs.h:2071
@ WCSERR_BAD_WORLD_COORD
Definition wcs.h:2084
@ WCSERR_BAD_COORD_TRANS
Definition wcs.h:2076
@ WCSERR_NO_SOLUTION
Definition wcs.h:2085
@ WCSERR_BAD_SUBIMAGE
Definition wcs.h:2087
@ WCSERR_SUCCESS
Definition wcs.h:2069
@ WCSERR_UNSET
Definition wcs.h:2089
@ WCSERR_NULL_POINTER
Definition wcs.h:2070
@ WCSERR_ILL_COORD_TRANS
Definition wcs.h:2078
@ WCSERR_BAD_PARAM
Definition wcs.h:2075
int wcsmix(struct wcsprm *wcs, int mixpix, int mixcel, const double vspan[2], double vstep, int viter, double world[], double phi[], double theta[], double imgcrd[], double pixcrd[])
Hybrid coordinate transformation.
int wcsp2s(struct wcsprm *wcs, int ncoord, int nelem, const double pixcrd[], double imgcrd[], double phi[], double theta[], double world[], int stat[])
Pixel-to-world transformation.
int wcsini(int alloc, int naxis, struct wcsprm *wcs)
Default constructor for the wcsprm struct.
int wcssize(const struct wcsprm *wcs, int sizes[2])
Compute the size of a wcsprm struct.
int wcsnpv(int n)
Memory allocation for PVi_ma.
int wcsfree(struct wcsprm *wcs)
Destructor for the wcsprm struct.
int auxsize(const struct auxprm *aux, int sizes[2])
Compute the size of a auxprm struct.
int wcssptr(struct wcsprm *wcs, int *i, char ctype[9])
Spectral axis translation.
int wcss2p(struct wcsprm *wcs, int ncoord, int nelem, const double world[], double phi[], double theta[], double imgcrd[], double pixcrd[], int stat[])
World-to-pixel transformation.
int wcscompare(int cmp, double tol, const struct wcsprm *wcs1, const struct wcsprm *wcs2, int *equal)
Compare two wcsprm structs for equality.
int wcsccs(struct wcsprm *wcs, double lng2p1, double lat2p1, double lng1p2, const char *clng, const char *clat, const char *radesys, double equinox, const char *alt)
Change celestial coordinate system.
int wcssub(int alloc, const struct wcsprm *wcssrc, int *nsub, int axes[], struct wcsprm *wcsdst)
Subimage extraction routine for the wcsprm struct.
int wcsperr(const struct wcsprm *wcs, const char *prefix)
Print error messages from a wcsprm struct.
int wcsauxi(int alloc, struct wcsprm *wcs)
Default constructor for the auxprm struct.
int wcsinit(int alloc, int naxis, struct wcsprm *wcs, int npvmax, int npsmax, int ndpmax)
Default constructor for the wcsprm struct.
int wcstrim(struct wcsprm *wcs)
Free unused arrays in the wcsprm struct.
int wcsprt(const struct wcsprm *wcs)
Print routine for the wcsprm struct.
const char * wcs_errmsg[]
Status return messages.
const char * wcslib_version(int vers[3])
wcsenq_enum
Definition wcs.h:2045
@ WCSENQ_SET
Definition wcs.h:2047
@ WCSENQ_CHK
Definition wcs.h:2049
@ WCSENQ_MEM
Definition wcs.h:2046
@ WCSENQ_BYP
Definition wcs.h:2048
int wcsset(struct wcsprm *wcs)
Setup routine for the wcsprm struct.
int wcsnps(int n)
Memory allocation for PSi_ma.
int wcsenq(const struct wcsprm *wcs, int enquiry)
enquire about the state of a wcsprm struct.
int wcsbchk(struct wcsprm *wcs, int bounds)
Enable/disable bounds checking.