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
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
|
subroutine readobs (line)
c******************************************************************************
c this routine reads in observed spectra. the format is indicated by
c the value of "specfileopt"
c******************************************************************************
implicit real*8 (a-h,o-z)
include 'Atmos.com'
include 'Linex.com'
include 'Pstuff.com'
include 'Equivs.com'
include 'Obspars.com'
byte int1(2880), onebyte
integer*2 int2(1440)
integer*4 int4(720)
real*4 real4(720)
real*4 xtemp, ytemp
equivalence (int1,int2,int4,real4)
character head*2880
c*****first set some initial parameters
write (obsitle,1021)
iend = 0
naxis = 1
naxis1 = 0
bscale = 0.
bzero = 0.
do j=1,9
disp(j) = 0.
enddo
c*****branch to the appropriate data file type
go to (100,200,300,400,500), specfileopt
c*****here is a pure FITS file type. Reading is done with ordinary
c FORTRAN read statements
c first get the header records and search for key parameters
100 irec = 1
101 read (unit=nfobs,rec=irec,err=1002,iostat=ierr) head
call obshead (head,iend,line)
if (lount == -1) return
if (iend == 0) then
irec = irec + 1
go to 101
endif
c next read the flux array from the file
nrec = lount/nblock
ipt = 0
if (mod(lount,nblock) /= 0) nrec = nrec + 1
do j=1,nrec
irec = irec + 1
jpt = min0(nblock,lount-ipt)
if (ibits == 16) then
read (unit=nfobs,rec=irec,err=1006,iostat=ierr) int1
if (byteswap == 1) then
do k=2,2880,2
onebyte = int1(k)
int1(k) = int1(k-1)
int1(k-1) = onebyte
enddo
endif
do k=1,jpt
yobs(ipt+k) = bzero + bscale*real(int2(k))
enddo
elseif (ibits == 32) then
read (unit=nfobs,rec=irec,err=1006,iostat=ierr) int1
if (byteswap == 1) then
do k=4,2880,4
onebyte = int1(k)
int1(k) = int1(k-3)
int1(k-3) = onebyte
onebyte = int1(k-1)
int1(k-1) = int1(k-2)
int1(k-2) = onebyte
enddo
endif
do k=1,jpt
yobs(ipt+k) = bzero + bscale*real(int4(k))
enddo
elseif(ibits == -32) then
read (unit=nfobs,rec=irec,err=1006,iostat=ierr) int1
if (byteswap == 1) then
do k=4,2880,4
onebyte = int1(k)
int1(k) = int1(k-3)
int1(k-3) = onebyte
onebyte = int1(k-1)
int1(k-1) = int1(k-2)
int1(k-2) = onebyte
enddo
endif
do k=1,jpt
yobs(ipt+k) = real4(k)
enddo
endif
ipt = ipt + jpt
enddo
c now adjust the continuum, if desired
do j=1,lount
yobs(j) = yobs(j)*contnorm
enddo
c now fill in the wavelength array
do j=1,lount
xobs(j) = sngl(wavecalc(sngl(dfloat(j)),lount,disp))
enddo
lim1obs = 1
lim2obs = lount
if(xobs(2) < xobs(1)) then
do j=1,lount/2
xtemp = xobs(j)
ytemp = yobs(j)
xobs(j) = xobs(lount+1-j)
yobs(j) = yobs(lount+1-j)
xobs(lount+1-j) = xtemp
yobs(lount+1-j) = ytemp
enddo
endif
return
c or if trouble, give error messages, and abort
1002 write(array,1003) ierr
istat = ivwrite(line+2,3,array,39)
go to 1007
1006 write(array,1004) ierr
istat = ivwrite(line+2,3,array,39)
1007 lount = -1
return
c*****this is a dummy for future file use
200 return
c*****this is a dummy for future file use
300 return
c*****this is a dummy for future file use
400 return
c*****here is a MONGO-style input array
500 rewind nfobs
read (nfobs,5001) obsitle
i = 1
501 read (nfobs,*,end=525) xobs(i),yobs(i)
i = i + 1
go to 501
525 lount = i - 1
lim1obs = 1
lim2obs = lount
if(xobs(2) < xobs(1)) then
do j=1,lount/2
xtemp = xobs(j)
ytemp = yobs(j)
xobs(j) = xobs(lount+1-j)
yobs(j) = yobs(lount+1-j)
xobs(lount+1-j) = xtemp
yobs(lount+1-j) = ytemp
enddo
endif
c now adjust the continuum, if desired
do j=1,lount
yobs(j) = yobs(j)*contnorm
enddo
return
c*****format statements
1003 format('ERROR IN FITS HEADER READ: ERROR=',i4)
1004 format('ERROR IN DATA READ: ERROR=',i4,'/')
1021 format (80(' '))
5001 format (a80)
end
|
