← Back to branch summary

~ubuntu-branches/debian/jessie/eso-midas/jessie

« back to all changes in this revision

Viewing changes to stdred/spec/src/maktab.for

  • Committer: Package Import Robot
  • Author(s): Ole Streicher
  • Date: 2014-04-22 14:44:58 UTC
  • Revision ID: package-import@ubuntu.com-20140422144458-okiwi1assxkkiz39
Tags: upstream-13.09pl1.2+dfsg
Import upstream version 13.09pl1.2+dfsg

Show diffs side-by-side

added added

removed removed

Lines of Context:
stdred/spec/src/maktab.for
 
1
C @(#)maktab.for        19.1 (ES0-DMD) 02/25/03 14:29:42
 
2
C===========================================================================
 
3
C Copyright (C) 1995 European Southern Observatory (ESO)
 
4
C
 
5
C This program is free software; you can redistribute it and/or 
 
6
C modify it under the terms of the GNU General Public License as 
 
7
C published by the Free Software Foundation; either version 2 of 
 
8
C the License, or (at your option) any later version.
 
9
C
 
10
C This program is distributed in the hope that it will be useful,
 
11
C but WITHOUT ANY WARRANTY; without even the implied warranty of
 
12
C MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 
13
C GNU General Public License for more details.
 
14
C
 
15
C You should have received a copy of the GNU General Public 
 
16
C License along with this program; if not, write to the Free 
 
17
C Software Foundation, Inc., 675 Massachusetss Ave, Cambridge, 
 
18
C MA 02139, USA.
 
19
C
 
20
C Corresponding concerning ESO-MIDAS should be addressed as follows:
 
21
C       Internet e-mail: midas@eso.org
 
22
C       Postal address: European Southern Observatory
 
23
C                       Data Management Division 
 
24
C                       Karl-Schwarzschild-Strasse 2
 
25
C                       D 85748 Garching bei Muenchen 
 
26
C                       GERMANY
 
27
C===========================================================================
 
28
C
 
29
        program MAKTAB
 
30
C------------------------------------------------------------------------------
 
31
C  Makes input table in proper format for DBLEND when working with combined
 
32
C   input mode
 
33
C
 
34
C  Possible inputs:
 
35
C    fit option W, input table BTAB with center positions in column :X, input
 
36
C      table ttemp with other info in normal form.
 
37
C    fit option S, FWHM value to use in BTAB, input table ttemp with other info
 
38
C      in normal form.
 
39
C    fit option A, centers in input table ttemp, FWHM's in keyword FWHM(6)
 
40
C
 
41
C  Output:
 
42
C    table ttemp with all info in proper form for DBLEND
 
43
C------------------------------------------------------------------------------
 
44
C
 
45
C  variable declarations
 
46
C
 
47
        implicit none
 
48
 
 
49
        integer MADRID
 
50
 
 
51
        INTEGER     ACOL,AROW
 
52
        INTEGER     ACTVAL
 
53
        REAL        BDATA(6,4)
 
54
        REAL        CENT(6)
 
55
        INTEGER     COL(4)
 
56
        INTEGER     CPIX(6)
 
57
        character*16  CUNIT
 
58
        REAL        DATA(50000)
 
59
        INTEGER     DTID,BTID,TID
 
60
        character*1   FOPT
 
61
        REAL        FWHM(6)
 
62
        character*72  IDENT
 
63
        character*60  IMAGE
 
64
        character*12  INTAB
 
65
        INTEGER     KUN,KNUL
 
66
        REAL        MAXVAL(6)
 
67
        INTEGER     NAXIS
 
68
        INTEGER     NCOLD,NROWD,NCOLB,NROWB,NROWT
 
69
        INTEGER     NDATA
 
70
        INTEGER     NN
 
71
        INTEGER     NO
 
72
        INTEGER     NPIX(2)
 
73
        INTEGER     NSC
 
74
        logical       NULL(4)
 
75
        REAL        PAR(4)
 
76
        INTEGER*8   PNTR
 
77
        REAL        SIDE(2)
 
78
        DOUBLE PRECISION        START(2),STEP(2)
 
79
        INTEGER     STAT
 
80
        character*72  STRING
 
81
        REAL        TABLE(22,4)
 
82
        REAL        TDATA(16,4)
 
83
        INTEGER     XCOL
 
84
        REAL        YVAL
 
85
 
 
86
        INTEGER       II,JJ,KK,MM
 
87
 
 
88
        common /VMR/MADRID(1)
 
89
 
 
90
        data COL/1,2,3,4/
 
91
 
 
92
C  connect to MIDAS environment
 
93
 
 
94
        call stspro('MAKTAB')
 
95
 
 
96
C  get inputs from keywords
 
97
 
 
98
        call stkrdc('BTAB',1,1,12,ACTVAL,INTAB,KUN,KNUL,STAT)
 
99
        call stkrdc('FITOPT',1,1,1,ACTVAL,FOPT,KUN,KNUL,STAT)
 
100
        call stkrdc('IN_A',1,1,60,ACTVAL,IMAGE,KUN,KNUL,STAT)
 
101
 
 
102
        call upcas(FOPT,FOPT)
 
103
 
 
104
C  rename ttemp to dummy name
 
105
 
 
106
        call stfrnm('ttemp.tbl','ttemp2.tbl',STAT)
 
107
 
 
108
C  read in values to be saved
 
109
 
 
110
        call tbtopn('ttemp2.tbl',2,DTID,STAT)
 
111
        call tbiget(DTID,NCOLD,NROWD,NSC,ACOL,AROW,STAT)
 
112
 
 
113
        do  II = 1,NROWD
 
114
          call tbrrdr(DTID,II,4,COL,PAR,NULL,STAT)
 
115
          do  JJ = 1,4
 
116
            TDATA(II,JJ) = PAR(JJ)
 
117
          end do
 
118
        end do
 
119
 
 
120
C  close old table
 
121
 
 
122
        call tbtclo(DTID,STAT)
 
123
 
 
124
C  open image file
 
125
 
 
126
        call stiget(IMAGE,10,0,1,2,NAXIS,NPIX,START,STEP,
 
127
     &                IDENT,CUNIT,PNTR,NO,STAT)
 
128
 
 
129
C  for option W, open table with center positions and read them in
 
130
 
 
131
        if (FOPT.EQ.'W') then
 
132
 
 
133
          call tbtopn(INTAB,0,BTID,STAT)
 
134
          call tbiget(BTID,NCOLB,NROWB,NSC,ACOL,AROW,STAT)
 
135
          call tblser(BTID,'X',XCOL,STAT)  !could be ':X' instead
 
136
 
 
137
          do  II = 1,NROWB
 
138
            call tberdr(BTID,II,XCOL,CENT(II),NULL,STAT)
 
139
          end do
 
140
 
 
141
C  read in the image
 
142
 
 
143
          if (NAXIS.eq.1) then
 
144
            NDATA = NPIX(1)
 
145
            call read1d(MADRID(PNTR),NPIX(1),DATA)
 
146
          else if (NAXIS.eq.2.and.(NPIX(1).eq.1.or.NPIX(2).eq.1)) then
 
147
            if (NPIX(1).eq.1) then
 
148
              NDATA = NPIX(2)
 
149
              START(1) = START(2)
 
150
              STEP(1) = STEP(2)
 
151
            else
 
152
              NDATA = NPIX(1)
 
153
            end if
 
154
            call read2d(MADRID(PNTR),NPIX,DATA,NDATA)
 
155
          else
 
156
            STRING = ' Input image is not a suitable type'
 
157
            call sttput(STRING,STAT)
 
158
            go to 9999
 
159
          end if
 
160
 
 
161
C  determine pixel positions of centers
 
162
 
 
163
          do  II = 1,NROWB
 
164
            CPIX(II)=nint((CENT(II)-real(START(1)))/real(STEP(1)))+1
 
165
          end do
 
166
 
 
167
C  get data values at centers
 
168
 
 
169
          do  II = 1,NROWB
 
170
            MAXVAL(II) = DATA(CPIX(II))
 
171
          end do
 
172
 
 
173
C  construct new table rows
 
174
 
 
175
          do  II = 1,NROWB
 
176
            BDATA(II,1) = CENT(II)
 
177
            BDATA(II,2) = MAXVAL(II)
 
178
            BDATA(II,3) = float(CPIX(II))
 
179
            BDATA(II,4) = 1.0
 
180
          end do
 
181
 
 
182
C  construct table
 
183
 
 
184
          NROWT = NROWD + NROWB
 
185
 
 
186
          do  JJ = 1,4
 
187
 
 
188
            do  II = 1,4
 
189
              TABLE(II,JJ) = TDATA(II,JJ)
 
190
            end do
 
191
 
 
192
            KK = 4
 
193
 
 
194
            do  II = 5,NROWT
 
195
              if (mod(II,3).eq.0) then
 
196
                TABLE(II,JJ) = BDATA((II/3)-1,JJ)
 
197
              else
 
198
                KK = KK+1
 
199
                TABLE(II,JJ) = TDATA(KK,JJ)
 
200
              end if
 
201
            end do
 
202
          end do
 
203
 
 
204
C  for option S, interpret number in INTAB
 
205
 
 
206
        else if (FOPT.eq.'S') then
 
207
 
 
208
          read (INTAB,*) FWHM(1)                !test this
 
209
 
 
210
C  construct new table rows
 
211
 
 
212
          NROWB = 0
 
213
 
 
214
          do  II = 5,NROWD
 
215
            SIDE(1) = TDATA(II,1) - (FWHM(1)/2)
 
216
            SIDE(2) = SIDE(1) + FWHM(1)
 
217
            YVAL = TDATA(II,2) / 2
 
218
            do  JJ = 1,2
 
219
              NROWB = NROWB + 1
 
220
              BDATA(NROWB,1) = SIDE(JJ)
 
221
              BDATA(NROWB,2) = YVAL
 
222
              BDATA(NROWB,3) = anint((SIDE(JJ)-real(START(1)))
 
223
     &                         /real(STEP(1)))+1.
 
224
              BDATA(NROWB,4) = 1.0
 
225
            end do
 
226
          end do
 
227
 
 
228
C  construct new table
 
229
 
 
230
          NROWT = NROWD + NROWB
 
231
 
 
232
          do  JJ = 1,4
 
233
 
 
234
            do  II = 1,4
 
235
              TABLE(II,JJ) = TDATA(II,JJ)
 
236
            end do
 
237
 
 
238
            KK = 4
 
239
            MM = 0
 
240
 
 
241
            do  II = 5,NROWT
 
242
              if (mod(II,3).eq.0) then
 
243
                KK = KK+1
 
244
                TABLE(II,JJ) = TDATA(KK,JJ)
 
245
              else
 
246
                MM = MM+1
 
247
                TABLE(II,JJ) = BDATA(MM,JJ)
 
248
              end if
 
249
            end do
 
250
          end do
 
251
 
 
252
C  for option A, get FWHM from keyword
 
253
 
 
254
        else if (FOPT.eq.'A') then
 
255
 
 
256
          call stkrdr('FWHM',1,6,ACTVAL,FWHM,KUN,KNUL,STAT)
 
257
 
 
258
C  construct new table rows
 
259
 
 
260
          NROWB = 0
 
261
 
 
262
          do  II = 5,NROWD
 
263
            SIDE(1) = TDATA(II,1) - (FWHM(II-4)/2)
 
264
            SIDE(2) = SIDE(1) + FWHM(II-4)
 
265
            YVAL = TDATA(II,2) / 2
 
266
            do  JJ = 1,2
 
267
              NROWB = NROWB + 1
 
268
              BDATA(NROWB,1) = SIDE(JJ)
 
269
              BDATA(NROWB,2) = YVAL
 
270
              BDATA(NROWB,3) = anint((SIDE(JJ)-real(START(1)))
 
271
     &                         /real(STEP(1)))+1.
 
272
              BDATA(NROWB,4) = 1.0
 
273
            end do
 
274
          end do
 
275
 
 
276
C  construct new table
 
277
 
 
278
          NROWT = NROWD + NROWB
 
279
 
 
280
          do  JJ = 1,4
 
281
 
 
282
            do  II = 1,4
 
283
              TABLE(II,JJ) = TDATA(II,JJ)
 
284
            end do
 
285
 
 
286
            KK = 4
 
287
            MM = 0
 
288
 
 
289
            do  II = 5,NROWT
 
290
              if (mod(II,3).eq.0) then
 
291
                KK = KK+1
 
292
                TABLE(II,JJ) = TDATA(KK,JJ)
 
293
              else
 
294
                MM = MM+1
 
295
                TABLE(II,JJ) = BDATA(MM,JJ)
 
296
              end if
 
297
            end do
 
298
          end do
 
299
 
 
300
C  invalid fit option used:
 
301
 
 
302
        else
 
303
 
 
304
 
 
305
C  write out error message
 
306
 
 
307
          STRING = 'fit option '// FOPT //
 
308
     +             ' incompatible with B input option'
 
309
          call sttput(STRING,STAT)
 
310
 
 
311
C  deposit status parameter in keyword
 
312
 
 
313
          PAR(1) = 1
 
314
          call stkwri('STATUS',PAR(1),1,1,KUN,STAT)   !takes KNUL par also?
 
315
 
 
316
C  put ttemp back where it was and exit
 
317
 
 
318
          call stfrnm('ttemp2.tbl','ttemp.tbl',STAT)
 
319
 
 
320
          call stsepi
 
321
 
 
322
        end if
 
323
 
 
324
C  write the table out
 
325
 
 
326
        call tbtini('ttemp',1,1,4,NROWT,TID,STAT)
 
327
        call tbcini(TID,10,1,'G13.6',CUNIT,'X_AXIS',NN,STAT)
 
328
        call tbcini(TID,10,1,'G13.6','COUNTS','Y_AXIS',NN,STAT)
 
329
        call tbcini(TID,10,1,'I5','PIXEL','LINE_NO',NN,STAT)
 
330
        call tbcini(TID,10,1,'I5','PIXEL','PIXEL_NO',NN,STAT)
 
331
 
 
332
        do  II = 1,NROWT
 
333
          do  JJ = 1,4
 
334
            PAR(JJ) = TABLE(II,JJ)
 
335
          end do
 
336
          call tbrwrr(TID,II,4,COL,PAR,STAT)
 
337
        end do
 
338
 
 
339
        call tbtclo(TID,STAT)
 
340
 
 
341
C  delete old table
 
342
 
 
343
        call STFDEL('ttemp2.tbl',STAT)
 
344
 
 
345
C  finished
 
346
 
 
347
 9999        call stsepi
 
348
             end