SUBROUTINE VEDG3L ( * * inputs * : NPTS, V2, V3, COUNT, DCOUNT, EPOCH, LEVEL, * * output * : V2CEN, V3CEN, DV2CEN, DV3CEN, EPOMIN, EPOMAX, : STATUS) * * Module number: 15.11.2.3.2.2 * * Module name: lgaperloc * * Keyphrase: * ---------- * determine the edge and the center coordinates for large aperture of Phase III. * * Description: * ------------ * * FORTRAN name: VEDG3L.FOR * * Keywords of accessed files and tables: * -------------------------------------- * Name I/O Description / Comments * None * * Subroutines Called: * ------------------- * CDBS: * None * SDAS: * 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 SDAS * *------------------------------------------------------------------------------- * *== input: * --total number of data points of the * --input time series INTEGER NPTS * --input arrays of count rate, its standard * --deviation, V2 and V3 coordinates REAL COUNT(1), DCOUNT(1), V2(1), V3(1), * --ratio between count levels of edge * --and the maximum (center) : LEVEL * --epoch of observation DOUBLE PRECISION EPOCH(1) * *== output * --center coordinates and their standard * --errors REAL V2CEN, V3CEN, DV2CEN, DV3CEN * --epoch range DOUBLE PRECISION EPOMIN, EPOMAX * --error status INTEGER STATUS * *== local: * --error number and loop indices INTEGER STATOK, I, J, K, L, M, N, JM, JN, * --dummies : NP(4), INDX(4, 100), JPOS(4), JNEG(4) * --V2, V3 limits REAL V2MIN, V2MAX, V3MIN, V3MAX, : V2MID, V3MID, * --approximate radius : R, * --dummies : EDGECT, DEDGE, ARR(4, 100), POS, NEG, : LOC(4), DLOC(4) * --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=========================================== *------------------------------------------------------------------------------ * * look for minimum, maximum, and the middle values of V2 and V3 * V2MIN = V2(1) V2MAX = V2(1) V2MID = V2(1) V3MIN = V3(1) V3MAX = V3(1) V3MID = V3(1) * DO 10 I = 2, NPTS IF (V2(I) .LT. V2MIN) THEN V2MIN = V2(I) V3MID = V3(I) ELSE IF (V2(I) .GT. V2MAX) THEN V2MAX = V2(I) END IF * IF (V3(I) .LT. V3MIN) THEN V3MIN = V3(I) V2MID = V2(I) ELSE IF (V3(I) .GT. V3MAX) THEN V3MAX = V3(I) END IF 10 CONTINUE * * locate the center measurement by requiring its coordinates to be closest to * the middle values * R = (V2MAX - V2MIN) / 2. * DO 20 I = 1, NPTS IF (ABS(V2(I) - V2MID) .LT. 0.1 * R .AND. : ABS(V3(I) - V3MID) .LT. 0.1 * R) THEN GO TO 30 END IF 20 CONTINUE * * if no center measurement, flag error and exit * STATUS = ERRNUM(1) CONTXT = 'no center measurement' GO TO 999 * 30 CONTINUE EDGECT = LEVEL * COUNT(I) DEDGE = LEVEL * DCOUNT(I) * * sort the data points into 4 groups, corresponding to the four segments * of the cross * DO 40 L = 1, 4 NP(L) = 0 40 CONTINUE * DO 50 I = 1, NPTS IF ((V2(I) - V2MID) .GT. 0.5 * R) THEN NP(1) = NP(1) + 1 ARR(1, NP(1)) = V2(I) INDX(1, NP(1)) = I ELSE IF ((V2(I) - V2MID) .LT. -0.5 * R) THEN NP(3) = NP(3) + 1 ARR(3, NP(3)) = V2(I) INDX(3, NP(3)) = I ELSE IF ((V3(I) - V3MID) .GT. 0.5 * R) THEN NP(2) = NP(2) + 1 ARR(2, NP(2)) = V3(I) INDX(2, NP(2)) = I ELSE IF ((V3(I) - V3MID) .LT. -0.5 * R) THEN NP(4) = NP(4) + 1 ARR(4, NP(4)) = V3(I) INDX(4, NP(4)) = I END IF 50 CONTINUE * * make sure each group has more than one points * DO 60 L = 1, 4 IF (NP(L) .LT. 2) THEN STATUS = ERRNUM(2) CONTXT = 'not enough points to determine the center' GO TO 999 END IF 60 CONTINUE * * determine edge coordinates * DO 80 L = 1, 4 POS = 1.E30 NEG = -1.E30 * DO 70 J = 1, NP(L) K = INDX(L, J) IF ((COUNT(K) - EDGECT) .LT. 0. .AND. : (COUNT(K) - EDGECT) .GT. NEG) THEN NEG = COUNT(K) JNEG(L) = J ELSE IF ((COUNT(K) - EDGECT) .GT. 0. .AND. : (COUNT(K) - EDGECT) .LT. POS) THEN POS = COUNT(K) JPOS(L) = J END IF 70 CONTINUE * JM = JPOS(L) JN = JNEG(L) M = INDX(L, JM) N = INDX(L, JN) LOC(L) = ARR(L, JN) + (EDGECT - COUNT(N)) * : (ARR(L, JM) - ARR(L, JN)) / (COUNT(M) - COUNT(N)) DLOC(L) = ABS(ARR(L, JM) - ARR(L, JN)) * SQRT((DCOUNT(N) * : (EDGECT - COUNT(M))) ** 2 + (DCOUNT(M) * : (EDGECT - COUNT(N))) ** 2 + (DEDGE * : (COUNT(M) - COUNT(N))) ** 2) / : (COUNT(M) - COUNT(N)) ** 2 80 CONTINUE * * calculate the center coordinate * V2CEN = (LOC(1) + LOC(3)) / 2. V3CEN = (LOC(2) + LOC(4)) / 2. DV2CEN = SQRT(DLOC(1) ** 2 + DLOC(3) ** 2) DV3CEN = SQRT(DLOC(2) ** 2 + DLOC(4) ** 2) * * find the range of epoch * EPOMIN = EPOCH(1) EPOMAX = EPOCH(1) DO 90 I = 2, NPTS IF (EPOCH(I) .LT. EPOMIN) THEN EPOMIN = EPOCH(I) ELSE IF (EPOCH(I) .GT. EPOMAX) THEN EPOMAX = EPOCH(I) END IF 90 CONTINUE * STATUS = OK GO TO 1000 * * write error message * 999 CALL UMSPUT ('VEDG3L: ' // CONTXT, DEST, PRIO, STATOK) * 1000 RETURN END