C @(#)lincol.for 17.1.1.1 (ES0-DMD) 01/25/02 17:40:35 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 Correspondence 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 SUBROUTINE LINCOL(A,NPIX,SIZE,B) C C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C C.IDENTIFICATION C subroutine LINCOL version 1.0 870109 C K. Banse ESO - Garching C C.KEYWORDS C lines to columns C C.PURPOSE C transposition of a matrix around minor diagonal C that is, the lines of an image are changed to columns C C.ALGORITHM C the input + result image must be completely mapped into memory C in order to minimize paging, this transposition is done in pieces C C.INPUT/OUTPUT C call as LINCOL(A,NPIX,SIZE,B) C C input par: C A: R*4 array input image C NPIX: I*4 array no. of pixels per line, no. of lines in A C SIZE: I*4 array size of chunk we work with in x,y C C output par: C B: R*4 array output image C C ------------------------------------------------------------------- C IMPLICIT NONE C INTEGER NPIX(2),SIZE(2) INTEGER NX,NY,N,K,M INTEGER XDIM,YDIM INTEGER XNO1,YNO1,NKK,YMULT,XMULT INTEGER OFFA,OFFB,AUX1 C REAL A(*),B(*) C C XDIM = NPIX(1) YDIM = NPIX(2) AUX1 = SIZE(2) - 1 YMULT = SIZE(1) * YDIM XMULT = SIZE(2) * XDIM C C loop through complete frame in chunks C NKK = 0 XNO1 = SIZE(1) - 1 DO 800,NX=1,XDIM,SIZE(1) YNO1 = AUX1 N = XDIM - NX IF (XNO1.GT.N) XNO1 = N C M = NX DO 600,NY=1,YDIM,SIZE(2) N = YDIM - NY IF (YNO1.GT.N) YNO1 = N C C was N = (NX-1)*YDIM + NY C before M = (NY-1)*XDIM + NX C N = NKK + NY C C loop through pixels in input image DO 400,K=M,M+XNO1 OFFA = K C DO 200,OFFB=N,N+YNO1 B(OFFB) = A(OFFA) OFFA = OFFA + XDIM 200 CONTINUE C N = N + YDIM 400 CONTINUE C M = M + XMULT 600 CONTINUE C NKK = NKK + YMULT 800 CONTINUE C RETURN END SUBROUTINE ROTM90(INTPOL,A,B,IMNI,IMNO,NPIX,STRIP) C C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C C.IDENTIFICATION C subroutine ROTM90 version 1.0 920311 C K. Banse ESO - Garching C 1.1 930902 C C.KEYWORDS C rotation, transposition C C.PURPOSE C rotation by multiples of 90. degrees of 2-D frames C C.ALGORITHM C chunks of input + result frame are processed, so this routine works for C any size frames C C.INPUT/OUTPUT C call as ROTM90(INTPOL,A,B,IMNI,IMNO,NPIX,STRIP) C C input par: C INTPOL: I*4 = -1, for transposing lines to columns C (around minor diagonal) C 1, for 90 degrees counter clockwise C 2, for 180 degrees counter clockwise C 3, for 270 degrees counter clockwise C 4, for 360 degrees counter clockwise C A: R*4 array working buffer for input C B: R*4 array working buffer for output C IMNI: I*4 file no. of input image C IMNO: I*4 file no. of output image C NPIX: I*4 array x-, y-dim of input frame C STRIP: I*4 array chunk size for input + output frame C C.VERSIONS C 1.1 fix 180 degree rotation for non-square files C C ------------------------------------------------------------------- C IMPLICIT NONE C INTEGER INTPOL INTEGER IMNI,IMNO,NPIX(2),STRIP(2) INTEGER INDXA,INDXB,KY,KY1,KYY,ISTRIP,OSTRIP INTEGER IOFFA,IOFFB,NN,DATIN,DATOUT INTEGER HYI,HYO,INOFF,OUTOFF INTEGER KOUT,IOUT,KSIZE,IAV,STAT C REAL A(*) !size = NPIX(1) * STRIP(1) (input) REAL B(*) !size = NPIX(2) * STRIP(2) (output) C C init ISTRIP = STRIP(1) OSTRIP = STRIP(2) IOUT = 1 !x-start for input C C branch on INTPOL flag IF (INTPOL.LE.0) GOTO 8000 GOTO (100,4000,6000,7000),INTPOL C C handle 90. degrees C ------------------ C 100 OUTOFF = 0 !y-offset of output C 1100 HYI = NPIX(2) !high y-index of input INOFF = HYI - ISTRIP !y-offset of input C HYO = OUTOFF + OSTRIP !high y-index of output IF (HYO.GT.NPIX(1)) HYO = NPIX(1) KYY = HYO - OUTOFF !no. of output lines KOUT = 1 !x-start for output C 1200 IF (INOFF.LT.0) INOFF = 0 DATIN = (INOFF*NPIX(1)) + 1 !offset for input read KY = HYI - INOFF !no. of input lines KY1 = KY - 1 KSIZE = KY * NPIX(1) CALL STFGET(IMNI,DATIN,KSIZE,IAV,A,STAT) C IOFFA = KY1*NPIX(1) + IOUT IOFFB = KOUT C DO 2000, NN=1,KYY INDXA = IOFFA C DO 1500, INDXB=IOFFB,IOFFB+KY1 B(INDXB) = A(INDXA) INDXA = INDXA - NPIX(1) 1500 CONTINUE C IOFFA = IOFFA + 1 IOFFB = IOFFB + NPIX(2) 2000 CONTINUE C C we've done a KY * KYY piece - check further progress HYI = HYI - KY IF (HYI.GT.0) THEN KOUT = KOUT + KY INOFF = HYI - ISTRIP GOTO 1200 !get next y-strip of input frame ENDIF C C so HYI must be 0 - write output strip to disk C KSIZE = KYY * NPIX(2) DATOUT = (OUTOFF*NPIX(2)) + 1 CALL STFPUT(IMNO,DATOUT,KSIZE,B,STAT) OUTOFF = OUTOFF + KYY IF (OUTOFF.LT.NPIX(1)) THEN IOUT = IOUT + KYY GOTO 1100 !set up things for next output strip ENDIF C RETURN C C handle 180. degrees C ------------------- C 4000 INOFF = 0 !y-offset for in/output HYO = NPIX(2) KYY = NPIX(1) - 1 C 4100 HYI = INOFF + ISTRIP !high y-index of in/output IF (HYI.GT.NPIX(2)) HYI = NPIX(2) KY = HYI - INOFF KY1 = KY - 1 KSIZE = KY * NPIX(1) DATIN = (INOFF*NPIX(1)) + 1 !disk pointer for in/output CALL STFGET(IMNI,DATIN,KSIZE,IAV,A,STAT) OUTOFF = HYO - KY C IOFFA = 1 IOFFB = (KY1*NPIX(1)) + 1 C DO 4600, NN=1,KY INDXA = IOFFA + KYY C DO 4500, INDXB=IOFFB,IOFFB+KYY B(INDXB) = A(INDXA) INDXA = INDXA - 1 4500 CONTINUE C IOFFA = IOFFA + NPIX(1) IOFFB = IOFFB - NPIX(1) 4600 CONTINUE C C we've done a KY strip DATOUT = (OUTOFF*NPIX(1)) + 1 CALL STFPUT(IMNO,DATOUT,KSIZE,B,STAT) INOFF = HYI IF (INOFF.LT.NPIX(2)) THEN HYO = HYO - KY GOTO 4100 ENDIF C RETURN C C handle 270. degrees C ------------------- C 6000 HYO = NPIX(1) !high y-index of output C 6100 INOFF = 0 !y-offset of input OUTOFF = HYO - OSTRIP !y-offset of output IF (OUTOFF.LT.0) OUTOFF = 0 KYY = HYO - OUTOFF !no. of output lines KOUT = 1 !x-start for output C 6200 HYI = INOFF + ISTRIP !high y-index of input IF (HYI.GT.NPIX(2)) HYI = NPIX(2) DATIN = (INOFF*NPIX(1)) + 1 !offset for input read KY = HYI - INOFF !no. of input lines KY1 = KY - 1 KSIZE = KY * NPIX(1) CALL STFGET(IMNI,DATIN,KSIZE,IAV,A,STAT) C IOFFA = IOUT IOFFB = (KYY-1)*NPIX(2) + KOUT C DO 6800, NN=1,KYY INDXA = IOFFA C DO 6500, INDXB=IOFFB,IOFFB+KY1 B(INDXB) = A(INDXA) INDXA = INDXA + NPIX(1) 6500 CONTINUE C IOFFA = IOFFA + 1 IOFFB = IOFFB - NPIX(2) 6800 CONTINUE C C we've done a KY * KYY piece - check further progress INOFF = INOFF + KY IF (INOFF.LT.NPIX(2)) THEN KOUT = KOUT + KY GOTO 6200 !get next y-strip of input frame ENDIF C C so INOFF must be >= NPIX(2) - write output strip to disk C KSIZE = KYY * NPIX(2) DATOUT = (OUTOFF*NPIX(2)) + 1 CALL STFPUT(IMNO,DATOUT,KSIZE,B,STAT) HYO = HYO - KYY IF (HYO.GT.0) THEN IOUT = IOUT + KYY GOTO 6100 !set up things for next output strip ENDIF C RETURN C C C handle 360. degrees (= simple copy...) C -------------------------------------- C 7000 CALL STTPUT('No rotation - file is just copied...',STAT) CALL COPYFX(A,IMNI,IMNO,NPIX,STRIP(1)) RETURN C C handle transposition of lines + columns C --------------------------------------- C 8000 OUTOFF = 0 !y-offset of output C 8100 INOFF = 0 !y-offset of input HYO = OUTOFF + OSTRIP !high y-index of output IF (HYO.GT.NPIX(1)) HYO = NPIX(1) KYY = HYO - OUTOFF !no. of output lines KOUT = 1 !x-start for output C 8200 HYI = INOFF + ISTRIP IF (HYI.GT.NPIX(2)) HYI = NPIX(2) DATIN = (INOFF*NPIX(1)) + 1 !offset for input read KY = HYI - INOFF !no. of input lines KY1 = KY - 1 KSIZE = KY * NPIX(1) CALL STFGET(IMNI,DATIN,KSIZE,IAV,A,STAT) C IOFFA = IOUT IOFFB = KOUT C DO 8800, NN=1,KYY INDXA = IOFFA C DO 8500, INDXB=IOFFB,IOFFB+KY1 B(INDXB) = A(INDXA) INDXA = INDXA + NPIX(1) 8500 CONTINUE C IOFFA = IOFFA + 1 IOFFB = IOFFB + NPIX(2) 8800 CONTINUE C C we've done a KY * KYY piece - check further progress INOFF = INOFF + KY IF (INOFF.LT.NPIX(2)) THEN KOUT = KOUT + KY GOTO 8200 !get next y-strip of input frame ENDIF C C so INOFF must be >= NPIX(2) - write output strip to disk C KSIZE = KYY * NPIX(2) DATOUT = (OUTOFF*NPIX(2)) + 1 CALL STFPUT(IMNO,DATOUT,KSIZE,B,STAT) OUTOFF = OUTOFF + KYY IF (OUTOFF.LT.NPIX(1)) THEN IOUT = IOUT + KYY GOTO 8100 !set up things for next output strip ENDIF C RETURN END