SUBROUTINE ZVIGF * * Module number: * * Module name: ZVIGF * * Keyphrase: * ---------- * Fit smooth curve to vignetting values * * Description: * ------------ * This routine fits a least sqares cubic spline using user * specified node positions to the raw input vignetting function. * The function is then normalized to equal 1.0 at sample position * 280.0 * * FORTRAN name: zvigf.for * * Keywords of accessed files and tables: * -------------------------------------- * Name I/O Description / Comments * * input I input raw vignetting template * mask I input mask template * output1 O fit to input data * output2 O output vignetting function * outtable O output node table * intable I optional input node table * nodes I input number of nodes * niter I number of least sqares iterations * template I output file template * * Subroutines Called: * ------------------- * CDBS: * zspfit, zlintp, ztplat * SDAS: * uclgs* , umsput * uttinn, utppti, utcdef, utrpt*, utcpt*, uthad*, uttclo, utccre * uttopn, utpgti, utcfnd, utrgt*, utcgt*, uthgt*, uttclo * uimotp, uimxtp, uimctp, uimopn, uimgid, uimclo, uhdgs* uimclo * * History: * -------- * Version Date Author Description * 1 Oct 88 D. Lindler Designed and coded * 1.1 Jan 92 S. Hulbert New grating values *------------------------------------------------------------------------------- C C FILE I/O ACCESS MODES C INTEGER RDONLY PARAMETER (RDONLY = 1) INTEGER RDWRIT PARAMETER (RDWRIT = 2) C C CODES FOR DATA TYPES C INTEGER TYREAL PARAMETER (TYREAL = 6) INTEGER TYDOUB PARAMETER (TYDOUB = 7) C C UMSPUT DESTINATIONS -- CB, DAO, 4-SEP-87 C INTEGER STDOUT PARAMETER (STDOUT = 1) INTEGER STDERR PARAMETER (STDERR = 2) C C UHDAS HEADER PARM TYPES -- CB, DAO, 5-SEP-87 C INTEGER GENHDR PARAMETER (GENHDR = 0) INTEGER IMSPEC PARAMETER (IMSPEC = 1) C C THIS SECTION IS FOR PARAMETERS RELEVANT TO TABLE I/O. C C THESE MAY BE SET BY UTPPTI AND/OR READ BY UTPGTI: C C C THESE MAY BE READ BY UTPGTI BUT MAY NOT BE SET: C C NUMBER OF ROWS WRITTEN TO INTEGER TBNROW PARAMETER (TBNROW = 21) C END IRAF77.INC C C ERROR PROCESSING PARAMETERS C INTEGER STATUS,ISTAT,ISTATS(20) CHARACTER*130 CONTXT C C KEYWORD PARAMETERS C INTEGER NODES C --->NUMBER OF EQUALLY SPACED NODES INTEGER NITER C --->NUMBER OF ITERATIONS OF FIT ROUTINE C C INPUT FILE I/O C CHARACTER*130 INPUT,MASKT C --->TEMPLATES CHARACTER*64 NAME,MNAME,TEMPLT,OUTPT1,NAME1,OUTPT2 C --->FILE NAMES INTEGER IDTEMP,IDTMP2 C --->TEMPLATE IDS INTEGER IDIN,IDIN2 C --->FILE IDS INTEGER NAXIS,DTYPE,DIMEN(8),NS C --->DATA DESCRIPTIONS DOUBLE PRECISION DATA(2000),MASK(2000) C --->DATA BUFFERS LOGICAL FIRST C --->FIRST FILE IN TEMPLATE C C HEADER PARAMETERS C DOUBLE PRECISION SAMPLE,DELTAS,LINE CHARACTER*5 GRAT INTEGER CARPOS CHARACTER*8 HEADER(5) C C NODE TABLES C CHARACTER*64 OUTNOD,INNODE CHARACTER*10 COLNAM(2) DOUBLE PRECISION SNODES(30),VNODES(30) CHARACTER*8 CUNITS(2),CFORM(2) INTEGER IDOUT,CTYPE(2),COLIDS(2) LOGICAL NULLS(30) C C FIT ROUTINE PARAMETERS C DOUBLE PRECISION FIT(2000),TOTAL(2000),SAMP(2000),DIF,VIG(561) DOUBLE PRECISION VNORM INTEGER NADDS(2000),NOUT,I C C HEADER PARAMETER DATA C DATA HEADER/'GRATING','LINE','SAMPLE','DELTAS', * 'CARPOS'/ C C NODE TABLE DATA C DATA COLNAM/'SAMPLE','VALUE'/ DATA CUNITS/2*' '/ DATA CFORM/2*' '/ DATA CTYPE/2*TYDOUB/ C C-----------------------------------------------------------------INPUT CL C C READ INPUT CL PARAMETERS C CALL UCLGST('input',INPUT,ISTATS(1)) CALL UCLGST('mask',MASKT,ISTATS(2)) CALL UCLGST('output1',OUTPT1,ISTATS(3)) CALL UCLGST('output2',OUTPT2,ISTATS(3)) CALL UCLGST('outtable',OUTNOD,ISTATS(4)) CALL UCLGSI('nodes',NODES,ISTATS(5)) CALL UCLGSI('niter',NITER,ISTATS(6)) CALL UCLGST('intable',INNODE,ISTATS(7)) CALL UCLGST('template',TEMPLT,ISTATS(8)) DO 10 I=1,8 IF(ISTATS(I).NE.0)THEN CONTXT='ERROR getting CL parameter value' GO TO 999 ENDIF 10 CONTINUE C C INITIALIZATION C DO 15 I=1,2000 NADDS(I)=0 TOTAL(I)=0.0 15 MASK(I)=1.0 C----------------------------------------------------------INPUT OBSERVATIONS C C OPEN TEMPLATES C CALL UIMOTP(INPUT,IDTEMP,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error opening input filename template '//INPUT GO TO 999 ENDIF IF(MASKT.NE.' ')THEN CALL UIMOTP(MASKT,IDTMP2,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error opening mask file template '//MASKT GO TO 999 ENDIF ENDIF C C GET NEXT FILE NAME C FIRST=.TRUE. 20 CALL UIMXTP(IDTEMP,NAME,ISTAT) IF((ISTAT.LT.0).AND.FIRST)THEN CONTXT='NO DATA FOUND IN INPUT TEMPLATE' GO TO 999 ENDIF IF (ISTAT.LT.0) GO TO 100 IF(ISTAT.NE.0)THEN CONTXT='Error getting filename from template '//INPUT GO TO 999 ENDIF C C GET NEXT MASK FILE NAME C IF(MASKT.NE.' ')THEN CALL UIMXTP(IDTMP2,MNAME,ISTAT) IF(ISTAT.LT.0)THEN CONTXT='Insufficient Mask files supplied' GO TO 999 ENDIF IF(ISTAT.NE.0)THEN CONTXT='Error getting filename from template '//MASKT GO TO 999 ENDIF ENDIF C C OPEN INPUT FILES C CALL UIMOPN(NAME,RDONLY,IDIN,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error opening input file '//NAME GO TO 999 ENDIF IF(MASKT.NE.' ')THEN CALL UIMOPN(MNAME,RDONLY,IDIN2,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error opening mask file '//MNAME GO TO 999 ENDIF ENDIF C C READ IMAGE INFO C CALL UIMGID(IDIN,DTYPE,NAXIS,DIMEN,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error reading input file '//NAME GO TO 999 ENDIF C C CHECK FOR VALID DATA C IF((NAXIS.NE.1).OR.(DIMEN(1).GT.2000))THEN CONTXT='Input data must vector with max of 2000 points' GO TO 999 ENDIF C C IF FIRST FILE THEN GET HEADER INFO C IF(FIRST)THEN FIRST=.FALSE. NAME1=NAME NS=DIMEN(1) CALL UHDGST(IDIN,'GRATING',GRAT,ISTATS(1)) CALL UHDGSD(IDIN,'LINE',LINE,ISTATS(2)) CALL UHDGSD(IDIN,'SAMPLE',SAMPLE,ISTATS(3)) CALL UHDGSD(IDIN,'DELTAS',DELTAS,ISTATS(4)) CALL UHDGSI(IDIN,'CARPOS',CARPOS,ISTATS(5)) DO 55 I=1,5 IF(ISTATS(I).NE.0)THEN CONTXT='ERROR getting input header parameter '// * HEADER(I) GO TO 999 ENDIF 55 CONTINUE ENDIF IF(DIMEN(1).NE.NS)THEN CONTXT='ERROR: ALL OBSERVATIONS MUST BE SAME LENGTH' GO TO 999 ENDIF C C READ MASK FILE INFO C IF(MASKT.NE.' ')THEN CALL UIMGID(IDIN2,DTYPE,NAXIS,DIMEN,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error reading mask file '//MNAME GO TO 999 ENDIF C C CHECK FOR VALID DATA C IF((NAXIS.NE.1).OR.(DIMEN(1).NE.NS))THEN CONTXT='Mask vector size must match input file' GO TO 999 ENDIF ENDIF C C READ DATA C CALL UIGL1D(IDIN,DATA,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error reading input file '//NAME GO TO 999 ENDIF IF(MASKT.NE.' ')THEN CALL UIGL1D(IDIN2,MASK,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error reading input file '//MNAME GO TO 999 ENDIF ENDIF C C CLOSE IMAGES C CALL UIMCLO(IDIN,ISTAT) IF(MASKT.NE.' ')CALL UIMCLO(IDIN2,ISTAT) C C ACCUMULATE AVERAGES C DO 90 I=1,NS IF(MASK(I).NE.0)THEN NADDS(I)=NADDS(I)+1 TOTAL(I)=TOTAL(I)+DATA(I) ENDIF 90 CONTINUE C C GO GET NEXT IMAGE C GO TO 20 100 CONTINUE C C DONE READING INPUT DATA ---------------------------------------------------- C C C CONSTRUCT SAMPLE VECTOR C DO 110 I=1,NS 110 SAMP(I)=(I-1)*DELTAS+SAMPLE C C FIND VALID DATA POINTS C NOUT=0 DO 120 I=1,NS IF(NADDS(I).GT.0)THEN NOUT=NOUT+1 TOTAL(NOUT)=TOTAL(I)/NADDS(I) SAMP(NOUT)=SAMP(I) ENDIF 120 CONTINUE C C READ INPUT NODE TABLE (IF NOT SUPPLIED) USE NUMBER OF NODES TO COMPUTE---- C EQUALLY SPACED NODE POSITIONS C IF(INNODE.NE.' ')THEN CALL UTTOPN(INNODE,RDONLY,IDIN,ISTAT) C --->open table IF(ISTAT.NE.0)THEN CONTXT='Error opening input table '//INNODE GO TO 999 ENDIF CALL UTPGTI(IDIN,TBNROW,NODES,ISTAT) C --->get number of rows IF(ISTAT.NE.0)THEN CONTXT='Error reading input table '//INNODE GO TO 999 ENDIF IF((NODES.GT.30).OR.(NODES.LT.2))THEN CONTXT='Error: '//INNODE//' Must have between 2 and '// * '30 node positions' GO TO 999 ENDIF CALL UTCFND(IDIN,COLNAM,1,COLIDS,ISTAT) C --->locate column IF(ISTAT.NE.0)THEN CONTXT='Error locating sample column in input table ' * //INNODE GO TO 999 ENDIF C C READ TABLE DATA C CALL UTCGTD(IDIN,COLIDS(1),1,NODES,SNODES,NULLS,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error reading input table '//INNODE GO TO 999 ENDIF CALL UTTCLO(IDIN,ISTAT) CALL UMSPUT('Input node positions taken from table '// * INNODE,STDOUT,0,ISTAT) ELSE C C USE EQUALLY SPACED NODES C IF(NODES.GT.1)THEN DIF=SAMP(NOUT)-SAMP(1) DIF=DIF/(NODES-1) DO 220 I=1,NODES SNODES(I)=SAMP(1)+DIF*(I-1) 220 CONTINUE ENDIF ENDIF C C-------------------------------------------------------------------PERFORM FIT C CALL ZSPFIT(SAMP,TOTAL,NOUT,SNODES,VNODES,NODES,NITER, * FIT,STATUS) IF(STATUS.NE.0)THEN CONTXT='ERROR DOING LEAST SQUARES SPLINE FIT' GO TO 999 ENDIF C C COMPUTE SPLINE VALUES AT SAMPLE POSITIONS OF INPUT DATA C DO 300 I=1,NS 300 SAMP(I)=(I-1)*DELTAS+SAMPLE CALL ZSPLIN(SAMP,NS,SNODES,VNODES,NODES,FIT,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error computing spline curve for output1' GO TO 999 ENDIF C C COMPUTE SPLINE VALUES AT SAMPLE POSITIONS 0.0 TO 560.0 C DO 305 I=1,561 305 SAMP(I)=I-1.0 CALL ZSPLIN(SAMP,561,SNODES,VNODES,NODES,VIG,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error computing spline curve for output2' GO TO 999 ENDIF C C NORMALIZE TO UNITY AT SAMPLE POSITION 280.0 C VNORM=VIG(281) DO 310 I=1,561 310 VIG(I)=VIG(I)/VNORM C C WRITE OUTPUT TABLE OF NODES -------------------------------------OUTTABLE C CALL UTTINN(OUTNOD,IDOUT,ISTATS(1)) CALL UTCDEF(IDOUT,COLNAM,CUNITS,CFORM,CTYPE,2,COLIDS,ISTATS(2)) CALL UTTCRE(IDOUT,ISTATS(3)) DO 200 I=1,3 IF(ISTATS(I).NE.0)THEN CONTXT='Error creating output table '//OUTNOD GO TO 999 ENDIF 200 CONTINUE C C COPY RESULTS TO TABLE C CALL UTCPTD(IDOUT,COLIDS(1),1,NODES,SNODES,ISTATS(1)) CALL UTCPTD(IDOUT,COLIDS(2),1,NODES,VNODES,ISTATS(2)) DO 210 I=1,2 IF(ISTATS(I).NE.0)THEN CONTXT='Error writing to output node table' GO TO 999 ENDIF 210 CONTINUE CALL UTTCLO(IDOUT,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error closing output table '//OUTNOD GO TO 999 ENDIF C C WRITE FIT AT INPUT SAMPLE POSITIONS ----------------------------- OUTPUT1 C CALL ZTPLAT(NAME1,OUTPT1,FIT,NS,IDOUT,ISTAT) IF (ISTAT.NE.0)THEN CONTXT='ERROR WRITING FILE OUTPUT1 '//OUTPT1 GO TO 999 ENDIF CALL UIMCLO(IDOUT,ISTAT) IF (ISTAT.NE.0)THEN CONTXT='ERROR CLOSING FILE OUTPUT1 '//OUTPT1 GO TO 999 ENDIF C C WRITE OUTPUT REFERENCE FILES USING INPUT TEMPLATE FILES C CALL ZTPLAT(TEMPLT,OUTPT2,VIG,561,IDOUT,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Unable to write output vignetting file' GO TO 999 ENDIF CALL UHDPST(IDOUT,'GRATING',GRAT,ISTATS(1)) CALL UHDPSD(IDOUT,'LINE_POS',LINE,ISTATS(2)) CALL UHDPSI(IDOUT,'CAR_POS',CARPOS,ISTATS(3)) CALL UHDPSR(IDOUT,'SAMPBEG',0.0,ISTATS(4)) CALL UHDPSR(IDOUT,'SAMPOFF',1.0,ISTATS(5)) DO 500 I=1,5 IF(ISTATS(I).NE.0)THEN CONTXT='ERROR WRITING TO OUTPUT2 HEADER' GO TO 999 ENDIF 500 CONTINUE CALL UIMCLO(IDOUT,ISTAT) IF(ISTAT.NE.0)THEN CONTXT='Error closing file '//OUTPT2 GO TO 999 ENDIF C C DONE C GO TO 1000 999 CALL UMSPUT(CONTXT,STDOUT+STDERR,0,ISTAT) 1000 CALL UIMCTP(IDTEMP,ISTAT) IF (MASKT.NE.' ') CALL UIMCTP(IDTMP2,ISTAT) END