Attachment 'stcschan-demo3.c'
Download 1 /* Name:
2 stcschan-demo3.c
3
4 Purpose:
5 A demonstration of the facilities provided by the AST library
6 for reading STC metadata encoded using the STC-S linear string
7 format.
8
9 Description:
10 This program reads an STC-S description from a text file, and also
11 reads a set of 2-D spatial FITS-WCS headers from another (text) file.
12 It then opens a specified graphics device, and displays an annotated
13 coordinate grid covering the region described by the FITS headers.
14 Finally, it draws the outline of the spatial extent of the STC-S
15 description over the top of the annotated coordinate grid.
16
17 Usage:
18 % stcschan-demo3 <stcs-file> <header-file> <device>
19
20 <stcs-file>: The path to a text file containing the STC-S
21 description.
22
23 <header-file>: The path to a text file containing a set of FITS-WCS
24 headers.
25
26 <device>: The name of an available PGPLOT graphics device. If not
27 supplied, the available device names will be listed and the program
28 will then exit.
29
30 Example:
31 % stcschan-demo3 m31.stcs andromeda.head /xserve
32
33 To compile and link:
34 Assuming your starlink distribution is in "/star":
35
36 % gcc -o stcschan-demo3 stcschan-demo3.c -g -L/star/lib \
37 -I/star/include `ast_link -pgplot` -lcpgplot
38
39 */
40
41 /* Include system headers. */
42 #include <stdio.h>
43 #include <string.h>
44
45 /* Include the PGPLOT header. */
46 #include "cpgplot.h"
47
48 /* Include the AST library header. */
49 #include "ast.h"
50
51 /* Maximum allowed length for a single line of text from the disk file. */
52 #define MAX_LINE_LEN 100
53
54 /* Prototypes */
55 const char *source( void );
56 AstFrameSet *ReadFitsHeaders( const char *, int *, int * );
57 AstRegion *ReadStcs( const char * );
58
59
60 int main( int argc, char **argv ){
61
62 /* Local variables: */
63 AstBox *pixbox;
64 AstFrame *pixfrm;
65 AstFrameSet *align_fs = NULL;
66 AstFrameSet *fset = NULL;
67 AstPlot *plot;
68 AstRegion *reg = NULL;
69 AstRegion *wcsbox = NULL;
70 AstRegion *wcsreg = NULL;
71 const char *dev;
72 double bbox[ 4 ];
73 float gbox[ 4 ];
74 int overlap_flag, status, naxis1, naxis2, ibase;
75
76 /* Initialise the returned system status to indicate failure. */
77 status = 1;
78
79 /* Start an AST object context. This means we do not need to annull
80 each AST Object individually. Instead, all Objects created within
81 this context will be annulled automatically by the corresponding
82 invocation of astEnd. */
83 astBegin;
84
85 /* Check there are enough command line arguments. */
86 if( argc < 3 ) {
87 printf( "Usage: stcschan-demo3 <stcs-file> <header-file> <device>\n" );
88
89 /* If so, attempt to read the STC-S description, creating a corresponding
90 AST Region. If this is successful, attempt to read the FITS header file
91 and create an equivalent AST FrameSet. */
92 } else {
93 reg = ReadStcs( argv[ 1 ] );
94 if( reg ) fset = ReadFitsHeaders( argv[ 2 ], &naxis1, &naxis2 );
95 }
96
97 /* Check we obtained a Region and a FrameSet successfully. */
98 if( reg && fset ){
99
100 /* Check that we can align the FITS WCS with the spatial axes in the
101 STC-S description. AST contains various built-in conversions between
102 standard celestial coordinate system. The necessary conversion will be
103 identified and used automatically if necessary to achieve alignment.
104 The returned object "align_fs" is a FrameSet that encapsulates the
105 the FITS WCS coordinate system, the STC-S spatial coordinate system and
106 the Mapping between them. A NULL pointer is returned if alignment is
107 not possible. Note, the astConvert method changes the base Frame in
108 any supplied FrameSet to indicate which coordinate frame was used for
109 alignment. So we first note the original base Frame index and then
110 re-instate it afterwards. */
111 ibase = astGetI( fset, "Base" );
112 align_fs = astConvert( reg, fset, " " );
113 astSetI( fset, "Base", ibase );
114
115 if( !align_fs ) {
116 printf( "Could not align the FITS WCS with the spatial axes "
117 "in the STC-S\n" );
118
119 /* If alignment was possible, use the alignment FrameSet to create a new
120 Region representing the same region as the supplied STC-S, but
121 expressed in the coordinate system of the FITS WCS. The FrameSet class
122 inherits from both Frame and Mapping, and so the "align_fs" FrameSet can
123 be used both as the Mapping in this call, and as the Frame. When used
124 as a Mapping, a FrameSet represents the transformation between its base
125 and current Frame. When used as a Frame, a FrameSet represents its
126 current Frame. */
127 } else {
128 wcsreg = astMapRegion( reg, align_fs, align_fs );
129
130 /* It would be nice to warn the user if the STC-S AstroCoordsArea does
131 not overlap the FITS grid. To do so we need a Region representing the
132 FITS grid. Create one now (an AST Box). First store the bounds of the
133 FITS grid, in pixel coordinates (i.e. a system in which the centre of
134 the bottom left pixel is at (1.0,1.0) ). */
135 bbox[ 0 ] = 0.5;
136 bbox[ 1 ] = 0.5;
137 bbox[ 2 ] = (double) naxis1 + 0.5;
138 bbox[ 3 ] = (double) naxis2 + 0.5;
139
140 /* Get a pointer to the Frame describing the FITS pixel coordinate system
141 (the base Frame in the FITS FrameSet). */
142 pixfrm = astGetFrame( fset, AST__BASE );
143
144 /* Create a Box that encompasses the required range of axis values within
145 the pixel coordinate Frame. */
146 pixbox = astBox( pixfrm, 1, bbox, bbox + 2, NULL, " " );
147
148 /* Create another Region that represents the same area, but in the FITS
149 WCS. */
150 wcsbox = astMapRegion( pixbox, fset, fset );
151
152 /* If the previous step failed, it probably means the FITS header covers
153 the entire sky, resulting the corners of the pixel grid having invalid
154 sky positions. So cancel the error and omit the overlap test. */
155 if( !wcsbox ) {
156 astClearStatus;
157 printf("\nContinuing, but omitting overlap test...\n\n");
158
159 /* Now see if the Region representing the FITS grid overlaps the region
160 read from the STC-S description.*/
161 } else {
162 overlap_flag = astOverlap( wcsreg, wcsbox );
163
164 if( overlap_flag == 1 || overlap_flag == 6 ) {
165 printf( "\nThere is no overlap between the FITS grid and the "
166 "STC-S AstroCoordsArea\n\n" );
167
168 } else if( overlap_flag == 3 || overlap_flag == 5 ) {
169 printf( "\nThe FITS grid is completely contained within the "
170 "STC-S AstroCoordsArea\n\n" );
171 }
172
173 }
174 }
175 }
176
177 /* Check we obtained a FITS-WCS to STC-S Mapping, and that no error has
178 occurred in AST. */
179 if( wcsreg && astOK ){
180
181 /* Open PGPLOT using the specified device. Prompt the user for a device if
182 none was supplied on the command line. */
183 if( argc < 4 ) {
184 dev = "?";
185 } else {
186 dev = argv[ 3 ];
187 }
188 if( cpgbeg( 0, dev, 1, 1 ) == 1 ) {
189
190 /* Clear the screen. */
191 cpgpage();
192
193 /* Ensure the graphics window has equal scales on both axes. */
194 cpgwnad( 0.0f, 1.0f, 0.0f, 1.0f );
195
196 /* Find the extent of the graphics window, and store in an array suitable
197 for passing to the astPlot function. */
198 cpgqwin( gbox, gbox + 2, gbox + 1, gbox + 3 );
199
200 /* Create an AST Plot. This is a special sort of FrameSet in which the
201 base Frame corresponds to graphics coordinates. All the coordinate
202 Frames and Mappings read from the FITS-WCS headers are added into the
203 Plot so that graphics can be drawn in any coordinate system. The extent
204 of the FITS array in pixel coordinates is mapped onto the extent of the
205 graphics device as returned above by cpgqwin. The AST library comes with
206 a driver module that provides primitive drawing functions by calling
207 appropriate PGFPLOT functions. It is simple to write driver modules
208 for other graphics systems such as Tcl/Tk, Java/Swing, etc. Set a few
209 graphics attributes to show the sort of thing that can be done. */
210 plot = astPlot( fset, gbox, bbox, "Colour(border)=2,Colour(ticks)=2,"
211 "Colour(axes)=2,Grid=1,Colour(grid)=3,"
212 "Style(grid)=4" );
213
214 /* Draw a set of annotated coordinate axes labelling the FITS WCS axes. */
215 astGrid( plot );
216
217 /* If there is any overlap (or if the overlap test could not be performed),
218 add the STC-S Region into the Plot. We use a UnitMap to connect it to
219 the current Frame (the FITS WCS frame). */
220 if( 1 || overlap_flag == 2 || overlap_flag == 4 || !wcsbox ) {
221 astAddFrame( plot, AST__CURRENT, astUnitMap( 2, " " ), wcsreg );
222
223 /* Now draw the border round the STC-S Region. A Region is a sub-class of
224 Mapping and so can be used to transform positions. When a Region is used
225 as a Mapping, positions that are inside the Region are left unchanged
226 by the transformation, and positions that are outside the Region are
227 transformed into "bad" positions (i.e. every axis value has the nmagic
228 value AST__BAD indicating that the axis value is undefined). The
229 astBorder method is a generic function that will outline the areas
230 within the current coordinate Frame of the Plot that correspond to
231 valid (i.e. non-bad) positions. Set the colour and line thickness first
232 to emphasise the border. */
233 astSet( plot, "Colour(border)=4,Width(border)=8" );
234 (void) astBorder( plot );
235 }
236
237 /* Set the returned system status to indicate success. */
238 status = 0;
239
240 /* Close down PGPLOT. */
241 cpgend();
242 }
243 }
244
245 /* End the AST Object context. All Objects created since the
246 corresponding invocation of astbegin will be annulled automatically. */
247 astEnd;
248
249 /* If an error occurred in the AST library, set the returned system
250 status non-zero. */
251 if( !astOK ) status = 1;
252 return status;
253 }
254
255
256
257
258
259
260
261 /* -------------------------------------------------------------------
262 * This is a function that reads a line of text from the disk file and
263 * returns it to the AST library. It is called from within the astRead
264 * function.
265 */
266
267 const char *source( void ){
268 static char buffer[ MAX_LINE_LEN + 2 ];
269 FILE *fd = astChannelData;
270 return fgets( buffer, MAX_LINE_LEN + 2, fd );
271 }
272
273
274
275 /* -------------------------------------------------------------------
276 * This function reads a set of FITS-WCS headers from a given text file,
277 * and attempts to convert them into an AST FrameSet. If successful, a
278 * pointer to the FrameSet is returned. A NULL pointer is returned if
279 * anything goes wrong, or if the WCS is not 2-dimensional. The values of
280 * the NAXIS1 and NAXIS2 headers are returned in "*naxis1" and "*naxis2".
281 */
282
283 AstFrameSet *ReadFitsHeaders( const char *file, int *naxis1, int *naxis2 ){
284 AstFitsChan *chan;
285 AstFrameSet *result;
286 AstObject *object;
287 FILE *fd;
288
289 /* Initialise the returned pointer to indicate that no FrameSet has yet
290 been read. */
291 result = NULL;
292
293 /* Attempt to open the FITS header file */
294 fd = fopen( file, "r" );
295 if( !fd ) {
296 printf("Failed to open FITS header file '%s'.\n", file );
297
298 /* If successful, create a FitsChan. This is the object that converts
299 external FITS headers into corresponding AST Objects. Tell it to use
300 the "source" function for obtaining lines of text from the disk file. */
301 } else {
302 chan = astFitsChan( source, NULL, " " );
303
304 /* Associate the descriptor for the input disk file with the StcsChan.
305 This makes it available to the "source" function. Since this
306 application is single threaded, we could instead have made "fd" a
307 global variable, but the ChannelData facility is used here to illustrate
308 how to pass data to a source or sink function safely in a multi-threaded
309 application. */
310 astPutChannelData( chan, fd );
311
312 /* Attempt to read the FITS heades and convert them into an AST FrameSet. */
313 object = astRead( chan );
314
315 /* The astRead function is a generic function and so returns a generic
316 AstObject pointer. Check an Object was created successfully. */
317 if( !object ) {
318 printf( "Failed to read an AST Object from FITS header file '%s'.\n",
319 file );
320
321 /* Now check that the object read is actually an AST FrameSet, rather than
322 some other class of AST Object. */
323 } else if( !astIsAFrameSet( object ) ) {
324 printf( "Expected a FrameSet but read a %s from FITS header "
325 "file '%s'.\n", astGetC( object, "Class" ), file );
326
327 /* If the Object is a FrameSet, return the FrameSet pointer. */
328 } else {
329 result = (AstFrameSet *) object;
330
331 /* Check the WCS is 2-dimensional. If not, report an error and set the
332 returned pointer to NULL. The memory used to store the FrameSet will
333 be released when the current AST object context is ended (by calling
334 astEnd). */
335 if( astGetI( result, "Naxes" ) != 2 ) {
336 printf( "The FITS WCS is not 2-dimensional.\n");
337 result = NULL;
338
339 /* If it is 2-dimensional, get the NAXIS1 and NAXIS2 keyword values. */
340 } else {
341
342 if( !astGetFitsI( chan, "NAXIS1", naxis1 ) ){
343 printf("Keyword 'NAXIS1' not found in header\n" );
344 result = NULL;
345 }
346
347 if( !astGetFitsI( chan, "NAXIS2", naxis2 ) ){
348 printf("Keyword 'NAXIS2' not found in header\n" );
349 result = NULL;
350 }
351
352 }
353 }
354
355 /* Close the file. */
356 fclose( fd );
357 }
358
359 return result;
360 }
361
362
363 /* -------------------------------------------------------------------
364 * This function reads an STC-S description from a given text file,
365 * and attempts to convert them into an AST Region. If successful, a
366 * pointer to the Region is returned. A NULL pointer is returned if
367 * anything goes wrong.
368 */
369
370 AstRegion *ReadStcs( const char *file ){
371 AstStcsChan *chan;
372 AstRegion *result;
373 AstObject *object;
374 FILE *fd;
375
376 /* Initialise the returned pointer to indicate that no Region has yet
377 been read. */
378 result = NULL;
379
380 /* Attempt to open the STC-S file */
381 fd = fopen( file, "r" );
382 if( !fd ) {
383 printf("Failed to open STC-S descrption file '%s'.\n", file );
384
385 /* If successful, create an StcsChan. This is the object that converts
386 external STC-S descriptions into corresponding AST Objects. Tell it to
387 use the "source" function for obtaining lines of text from the disk
388 file. */
389 } else {
390 chan = astStcsChan( source, NULL, " " );
391
392 /* Associate the descriptor for the input disk file with the StcsChan.
393 This makes it available to the "source" function. Since this
394 application is single threaded, we could instead have made "fd" a
395 global variable, but the ChannelData facility is used here to illustrate
396 how to pass data to a source or sink function safely in a multi-threaded
397 application. */
398 astPutChannelData( chan, fd );
399
400 /* The default behaviour of the astRead function when used on an StcsChan is
401 to read and return the AstroCoordArea as an AST Region. This behaviour
402 can be changed by assigning appropriate values to the StcsChan attributes
403 "StcsArea", "StcsCoords" and "StcsProps". Options exist to return the
404 AstroCoords as an AST PointList, and/or to return the individual
405 property values read from the STC-S text in the form of an AST KeyMap
406 (a sort of hashmap). For now, just take the default action of reading the
407 AstroCoordsArea. */
408 object = astRead( chan );
409
410 /* The astRead function is a generic function and so returns a generic
411 AstObject pointer. Check an Object was created successfully. */
412 if( !object ) {
413 printf( "Failed to read an AST Object from STC-S description "
414 "file '%s'.\n", file );
415
416 /* Now check that the object read is actually an AST Region, rather than
417 some other class of AST Object. */
418 } else if( !astIsARegion( object ) ) {
419 printf( "Expected a Region but read a %s from STC-S description "
420 "file '%s'.\n", astGetC( object, "Class" ), file );
421
422 /* If the Object is a Region, return the Region pointer. */
423 } else {
424 result = (AstRegion *) object;
425 }
426
427 /* Close the file. */
428 fclose( fd );
429 }
430
431 return result;
432 }
433
434
435
436
Attached Files
To refer to attachments on a page, use attachment:filename, as shown below in the list of files. Do NOT use the URL of the [get] link, since this is subject to change and can break easily.- [get | view] (2009-09-24 10:49:44, 501.2 KB) [[attachment:adass-2007.pdf]]
- [get | view] (2009-10-27 14:29:15, 709.4 KB) [[attachment:adass2009.pdf]]
- [get | view] (2011-11-22 16:37:35, 1053.8 KB) [[attachment:adass2011.pdf]]
- [get | view] (2009-09-24 10:56:59, 11508.3 KB) [[attachment:ast-5.2-0.tar.gz]]
- [get | view] (2009-10-30 13:30:20, 11553.7 KB) [[attachment:ast-5.3-1.tar.gz]]
- [get | view] (2011-04-05 12:39:19, 11994.6 KB) [[attachment:ast-5.6-0.tar.gz]]
- [get | view] (2011-05-24 14:49:42, 12040.2 KB) [[attachment:ast-5.7-0.tar.gz]]
- [get | view] (2011-06-07 20:58:18, 12075.5 KB) [[attachment:ast-5.7-1.tar.gz]]
- [get | view] (2011-07-18 18:52:08, 12141.4 KB) [[attachment:ast-5.7-2.tar.gz]]
- [get | view] (2011-11-22 16:40:11, 12192.1 KB) [[attachment:ast-6.0-1.tar.gz]]
- [get | view] (2011-10-31 11:00:18, 12189.3 KB) [[attachment:ast-6.0.tar.gz]]
- [get | view] (2011-11-22 16:37:59, 31.7 KB) [[attachment:ast.news]]
- [get | view] (2020-09-15 16:09:39, 4.7 KB) [[attachment:example.asdf]]
- [get | view] (2009-09-24 09:34:46, 3.1 KB) [[attachment:img1.gif]]
- [get | view] (2009-09-24 09:35:03, 3.3 KB) [[attachment:img2.gif]]
- [get | view] (2009-09-24 09:35:22, 3.2 KB) [[attachment:img3.gif]]
- [get | view] (2020-09-15 16:09:08, 4.1 KB) [[attachment:readasdf.c]]
- [get | view] (2009-09-25 16:24:18, 10.0 KB) [[attachment:stcschan-demo1.c]]
- [get | view] (2009-09-25 16:24:47, 8.7 KB) [[attachment:stcschan-demo2.c]]
- [get | view] (2009-09-25 16:25:12, 15.9 KB) [[attachment:stcschan-demo3.c]]
- [get | view] (2009-09-25 16:25:41, 8.5 KB) [[attachment:stcschan-demo4.c]]
- [get | view] (2009-09-25 16:26:06, 10.4 KB) [[attachment:stcschan-demo5.c]]
- [get | view] (2022-05-06 10:57:56, 2800.5 KB) [[attachment:sun210.pdf]]
- [get | view] (2022-05-06 10:58:45, 2804.1 KB) [[attachment:sun211.pdf]]
- [get | view] (2020-09-15 16:28:09, 2.8 KB) [[attachment:writeasdf.c]]
You are not allowed to attach a file to this page.