← Back to branch summary

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

« back to all changes in this revision

Viewing changes to stdred/echelle/src/necconv.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/echelle/src/necconv.for
 
1
C @(#)necconv.for       19.1 (ES0-DMD) 02/25/03 14:20:22
 
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
C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
30
C.COPYRIGHT: Copyright (c) 1987 European Southern Observatory,
 
31
C                                         all rights reserved
 
32
C
 
33
C.VERSION: 1.0  ESO-FORTRAN Conversion, AA  21:03 - 3 DEC 1987
 
34
C
 
35
C.LANGUAGE: F77+ESOext
 
36
C
 
37
C.AUTHOR: D.PONZ
 
38
C
 
39
C.IDENTIFICATION
 
40
C
 
41
C  program ECHCONV2.FOR
 
42
C
 
43
C.PURPOSE
 
44
C
 
45
C  Execute the command
 
46
C  CONVERT/ECHELLE input output sampling_domain function parms option
 
47
C
 
48
C  input            - input image with lines corresponding to 
 
49
C                       echelle orders
 
50
C  output           - output image with lines corresponding to 
 
51
C                       echelle orders
 
52
C  sampling_domain  - start,step or filename defining start,step
 
53
C  function         - rebinning function
 
54
C  params           - define the transformation coeffs (def 0.,1.)
 
55
C  option           - rebinning method as PIX (def)
 
56
C                                         LIN
 
57
C                                         SPG
 
58
C.ALGORITHM
 
59
C
 
60
C  Linear rebinning
 
61
C
 
62
C.KEYWORDS
 
63
C
 
64
C  echelle, flux calibration, rebin
 
65
C
 
66
C.INPUT/OUTPUT
 
67
C
 
68
C  P1 - P5    contain input parameters
 
69
C
 
70
C.MODIFS
 
71
C 901024, M Peron , remove list_directed internal I/O
 
72
C 901217, P.Ballester, Skip bad orders of extracted standard star spectrum
 
73
C-----------------------------------------------------------
 
74
C
 
75
C
 
76
      PROGRAM ECHCV2
 
77
      IMPLICIT NONE
 
78
C
 
79
      INTEGER     NFPAR
 
80
      PARAMETER   (NFPAR=12)
 
81
      INTEGER     MADRID
 
82
      INTEGER     NPTOT1(100),NPIX1(2),NPIXR(2),IDUM(2)
 
83
      INTEGER     NPTOT2(100),NPIX2(2),NORDER(100)
 
84
      INTEGER*8   PNTR1, PNTR2
 
85
      INTEGER     IMNO, NA, KUN, KNUL, ISTAT
 
86
      INTEGER     INOP, I, J, NAXIS, NFPACT, NN, NDIM 
 
87
      INTEGER     II, NPMAX, STATUS,NFUNC, IMNO1, IMR
 
88
      DOUBLE PRECISION STEP1(2),START1(2),WSTR1(100)
 
89
      DOUBLE PRECISION STEP2(2),START2(2),WSTR2(100)
 
90
      REAL        CUTS(4),RMIN,RMAX,RDUM(2)
 
91
      DOUBLE PRECISION RPAR(2)
 
92
      DOUBLE PRECISION DPARM(NFPAR),FPARM(NFPAR)
 
93
      CHARACTER*8  IMAGE1,IMAGE2,METHOD,OPTION
 
94
      CHARACTER*64 IMAGER
 
95
      CHARACTER*72 IDENT
 
96
      CHARACTER*80 CUNIT
 
97
      CHARACTER*3  FUN(9)
 
98
C
 
99
      INCLUDE 'MID_INCLUDE:ST_DEF.INC'
 
100
      COMMON/VMR/MADRID(1)
 
101
      INCLUDE 'MID_INCLUDE:ST_DAT.INC'
 
102
C
 
103
      DATA FUN/'LIN','POL','INV','EXP','DEX','LOG','DLG','IPO','U01'/
 
104
C
 
105
C ... get into MIDAS
 
106
C
 
107
      CALL STSPRO('ECHCV2')
 
108
C
 
109
C ... get input parameters + default
 
110
C
 
111
      CALL STKRDC('P1',1,1,8,NA,IMAGE1,KUN,KNUL,ISTAT)
 
112
      CALL STKRDC('P2',1,1,8,NA,IMAGE2,KUN,KNUL,ISTAT)
 
113
      CALL STKRDC('P3',1,1,64,NA,IMAGER,KUN,KNUL,ISTAT)
 
114
      CALL STKRDC('P4',1,1,8,NA,METHOD,KUN,KNUL,ISTAT)
 
115
      CALL STKRDC('P6',1,1,8,NA,OPTION,KUN,KNUL,ISTAT)
 
116
      CALL STKRDD('INPUTD',1,NFPAR,NA,FPARM,KUN,KNUL,ISTAT) 
 
117
      INOP   = 4
 
118
      IF (OPTION(1:3).EQ.'PIX') INOP   = 1
 
119
      IF (OPTION(1:3).EQ.'LIN') INOP   = 2
 
120
      IF (OPTION(1:3).EQ.'SPG') INOP   = 3
 
121
C
 
122
C ... translate FUNCTION into function number
 
123
C
 
124
      NFUNC  = 0
 
125
      DO 10 I = 1,9
 
126
          IF (METHOD(1:3).EQ.FUN(I)) NFUNC  = I
 
127
   10 CONTINUE
 
128
      IF (NFUNC.EQ.0) THEN
 
129
          CALL STTPUT(' Specified function non-existent...',ISTAT)
 
130
          GO TO 50
 
131
      END IF
 
132
C
 
133
C ... search for number of active parameters and DOUBLE
 
134
C
 
135
      DO 20 J = NFPAR,1,-1
 
136
          NFPACT = J
 
137
          IF (FPARM(J).NE.0.0D0) GO TO 30
 
138
   20 CONTINUE
 
139
   30 DO 40 J = 1,NFPACT
 
140
          DPARM(J) = FPARM(J)
 
141
   40 CONTINUE
 
142
C
 
143
C
 
144
C ... read input image
 
145
C
 
146
      CALL STIGET(IMAGE1,D_R4_FORMAT,F_I_MODE,F_IMA_TYPE,
 
147
     .     2,NAXIS,NPIX1,START1,STEP1,IDENT,CUNIT,
 
148
     .     PNTR1,IMNO1,ISTAT)
 
149
      IF (NAXIS.GT.1 .AND. NPIX1(2).GT.1) THEN
 
150
CCCCCCCCCCCCC          CALL SXDFND(IMAGE1,'WSTART',TYPE,NO,NB,ISTAT)
 
151
C          IF (TYPE(1:1).EQ.' ') THEN
 
152
C              CALL STTPUT(' Wrong input image ',ISTAT)
 
153
C              GO TO 50
 
154
C          END IF
 
155
          CALL STDRDD(IMNO1,'WSTART',1,NPIX1(2),NN,WSTR1,
 
156
     .                KUN,KNUL,ISTAT)
 
157
          CALL STDRDI(IMNO1,'NPTOT',1,NPIX1(2),NN,NPTOT1,
 
158
     .                KUN,KNUL,ISTAT)
 
159
          NDIM   = 2
 
160
      ELSE
 
161
          NPIX1(2) = 1
 
162
          NDIM   = 1
 
163
          NPTOT1(1) = NPIX1(1)
 
164
          WSTR1(1) = START1(1)
 
165
      END IF
 
166
      START2(1) = START1(1)
 
167
      START2(2) = START1(2)
 
168
      STEP2(2) = STEP1(2)
 
169
C
 
170
C ... find start, step, npix for each order in the output image
 
171
C
 
172
      II     = INDEX('-+.0123456789',IMAGER(1:1))
 
173
      IF (II.EQ.0) THEN
 
174
          CALL STFOPN(IMAGER,D_R4_FORMAT,0,F_IMA_TYPE,IMR,ISTAT)
 
175
          CALL STDRDI(IMR,'NPIX',1,2,NN,NPIXR,KUN,KNUL,ISTAT)
 
176
          IF (NPIXR(2).NE.NPIX1(2)) THEN
 
177
              CALL STTPUT(' Error in reference image',ISTAT)
 
178
              GO TO 50
 
179
          END IF
 
180
          CALL STDRDD(IMR,'STEP',1,1,NN,STEP2,KUN,KNUL,ISTAT)
 
181
          CALL STDRDD(IMR,'WSTART',1,NPIX1(2),NN,WSTR2,
 
182
     .                KUN,KNUL,ISTAT)
 
183
          CALL STDRDI(IMR,'NPTOT',1,NPIX1(2),NN,NPTOT2,
 
184
     .                KUN,KNUL,ISTAT)
 
185
          CALL STDRDI(IMR,'NORDER',1,NPIX1(2),NN,NORDER,
 
186
     .                KUN,KNUL,ISTAT)
 
187
          NPMAX  = NPIXR(1)
 
188
      ELSE
 
189
          CALL GENCNV(IMAGER,4,2,IDUM,RDUM,RPAR,ISTAT)
 
190
          CALL OUTIMA(NPIX1(2),STEP1,WSTR1,NPTOT1,RPAR(1),RPAR(2),
 
191
     +            START2,STEP2,WSTR2,NPTOT2,NPMAX)
 
192
      END IF
 
193
C
 
194
C ... map output image
 
195
C
 
196
      NPIX2(1) = NPMAX
 
197
      NPIX2(2) = NPIX1(2)
 
198
      CALL STIPUT(IMAGE2,D_R4_FORMAT,F_O_MODE,F_IMA_TYPE,NDIM,
 
199
     .            NPIX2,START2,STEP2,IDENT,CUNIT,
 
200
     +            PNTR2,IMNO,STATUS)
 
201
C
 
202
C ... rebin image
 
203
C
 
204
      CALL APPREB(MADRID(PNTR1),NPIX1(1),NPIX1(2),STEP1,WSTR1,NPTOT1,
 
205
     +            MADRID(PNTR2),NPIX2(1),NPIX2(2),STEP2,WSTR2,NPTOT2,
 
206
     +            NFUNC,NFPAR,NFPACT,DPARM,INOP,RMIN,RMAX)
 
207
C
 
208
C ... write non standard descriptors
 
209
C
 
210
CCCCCCCCCCCCCCCCC      CALL SXDCOP(IMAGE1,IMAGE2,3,' ',ISTAT)
 
211
      CUTS(1) = RMIN
 
212
      CUTS(2) = RMAX
 
213
      CUTS(3) = RMIN
 
214
      CUTS(4) = RMAX
 
215
      CALL STDWRR(IMNO,'LHCUTS',CUTS,1,4,KUN,ISTAT)
 
216
      IF (NPIX1(2).GT.1) THEN
 
217
          CALL STDWRD(IMNO,'WSTART',WSTR2,1,NPIX2(2),KUN,ISTAT)
 
218
          CALL STDWRI(IMNO,'NPTOT',NPTOT2,1,NPIX2(2),KUN,ISTAT)
 
219
          CALL STDWRI(IMNO,'NORDER',NORDER,1,NPIX2(2),KUN,ISTAT)
 
220
      END IF
 
221
C
 
222
C ... end
 
223
C
 
224
   50 IF (ISTAT.NE.0) THEN
 
225
          CALL STTPUT(' Error in ECHCONV',STATUS)
 
226
      END IF
 
227
      CALL STSEPI
 
228
      END
 
229
 
 
230
      SUBROUTINE OUTIMA(NORDER,STEP1,WSTR1,NPTOT1,STARTR,STEPR,
 
231
     +                  START2,STEP2,WSTR2,NPTOT2,NPMAX)
 
232
C
 
233
C FIND SAMPLING DOMAIN IN THE OUTPUT IMAGE BASED ON THE
 
234
C REFERENCE START AND STEP VALUES
 
235
C
 
236
C INPUT ARGUMENTS
 
237
C NORDER  NO OF ORDERS
 
238
C STEP1   STEP IN INPUT IMAGE
 
239
C WSTR1   START IN INPUT IMAGE
 
240
C NPTOT1  NO OF PIXELS/ORDER
 
241
C STARTR  START USED AS REFERENCE
 
242
C STEPR   STEP USED AS REFERENCE
 
243
C OUTPUT ARGUMENTS
 
244
C START2  START IN OUTPUT IMAGE IF NORDER = 1
 
245
C STEP2   STEP IN OUTPUT IMAGE
 
246
C WSTR2   START IN OUTPUT IMAGE IF NORDER > 1
 
247
C NPTOT2  NO OF PIXELS/ORDER IF NORDER > 1
 
248
C NPMAX   MAX(NPTOT)
 
249
C
 
250
      IMPLICIT NONE
 
251
      DOUBLE PRECISION WSTR1(1),WSTR2(1)
 
252
      DOUBLE PRECISION STARTR, STEPR, W, STEP1, STEP2, START2
 
253
      INTEGER NPTOT1(1),NPTOT2(1)
 
254
      INTEGER I, NPMAX, NORDER, NPOS
 
255
C
 
256
C ... iterate for each order
 
257
C
 
258
      NPMAX  = 0
 
259
      DO 10 I = 1,NORDER
 
260
C
 
261
C ...     start, end for each order
 
262
C
 
263
          NPOS      = NINT((WSTR1(I)-STARTR)/STEPR)
 
264
          WSTR2(I)  = STARTR + (NPOS+1)*STEPR
 
265
          NPOS      = INT((WSTR1(I)+ (NPTOT1(I)-1)*STEP1-STARTR)/STEPR)
 
266
          W         = STARTR + (NPOS-1)*STEPR
 
267
          NPTOT2(I) = (W-WSTR2(I))/STEPR + 1
 
268
          NPMAX     = MAX(NPMAX,NPTOT2(I))
 
269
C   TYPE *,WSTR2(I),NPTOT2(I)
 
270
   10 CONTINUE
 
271
      STEP2  = STEPR
 
272
      IF (NORDER.EQ.1) START2 = WSTR2(1)
 
273
      RETURN
 
274
      END
 
275
 
 
276
      SUBROUTINE APPREB(X,NPIX11,NPIX12,STEP1,WSTR1,NPTOT1,
 
277
     +                  Y,NPIX21,NPIX22,STEP2,WSTR2,NPTOT2,
 
278
     +                  NFUNC,NFPAR,NFPACT,DPARM,INOP,RMIN,RMAX)
 
279
C
 
280
C REBIN EACH ORDER
 
281
C
 
282
      IMPLICIT NONE
 
283
      INTEGER NPIX11,NPIX12,NPIX21,NPIX22,NPTOT1(1),NPTOT2(1)
 
284
      REAL    X(NPIX11,NPIX12)
 
285
      REAL    Y(NPIX21,NPIX22)
 
286
      DOUBLE PRECISION WSTR1(1),WSTR2(1), STEP1, STEP2, WX
 
287
      REAL    RMIN,RMAX,RMIN1,RMAX1
 
288
      DOUBLE PRECISION DPARM(1), X0
 
289
      INTEGER   MADRID, NI, NI1, IS
 
290
      INTEGER*8 IP1, IP2, IP3
 
291
      INTEGER   NO, NO1
 
292
      INTEGER*8 JP1, JP2
 
293
      INTEGER   I, N1, N2, NINT, NP1, NFUNC,NFPAR
 
294
      INTEGER   NFPACT, INOP, J1, J2, J
 
295
C
 
296
      COMMON/VMR/MADRID(1)
 
297
C
 
298
C ALLOCATE WORK SPACE
 
299
C
 
300
      NI     = NPIX11
 
301
      NI1    = 8*NI
 
302
      CALL TDMGET(NI1,IP1,IS)
 
303
      CALL TDMGET(NI1,IP2,IS)
 
304
      CALL TDMGET(NI1,IP3,IS)
 
305
      NO     = NPIX21
 
306
      NO1    = 8*NO
 
307
      CALL TDMGET(NO1,JP1,IS)
 
308
      CALL TDMGET(NO1,JP2,IS)
 
309
      RMIN   = 0.
 
310
      RMAX   = 0.
 
311
      DO 20 I = 1,NPIX12
 
312
C
 
313
C   INPUT ARRAYS
 
314
C
 
315
          X0     = WSTR1(I)
 
316
          WX     = STEP1
 
317
          N1     = NPTOT1(I)
 
318
          CALL IMVAL3(N1,X0,WX,X(1,I),MADRID(IP1),MADRID(IP2),
 
319
     +                 MADRID(IP3),NP1)
 
320
C
 
321
C   OUTPUT ARRAYS
 
322
C
 
323
          N2     = NPTOT2(I)
 
324
          X0     = WSTR2(I)
 
325
          WX     = STEP2
 
326
          CALL IMVAL2(N2,X0,WX,MADRID(JP1),MADRID(JP2))
 
327
C
 
328
C         CALL INTERPOLATING (AND INTEGRATING) ROUTINES:
 
329
C
 
330
C  !!!! Watch out for VECI,VECD arrays in INTERP, 
 
331
C       dimension => NINT !!!!!!
 
332
C
 
333
          NINT   = 8
 
334
          CALL REBMET(NP1,MADRID(IP1),MADRID(IP3),MADRID(IP2),N2,
 
335
     +                 MADRID(JP1),MADRID(JP2),NFUNC,NFPAR,NFPACT,DPARM,
 
336
     +                 INOP,NINT,Y(1,I),RMIN1,RMAX1)
 
337
          RMIN   = AMIN1(RMIN1,RMIN)
 
338
          RMAX   = AMAX1(RMAX1,RMAX)
 
339
          J1     = N2 + 1
 
340
          J2     = NPIX21
 
341
          DO 10 J = J1,J2
 
342
              Y(J,I) = 0.
 
343
   10     CONTINUE
 
344
   20 CONTINUE
 
345
C
 
346
C  FREE MEMORY
 
347
C
 
348
      CALL TDMFRE(NI1,IP1,IS)
 
349
      CALL TDMFRE(NI1,IP2,IS)
 
350
      CALL TDMFRE(NI1,IP3,IS)
 
351
      CALL TDMFRE(NO1,JP1,IS)
 
352
      CALL TDMFRE(NO1,JP2,IS)
 
353
      RETURN
 
354
      END
 
355
 
 
356
      SUBROUTINE IMVAL2(NP,XS,WS,X,W)
 
357
C +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
358
C
 
359
C  Subroutine IMVAL2
 
360
C
 
361
C
 
362
C  Fill arrays X,W for image world coordinates and pixel size data
 
363
C
 
364
      IMPLICIT NONE
 
365
      INTEGER  I, NP
 
366
      DOUBLE PRECISION X(NP),W(NP),WSS,XSS
 
367
      DOUBLE PRECISION WS, XS
 
368
C
 
369
      WSS    = WS
 
370
      XSS    = XS
 
371
      DO 10 I = 1,NP
 
372
          X(I)   = XSS + WSS* (I-1)
 
373
          W(I)   = WSS
 
374
   10 CONTINUE
 
375
      RETURN
 
376
      END
 
377
 
 
378
      SUBROUTINE IMVAL3(NP,XS,WS,YY,X,W,Y,NP1)
 
379
C +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
380
C
 
381
C  Subroutine IMVAL3
 
382
C
 
383
C
 
384
C  Fill arrays X,W for image world coordinates and pixel size data
 
385
C  Convert YY(real*4) into Y(real*8)
 
386
C
 
387
      IMPLICIT NONE
 
388
      INTEGER  NP
 
389
      DOUBLE PRECISION X(NP),W(NP),Y(NP),WSS,XSS
 
390
      DOUBLE PRECISION WS, XS
 
391
      REAL YY(NP)
 
392
      INTEGER NP1, I
 
393
C
 
394
      NP1    = 0
 
395
      WSS    = WS
 
396
      XSS    = XS
 
397
      DO 10 I = 1,NP
 
398
           IF (YY(I).GE.0.0) THEN         ! Modif. 901217 P.B (old .GT.)
 
399
              NP1    = NP1 + 1
 
400
              X(NP1) = XSS + WSS* (I-1)
 
401
              W(NP1) = WSS
 
402
              IF (YY(I).GT.1.E30) THEN
 
403
                 Y(NP1) = 0.D0
 
404
              ELSE 
 
405
                 Y(NP1) = YY(I)
 
406
              ENDIF
 
407
           END IF
 
408
   10 CONTINUE
 
409
      RETURN
 
410
      END