C @(#)rebirbr.for 17.1.1.1 (ES0-DMD) 01/25/02 17:53:58 C=========================================================================== C Copyright (C) 1995 European Southern Observatory (ESO) C C This program is free software; you can redistribute it and/or C modify it under the terms of the GNU General Public License as C published by the Free Software Foundation; either version 2 of C the License, or (at your option) any later version. C C This program is distributed in the hope that it will be useful, C but WITHOUT ANY WARRANTY; without even the implied warranty of C MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the C GNU General Public License for more details. C C You should have received a copy of the GNU General Public C License along with this program; if not, write to the Free C Software Foundation, Inc., 675 Massachusetss Ave, Cambridge, C MA 02139, USA. C C Corresponding concerning ESO-MIDAS should be addressed as follows: C Internet e-mail: midas@eso.org C Postal address: European Southern Observatory C Data Management Division C Karl-Schwarzschild-Strasse 2 C D 85748 Garching bei Muenchen C GERMANY C=========================================================================== C C @(#)rebirbr.for 17.1.1.1 (ESO) 01/25/02 17:53:58 C +++++++++++++++++++++++++++++++++++++++++++++++++++++++ C .COPYRIGHT (C) 1993 European Southern Observatory C .IDENT .prg C .AUTHORS Pascal Ballester (ESO/Garching) C Cristian Levin (ESO/La Silla) C .KEYWORDS Spectroscopy, Long-Slit C .PURPOSE C .VERSION 1.0 Package Creation 17-MAR-1993 C ------------------------------------------------------- C @(#)rebirbr.for 1.1 (ESO-IPG) 2/17/93 16:46:20 C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C.COPYRIGHT: Copyright (c) 1987 European Southern Observatory, C all rights reserved C C.VERSION: 1.0 ESO-FORTRAN Conversion, AA 8:58 - 9 DEC 1987 C C.LANGUAGE: F77+ESOext C C.IDENTIFICATION C C SUBROUTINES TDRBII MODIFIED TO REBIN A 2D FRAME ROW BY ROW C (M.Pierre Nov. 1989) C.KEYWORDS C C IMAGES, NONLINEAR REBIN, INTERPOLATION, INTEGRATION C C.PURPOSE C C NONLINEAR REBIN IMAGE TO IMAGE IN 2 DIMENSIONS C C.ALGORITHM C C Rebin input data YI sampled at XI with binwidth WI to C output data YO sampled at XO with binwidth WO. C C Restriction: Input and output independent variables are C assumed to be monotonic increasing or decreasing C C Several FUNCTal relations between input and output independent C variables are foreseen, additional ones may be added in dummy call C U01 as follows: C C SUBROUTINE REB_U01 (DV,P,NP,NPC) C REAL*8 DV,P(NP) C DV = ........... C RETURN C END C C where DV is input-output value, P(NP) the parameter array and C NPC the nuber of active parameters in p(np). C C According to the definition of IMAGES, pixels have constant C bin width (STEP) [WI;WO], are not overlapping and have filling C factors of exactly one. World coordinates [XI;XO] and flux values C [YI;YO] are pertinent C to the center of a pixel. In flux conservation mode (the default) C the flux is interpreted as the area under a smooth curve between C the boundaries of the pixel (i.e. a histogram representation of C a smooth curve). C C Detailed flux conservation obtained by REAL*8 spline interpolation C and integration of input data over the projected C bin WI for fractions of input pixels and summation of flux in C entire input pixels covered within the valid domain of input data C i.e. from x(1)-0.5*step to x(n)+0.5*step. Several options available C to handle undefined input at the extremes and for NULLVALUE assigned C parts of input array. C C Optimized for long arrays by search for valid range, optimum C direction or DO LOOPS and updating of loop boundaries - therefore C overdoing it for small arrays. C C.INPUT C C Keys: IN_A/C/1/8 input data array C IN_B/C/1/8 reference frame or dim of output C CFUNC/C/1/8 FUNCTal relation between indep. C INPUTD/D/1/12 parameters or FUNCT C COPT/C/1/8 options C C.OUTPUT C C Key: OUT_A/C/1/8 output data array C C.MODIFS C M Peron sep 90 : remove stupid limit to the size of the output frame!!!!! C-------------------------------------------------------------------------- C PROGRAMM REBIRBR IMPLICIT NONE INTEGER NFPAR PARAMETER (NFPAR = 12) C INTEGER MADRID INTEGER I,IACT INTEGER*8 IP1,IP2,IP3,IPNTR INTEGER IS,ISTAT INTEGER*8 JP1,JP2,JP3,JPNTR INTEGER NAXISI,NAXISO,KUN,KNUL INTEGER NFPACT,NFUNC INTEGER NI,NINT,NN,NO,NBY,NBO,NCOL,NTOT INTEGER NPIXI(3),NPIXO(3),INOP,IMNO,IMNI,IREF INTEGER JJ,TID,NROW,NSC,NAC,NAR,ICOL(12),STATUS REAL RMIN, RMAX, CUTS(4) DOUBLE PRECISION STEPI(3),STARTI(3) REAL X0O,WXO,X0I,WXI DOUBLE PRECISION DPARM(NFPAR) DOUBLE PRECISION DSTART(3), DSTEP(3),YPOS CHARACTER*3 FUN(9) CHARACTER*60 INPIMA,OUTIMA,REFIMA,TBCOEFF CHARACTER*80 FUNCT,OPTION CHARACTER*72 IDENTO,IDENTI CHARACTER*80 CUNITO,CUNITI LOGICAL NULL(12) C INCLUDE 'MID_INCLUDE:TABLES.INC' COMMON /VMR/MADRID(1) INCLUDE 'MID_INCLUDE:TABLED.INC' C DATA FUN/'LIN','POL','INV','EXP','DEX','LOG','DLG','IPO','U01'/ C C..... GET INTO MIDAS CALL STSPRO('REBLL') C ... get input parameters C CALL STKRDC('OUT_A',1,1,60,IACT,OUTIMA,KUN,KNUL,ISTAT) CALL STKRDC('IN_A', 1,1,60,IACT,INPIMA,KUN,KNUL,ISTAT) CALL STKRDC('IN_B', 1,1,60,IACT,REFIMA,KUN,KNUL,ISTAT) CALL STKRDC('CFUNC',1,1,8,IACT,FUNCT,KUN,KNUL,ISTAT) CALL STKRDC('IN_C', 1,1,60,IACT,TBCOEFF,KUN,KNUL,ISTAT) CALL STKRDC('COPT',1,1,8,IACT,OPTION,KUN,KNUL,ISTAT) C C ... det. interpolation/integration options C ! Not implemented INOP = 4 C IF (OPTION(1:3).EQ.'PIX') INOP = 1 C IF (OPTION(1:3).EQ.'LIN') INOP = 2 C IF (OPTION(1:3).EQ.'SPG') INOP = 3 C .......................................... the actual values are: CALL FORUPC(OPTION,OPTION) IF (OPTION(1:3).EQ.'PIX') INOP = 2 IF (OPTION(1:3).EQ.'LIN') INOP = 3 IF (OPTION(1:3).EQ.'SPG') INOP = 1 C .......................................... Bitte, bitte C C ... translate FUNCT into FUNCT number C CALL FORUPC(FUNCT,FUNCT) NFUNC = 0 DO 10 I = 1,9 IF (FUNCT(1:3).EQ.FUN(I)) NFUNC = I 10 CONTINUE IF (NFUNC.EQ.0) THEN CALL STTPUT(' Specified FUNCT non-existent...',ISTAT) GO TO 60 END IF C C..... init the input table of dispersion coefficients CALL TBTOPN(TBCOEFF,F_U_MODE,TID,ISTAT) CALL TBIGET(TID,NCOL,NROW,NSC,NAC,NAR,ISTAT) NFPACT = NCOL-1 DO 40 I = 1,NFPACT ICOL(I) = I+1 40 CONTINUE C C ... get parameters form reference image C CALL STFOPN(REFIMA,D_R4_FORMAT,0,F_IMA_TYPE,IREF,ISTAT) CALL STDRDI(IREF,'NAXIS',1,1,NN,NAXISO,KUN,KNUL,ISTAT) CALL STDRDI(IREF,'NPIX',1,NAXISO,NN,NPIXO,KUN,KNUL,ISTAT) CALL STDRDD(IREF,'START',1,NAXISO,NN,DSTART,KUN,KNUL,ISTAT) CALL STDRDD(IREF,'STEP',1,NAXISO,NN,DSTEP,KUN,KNUL,ISTAT) CALL STDRDC(IREF,'IDENT',1,1,72,NN,IDENTO,KUN,KNUL,ISTAT) CALL STDRDC(IREF,'CUNIT',1,1,80,NN,CUNITO,KUN,KNUL,ISTAT) CALL STFCLO(IREF,ISTAT) NTOT = NPIXO(1)*NPIXO(2) C C ... get input image C CALL STIGET(INPIMA,D_R4_FORMAT,F_I_MODE,F_IMA_TYPE, .2,NAXISI,NPIXI,STARTI,STEPI,IDENTI,CUNITI,IPNTR,IMNI,ISTAT) C..... check dimensions IF(NPIXI(2).NE.NROW) THEN CALL STTPUT('Input table and image have not + the same number of rows',ISTAT) GO TO 50 ELSE IF(NPIXI(2).NE.NPIXO(2)) THEN CALL STTPUT('Input and output images have not + the same number of rows',ISTAT) GO TO 50 ENDIF ENDIF CALL TBRRDD(TID,1,1,1,YPOS,NULL,STATUS) IF (ABS(YPOS-DSTART(2)).GT.0.001) THEN CALL STTPUT(' WARNING:',ISTAT) CALL STTPUT('Y_position recorded in the 1st line + of input table is not the same',ISTAT) CALL STTPUT('as Y_start of the output image',ISTAT) ENDIF C ... map output image C CALL STIPUT(OUTIMA,D_R4_FORMAT,F_O_MODE,F_IMA_TYPE, .NAXISO,NPIXO,DSTART,DSTEP,IDENTO,CUNITO,JPNTR,IMNO,ISTAT) C C ... fill the arrays XI,WI and XO,WO C NI = NPIXI(1) NBY = 8*NI CALL TDMGET(NBY,IP1,IS) CALL TDMGET(NBY,IP2,IS) CALL TDMGET(NBY,IP3,IS) X0O = STARTI(1) WXO = STEPI(1) C NO = NPIXO(1) NBO = 8*NO CALL TDMGET(NBO,JP1,IS) CALL TDMGET(NBO,JP2,IS) CALL TDMGET(NBO,JP3,IS) X0I = DSTART(1) WXI = DSTEP(1) C..... DO LOOP ON THE IMAGE ROWS DO 100 JJ = 1,NROW CALL IMVAL3(NI,X0O,WXO,MADRID(IPNTR),MADRID(IP1),MADRID(IP2), + MADRID(IP3),JJ,NTOT) CALL IMVAL2(NO,X0I,WXI,MADRID(JP1),MADRID(JP2)) C C...... Read dispersion coefficients for row JJ CALL TBRRDD(TID,JJ,NFPACT,ICOL,DPARM,NULL,STATUS) WRITE(6,*) 'ROW',JJ C TYPE *,'COEFF',DPARM(1),DPARM(2),DPARM(3),DPARM(4),DPARM(5) C ... Now, at last, we call the interpolating (and integrating) routines: C C !!!! Watch out for VECI,VECD arrays in INTERP, dimension => NINT !!!!!! C NINT = 8 CALL REBMET(NI,MADRID(IP1),MADRID(IP3),MADRID(IP2),NO, + MADRID(JP1),MADRID(JP2),NFUNC,NFPAR,NFPACT,DPARM, + INOP,NINT,MADRID(JP3),RMIN,RMAX) C CALL IMVAL5(NO,JJ,MADRID(JP3),MADRID(JPNTR),NTOT) C 100 CONTINUE CALL TDMFRE(NBY,IP1,IS) CALL TDMFRE(NBY,IP2,IS) CALL TDMFRE(NBY,IP3,IS) CALL TDMFRE(NBO,JP1,IS) CALL TDMFRE(NBO,JP2,IS) C ... write cuts C CUTS(1) = RMIN CUTS(2) = RMAX CUTS(3) = RMIN CUTS(4) = RMAX CALL STDWRR(IMNO,'LHCUTS',CUTS,1,4,KUN,ISTAT) C C ... free memory C C ... end C CALL STFCLO(IMNO,ISTAT) 50 CALL STFCLO(IMNI,ISTAT) CALL TBTCLO(TID,ISTAT) 60 CALL STSEPI END C************************************************************** SUBROUTINE IMVAL2(NP,XS,WS,X,W) C C Fill arrays X,W for image world coordinates and pixel size data C IMPLICIT NONE INTEGER I, NP REAL XS, WS DOUBLE PRECISION X(NP),W(NP),WSS,XSS C WSS = WS XSS = XS DO 10 I = 1,NP X(I) = XSS + WSS* (I-1) W(I) = WSS 10 CONTINUE RETURN END C------------------------------------------------------------------ SUBROUTINE IMVAL3(NP,XS,WS,YY,X,W,Y,NL,NTOT) C C Fill arrays X,W for image world coordinates and pixel size data C Convert YY(real*4) into Y(real*8) C IMPLICIT NONE INTEGER NP, I,NTOT INTEGER NL,IPT DOUBLE PRECISION X(NP),W(NP),Y(NP),WSS,XSS REAL XS, WS REAL YY(NTOT) C WSS = WS XSS = XS DO 10 I = 1,NP IPT = (NL-1)*NP+I X(I) = XSS + WSS* (I-1) W(I) = WSS Y(I) = YY(IPT) 10 CONTINUE RETURN END C---------------------------------------------------------- SUBROUTINE IMVAL5(N,NL,X,Y,NTOT) C FILLS THE BUFFER INTEGER I,N,NL,NTOT,NPT REAL*4 Y(NTOT),X(N) DO 10 I = 1,N NPT = (NL-1)*N+I Y(NPT) = X(I) 10 CONTINUE RETURN END