SUBROUTINE VAPER2 * * Module number: 15.11.2.2 * * Module name: imgscale * * Keyphrase: * ---------- * Aperture location calibration Phase II. * * Description: * ------------ * Maximum number of points on each side of the area scan is 50. * * FORTRAN name: VAPER2.FOR * * Keywords of accessed files and tables: * -------------------------------------- * Name I/O Description / Comments * 'infiles' I FITS files containing input data * 'inmasks' I mask files associated with files in INFILES * 'nfiles' I number of input files * 'npeak' I number of highest points to be selected * in each input file * 'outtable' I output table name * 'low_limit' I fraction of the highest count rate, * points under such level will be rejected * 'okval' I mask value of gook points * 'hitval' I mask value of bad or particle event points * * Subroutines Called: * ------------------- * CDBS: * VAP2GT, CDMAP, CDCTRD, CDPLLS, VAP2PT * SDAS: * UCLGSI, UCLGST, UCLGSR, UMSPUT * Others: * None * * History: * -------- * Version Date Author Description * 1 01-20-87 J.-C. HSU design and coding * 2 11-10-87 J.-C. HSU F77 standard *------------------------------------------------------------------------------- * *== local: * --x and y deflection of the starting point INTEGER H0, V0, * --number of columns and rows of the * --area scan : HPT, VPT, * --horizontal and vertical steps per * --point of the are scan : HSTPT, VSTPT, * --number of integration per point : INTPT, * --total number of data points in the * --input file : NPTS, * --image/mask descripters : IMGDSC, MSKDSC, * --number of input files : NFILES, * --number of highest count points to be * --selected : NPEAK, * --number of usable files : NGOOD, * --dimension of COUNT, FLAG, X, AND Y : DIM, * --detector ID of object : DETOBJ, * --loop index : I * --instrument mode (SCP or SSP) CHARACTER*3 MODE * --HSP data format (byte, ..., all) CHARACTER*4 FORMAT * --observation ID CHARACTER*10 ROOT * --names of input file and mask CHARACTER*128 INFILE, INMASK * --observed (average) number of counts, * --its standard error, its flag and its * --coordinates of each point of the scan REAL COUNT(50, 50), DCOUNT(50, 50), : FLAG(50, 50), X(50, 50), Y(50, 50), * --input data and associated mask : DATA(2500), MASK(2500), * --V2 and V3 coordinates : V2(40), V3(40), * --peak count level : PEAKCT(40), * --weighting schemes : SCHEME(3), * --centroid coordinates and their error : XCENT(40), YCENT(40), : DXCENT(40), DYCENT(40), * --fitting coefficients : A, B, THETA, V2TRGT, V3TRGT, * --mask value of good data points : OKVAL, * --mask value of particle events : HITVAL, * --fraction of the highest count rate, * --points under such level will be rejected : LEVEL, * --maximum value : MAX * --epoch of the observation DOUBLE PRECISION EPOCH, * --V2 and V3 coordinates and their errors : VV2(40), VV3(40), DV2(40), DV3(40), * --centroid coordinates and their error : XXCENT(40), YYCENT(40), : DXX(40), DYY(40), * --chi-squared : CHISQ, * --covariance matrix : MATRIX(3, 3), * --fitting coefficients : COEFF1(3), COEFF2(3) * --output file name CHARACTER*128 OFILE * --error status INTEGER STATUS, STATOK, STAT(20) * --error message context CHARACTER*130 CONTXT * *=========================begin hsp.inc========================================= * --status return code INTEGER OK, ERRNUM(20) * --message destination and priority INTEGER DEST, PRIO DATA OK /0/ DATA ERRNUM /701, 702, 703, 704, 705, 706, 707, 708, 709, 710, : 711, 712, 713, 714, 715, 716, 717, 718, 719, 720/ DATA DEST, PRIO /1, 0/ *=========================end hsp.inc=========================================== *------------------------------------------------------------------------------ * DIM = 50 * * input parameters from PAR file * CALL UCLGST ('infiles', INFILE, STAT(1)) CALL UCLGST ('inmasks', INMASK, STAT(2)) CALL UCLGSI ('nfiles', NFILES, STAT(3)) CALL UCLGST ('outtable', OFILE, STAT(4)) CALL UCLGSR ('low_limit', LEVEL, STAT(5)) CALL UCLGSI ('npeak', NPEAK, STAT(6)) CALL UCLGSR ('okval', OKVAL, STAT(7)) CALL UCLGSR ('hitval', HITVAL, STAT(8)) * DO 10 I = 1, 8 IF (STAT(I) .NE. OK) THEN STATUS = ERRNUM(1) CONTXT = 'cannot read parameter(s) from PAR file' GO TO 999 END IF 10 CONTINUE * * maximum number of input files is 40 * IF (NFILES .GT. 40) THEN STATUS = ERRNUM(2) CONTXT = 'too many input files' GO TO 999 END IF * * initialize template processing * CALL UIMOTP (INFILE, IMGDSC, STAT(1)) CALL UIMOTP (INMASK, MSKDSC, STAT(2)) IF (STAT(1) .NE. OK .OR. STAT(2) .NE. OK) THEN STATUS = ERRNUM(3) CONTXT = 'cannot initialize template processing' GO TO 999 END IF * * repeat for each file * DO 20 I = 1, NFILES * * read parameters and data from the input file * CALL VAP2GT (IMGDSC, MSKDSC, : H0, V0, HPT, VPT, HSTPT, VSTPT, INTPT, : V2(I), V3(I), : DATA, MASK, EPOCH, : DETOBJ, MODE, FORMAT, ROOT, NPTS, STATUS) * * check array size * IF (HPT .GT. DIM .OR. VPT .GT. DIM) THEN STATUS = ERRNUM(4) CONTXT = 'oversized 2-D array' GO TO 999 END IF * IF (STATUS .NE. OK) THEN CONTXT = 'cannot read data of input file/mask' GO TO 999 END IF * * find out the coordinate of each point of the area scan * CALL CDMAP (NPTS, DIM, DATA, MASK, : OKVAL, HITVAL, H0, V0, HPT, VPT, : HSTPT, VSTPT, INTPT, : COUNT, DCOUNT, FLAG, X, Y, STATUS) * IF (STATUS .NE. OK) THEN GO TO 999 END IF * * find highest count levels of each file and calculate the coordinate of the * centroid * CALL CDCTRD (COUNT, FLAG, X, Y, OKVAL, HITVAL, DIM, NPEAK, : HPT, VPT, HSTPT, VSTPT, : PEAKCT(I), XCENT(I), YCENT(I), : DXCENT(I), DYCENT(I), STATUS) * IF (STATUS .NE. OK) THEN GO TO 999 END IF 20 CONTINUE * * look for highest count of all files * MAX = PEAKCT(1) DO 30 I = 2, NFILES IF (PEAKCT(I) .GT. MAX) THEN MAX = PEAKCT(I) END IF 30 CONTINUE * * reject low count points * NGOOD = 0 MAX = LEVEL * MAX DO 40 I = 1, NFILES IF (PEAKCT(I) .GT. MAX) THEN NGOOD = NGOOD + 1 XXCENT(NGOOD) = XCENT(I) YYCENT(NGOOD) = YCENT(I) DXX(NGOOD) = DXCENT(I) DYY(NGOOD) = DYCENT(I) VV2(NGOOD) = V2(I) VV3(NGOOD) = V3(I) END IF 40 CONTINUE * DO 50 I = 1, NGOOD DV2(I) = .01D0 DV3(I) = .01D0 50 CONTINUE * * weighting scheme * SCHEME(1) = -2. SCHEME(2) = -2. SCHEME(3) = -2. * * initial values of coefficients * DO 60 I = 1, 3 COEFF1(I) = 0.D0 COEFF2(I) = 0.D0 60 CONTINUE * IF (NGOOD. LE. 3) THEN STATUS = ERRNUM(5) CONTXT = 'too few points to fit transformation function' GO TO 999 END IF * * fit the transformation function * CALL CDPLLS (VV2, VV3, XXCENT, DV2, DV3, DXX, NGOOD, : 10, SCHEME, 1.E-12, 1., 3, 3, 3, : COEFF1, : CHISQ, MATRIX, STATUS) CALL CDPLLS (VV2, VV3, YYCENT, DV2, DV3, DYY, NGOOD, : 10, SCHEME, 1.E-12, 1., 3, 3, 3, : COEFF2, : CHISQ, MATRIX, STATUS) * * calculate the transformation coeficients * THETA = DATAN2(COEFF1(2), -COEFF1(1)) IF (COEFF1(1) .EQ. 0) THEN A = SIN(THETA) / COEFF1(2) ELSE A = -COS(THETA) / COEFF1(1) END IF * IF (COEFF2(2) .EQ. 0) THEN B = -SIN(THETA) / COEFF2(1) ELSE B = -COS(THETA) / COEFF2(2) END IF * V2TRGT = A * COEFF1(3) * COS(THETA) + B * COEFF2(3) : * SIN(THETA) V3TRGT = -A * COEFF1(3) * SIN(THETA) + B * COEFF2(3) : * COS(THETA) THETA = THETA * 57.2957795 * * write the result to output table * CALL VAP2PT (OFILE, A, B, THETA, V2TRGT, V3TRGT, NGOOD, NPEAK, : LEVEL, STATUS) * IF (STATUS .NE. OK) THEN GO TO 999 END IF * STATUS = OK GO TO 1000 * * write error message * 999 CALL UMSPUT ('VAPER2: ' // CONTXT, DEST, PRIO, STATOK) * 1000 RETURN END