1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
|
include <mach.h>
include "../lib/apphot.h"
include "../lib/center.h"
define RADIUS 15.0
define CRADIUS 5
# AP_QRAD -- Simple radial profile plotter.
procedure ap_qrad (ap, im, wx, wy, gd)
pointer ap # pointer to apphot structure
pointer im # pointero to the IRAF image
real wx, wy # cursor coordinates
pointer gd # pointer to graphics stream
real gwx, gwy, xcenter, ycenter, xc, yc, radius, rmin, rmax, imin, imax
real u1, u2, v1, v2, x1, x2, y1, y2, xold, yold
pointer gt, sp, pix, coords, index, r, cmd
int maxpix, npix, nx, ny, wcs, key, niter
real apstatr()
pointer ap_gtinit()
int ap_skypix(), clgcur(), apstati()
int nscan(), scan()
begin
# Check for open graphics stream.
if (gd == NULL)
return
call greactivate (gd, 0)
call gclear (gd)
call gflush (gd)
# Get the radius of the extraction region.
call printf ("Half width of extraction box (%4.1f) pixels:")
call pargr (RADIUS)
call flush (STDOUT)
if (scan () == EOF)
radius = RADIUS
else {
call gargr (radius)
if (nscan () < 1)
radius = RADIUS
}
maxpix = (2 * int (radius) + 1) ** 2
# Allocate temporary space.
call smark (sp)
call salloc (coords, maxpix, TY_INT)
call salloc (index, maxpix, TY_INT)
call salloc (pix, maxpix, TY_REAL)
call salloc (r, maxpix, TY_REAL)
call salloc (cmd, SZ_LINE, TY_CHAR)
# Fit the center using 3 iterations.
xold = wx
yold = wy
niter = 0
repeat {
call ap_ictr (im, xold, yold, CRADIUS, apstati (ap,
POSITIVE), xcenter, ycenter)
niter = niter + 1
if (abs (xcenter - xold) <= 1.0 && abs (ycenter - yold) <= 1.0)
break
xold = xcenter
yold = ycenter
} until (niter >= 3)
# Fetch the pixels for the radial profile.
npix = ap_skypix (im, xcenter, ycenter, 0.0, radius, Memr[pix],
Memi[coords], xc, yc, nx, ny)
if (npix <= 0) {
call gdeactivate (gd, 0)
call sfree (sp)
return
}
call ap_index (Memi[index], npix)
# Store old viewport and window coordinates.
call ggview (gd, u1, u2, v1, v2)
call ggwind (gd, x1, x2, y1, y2)
# Initialize the plot and store the viewport and window limits.
#call apstats (ap, IMNAME, Memc[cmd], SZ_FNAME)
call apstats (ap, IMROOT, Memc[cmd], SZ_FNAME)
call ap_ltov (im, xcenter, ycenter, xcenter, ycenter, 1)
gt = ap_gtinit (Memc[cmd], xcenter, ycenter)
# Compute the radius values.
call ap_xytor (Memi[coords], Memi[index], Memr[r], npix, xc, yc, nx)
call alimr (Memr[r], npix, rmin, rmax)
call alimr (Memr[pix], npix, imin, imax)
# Plot radial profiles.
call gclear (gd)
call ap_rset (gd, gt, 0.0, rmax, imin, imax, apstatr (ap, SCALE))
call ap_plotrad (gd, gt, Memr[r], Memr[pix], npix, "plus")
# Go into cursor mode.
call printf ("Waiting for cursor mode command: [:.help=help,q=quit]")
while (clgcur ("gcommands", gwx, gwy, wcs, key, Memc[cmd], SZ_LINE) !=
EOF) {
if (key == 'q')
break
call printf (
"Waiting for cursor mode command: [:.help=help,q=quit]")
}
# Restore old window and viewport coordinates.
call gsview (gd, u1, u2, v1, v2)
call gswind (gd, x1, x2, y1, y2)
# Free plots and working space.
call ap_gtfree (gt)
call sfree (sp)
call gdeactivate (gd, 0)
end
|
