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
120
|
include <mach.h>
include <ctype.h>
include <error.h>
include <gset.h>
include "pdm.h"
define MSIZE 2.0 # Mark size
define EDGESCL 5.0 # Percent extra space around plot data
# PDM_PPLOT -- Plot the Phase curve.
procedure pdm_pplot (pdmp, period, filename, flip)
pointer pdmp # pointer to PDM data structure
double period # period for which to plot the phase curve
char filename[SZ_FNAME] # name of the input data file
bool flip # flip the y-axis scale
int npt, i, index
double frequency
real x1, x2, y1, y2, scldif, sclspc
pointer gp, xtemp, ytemp, etemp
pointer title, system_id, sp
begin
npt = PDM_NPT(pdmp)
gp = PDM_GP(pdmp)
call gclear (PDM_GP(pdmp))
call smark (sp)
call salloc (system_id, SZ_LINE, TY_CHAR)
call salloc (title, PDM_SZ_TITLE, TY_CHAR)
call salloc (xtemp, 2*npt, TY_REAL)
call salloc (ytemp, 2*npt, TY_REAL)
call salloc (etemp, 2*npt, TY_REAL)
do i = 1, 2*npt {
Memr[xtemp+i-1] = PDM_XPH(pdmp,i)
Memr[ytemp+i-1] = PDM_YPH(pdmp,i)
Memr[etemp+i-1] = PDM_PHERR(pdmp,i)
}
# Scale the wcs.
call gascale (gp, Memr[xtemp], 2*npt, 1)
call gascale (gp, Memr[ytemp], 2*npt, 2)
# Get the boundaries in X and Y.
call ggwind (gp, x1, x2, y1, y2)
# Add boundry space.
scldif = x2 - x1
sclspc = scldif * (EDGESCL / 100.)
x1 = x1 - sclspc
x2 = x2 + sclspc
scldif = y2 - y1
sclspc = scldif * (EDGESCL / 100.)
y1 = y1 - sclspc
y2 = y2 + sclspc
# Flip the y-axis scale if flip = TRUE.
if (flip)
call gswind (gp, x1, x2, y2, y1)
else
call gswind (gp, x1, x2, y1, y2)
# Multiline title, save in an array and sprintf to it.
# Get the system identification.
call sysid (Memc[system_id], SZ_LINE)
# Calculate the frequency.
if (period != 0.0)
frequency = 1.0d+0/period
call sprintf (Memc[title], PDM_SZ_TITLE,
"%s\nPhase curve at period %12.12g, frequency %12.12g\nFile = %s, numpts = %d")
call pargstr (Memc[system_id])
call pargd (period)
call pargd (frequency)
call pargstr (filename)
call pargi (npt)
# Draw the axes.
call glabax (gp, Memc[title], "phase", "magnitude")
# Make the plot. If error bars are turned on, draw them.
if (PDM_EB(pdmp) == YES && Memr[etemp] > EPSILONR) {
call gpmark (gp, Memr[xtemp], Memr[ytemp],
2*npt, GM_CIRCLE, 1.0, 1.0)
call gpmark (gp, Memr[xtemp], Memr[ytemp],
2*npt, GM_CIRCLE+GM_FILL, 1.0, 1.0)
do i = 1, 2*npt
call gpmark (gp, Memr[xtemp+i-1], Memr[ytemp+i-1],
1, GM_VEBAR, MSIZE, -(2.0*Memr[etemp+i-1]))
} else {
if (npt <= PDM_PLUSPOINT(pdmp)) {
call gpmark (gp, Memr[xtemp], Memr[ytemp],
2*npt, GM_PLUS, MSIZE, MSIZE)
} else {
call gpmark (gp, Memr[xtemp], Memr[ytemp],
2*npt, GM_POINT, 1.0, 1.0)
}
}
# Draw an x over any deleted points.
do i = 1, npt {
index = PDM_SORT(pdmp,i)
if (PDM_INUSE(pdmp,index) == 0) {
call gscur (gp, Memr[xtemp+i-1]+1, Memr[ytemp+i-1])
call gmark (gp, Memr[xtemp+i-1]+1, Memr[ytemp+i-1], GM_CROSS,
MSIZE, MSIZE)
call gscur (gp, Memr[xtemp+i-1], Memr[ytemp+i-1])
call gmark (gp, Memr[xtemp+i-1], Memr[ytemp+i-1], GM_CROSS,
MSIZE, MSIZE)
}
}
call sfree (sp)
end
|
