C @(#)kolmogorov.for 17.1.1.1 (ES0-DMD) 01/25/02 17:12:45 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+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C.COPYRIGHT: Copyright (c) 1987 European Southern Observatory, C all rights reserved C C.VERSION: 1.0 ESO-FORTRAN Conversion, 29 JUN 89 C C.LANGUAGE: F77+ESOext C C.AUTHOR: M . PERON C C.IDENTIFICATION C C Program KOLMOGOROV C C.KEYWORDS C C C C.PURPOSE C C Test de Kolmogorov for one or two distributions. C C references: Numerical Recipes p472 C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ REAL FUNCTION PROBKS(ALAM) REAL EPS1,EPS2,TERM,FAC,TERMBF,ALAM REAL A2 INTEGER J PARAMETER (EPS1=0.001,EPS2=1.E-8) A2 = -2*ALAM**2 FAC = 2. PROBKS = 0. TERMBF = 0. DO 10 J=1,100 TERM = FAC*EXP(A2*J**2) PROBKS = PROBKS+TERM IF (ABS(TERM).LT.EPS1*TERMBF.OR.ABS(TERM).LE.EPS2*PROBKS) 1 RETURN FAC = -FAC TERMBF = ABS(TERM) 10 CONTINUE PROBKS = 1. RETURN END SUBROUTINE KOLM2DIS(X1,N1,X2,N2,D,PROB) INTEGER J1,J2,N1,N2 REALX1(N1),X2(N2),D,PROB,FN1,FN2,EN1,EN2 REAL D1,D2,DT,PROBKS J1 = 1 J2 = 1 FN1 = 0. FN2 = 0. EN1 = N1 EN2 = N2 D = 0. 10 IF(J1.LE.N1.AND.J2.LE.N2)THEN D1 = X1(J1) D2 = X2(J2) IF (D1.LE.D2) THEN FN1 = J1/EN1 J1 = J1+1 ENDIF IF (D2.LE.D1)THEN FN2 = J2/EN2 J2 = J2+1 ENDIF DT = ABS(FN2-FN1) IF (DT.GT.D) D = DT GOTO 10 ENDIF PROB = PROBKS(SQRT(EN1*EN2/(EN1+EN2))*D) RETURN END SUBROUTINE KOLM1D(DATA,N,DISTRI,PARA,D,PROB) INTEGER J,N REAL DATA(N),D,PROB,EN,F0,PARA(2),FN,FF,DT REAL ERFCC,GAMMAQ,PROBKS,UNIF,CONST CHARACTER*1 DISTRI EN = N D = 0. F0 = 0. DO 10 J=1,N FN = J/EN IF (DISTRI.EQ.'U') THEN FF = 1-UNIF(DATA(J),PARA) ELSE IF (DISTRI.EQ.'G') THEN CONST = 1./(SQRT(2.)*PARA(2)) FF = 1-0.5*ERFCC((DATA(J)-PARA(1))*CONST) ELSE IF (DISTRI.EQ.'P') THEN FF = GAMMAQ(DATA(J),PARA(1)) ELSE GOTO 80 ENDIF DT = MAX(ABS(F0-FF),ABS(FN-FF)) C D = MAX(ABS(FN-FF),D) IF (DT.GT.D) D = DT F0 = FN 10 CONTINUE PROB = PROBKS(SQRT(EN)*D) 80 RETURN END