C.........................................................Aug. 2018
C JPL maps 1 X 1 glat X glon 
C.........................................................Feb. 2015
C
c VERIFICATION OF MEDIAN ================
      subroutine smedstd3mm(CC,SMED,STD,STD1,STD2)
C
C Median smed, and standard deviation std, std1(neg),std2(pos) for -14+1 days at specified hour
C
C Input file: CC(ndnr,360) for ndnr days, 
C ncol=360 for longitudes -179.5:1:179.5 or ncol=360 for hours UT=0,1,...,23. MIN=00, 15, 30, 45
C Number days for median : ndnr=15	(14 preceding days+current day)
C Output arrays: SMED(360) for longitudes -179.5:1:179.5 
C                 STD1(360),STD2(360)
C
	DIMENSION CC(15,360),SMED(360),STD1(360),STD2(360)
	+,X1(15),X2(15),STD(360)
C--		COMMON /BL1/CC,SMED,STD,STD1,STD2
C "MEDIAN:"
	ncol=360
	ndnr=15
       DO 1030 M=1,ncol
      MM=0
C Add calculation of average:
	aver=0.
      DO N=1,ndnr
      X1(N)=CC(N,M)
	aver=aver+x1(n)
      ENDDO 
	aver=aver/ndnr
  880 XX=0.
CREM MAX X1(N)
  890 DO N=1,ndnr
  900 IF (XX.lt.X1(N)) THEN 
      XX=X1(N)
	ENDIF
      ENDDO      !N
  920 IF (XX.eq.0) GOTO 980
C
  930 DO 960 N=1,ndnr
      IF (X1(N).lt.XX) GOTO 960
      MM=MM+1 
      X2(MM)=XX 
	X1(N)=0.
  960 CONTINUE      ! N
  970 IF (MM.lt.NDNR) GOTO 880
  980 XMM=FLOAT(MM)
      XD1=XMM/2.
      ND1=INT(XD1)
	CD1=FLOAT(ND1) 
      ND2=ND1+1 
      IF (CD1.LT.XD1) THEN 
	XMD=X2(ND2)
	else
	XMD=(X2(ND1)+X2(ND2))/2.
	ENDIF
C
      SMED(M)=XMD
C Standard deviation:
	sumneg=0.
	sumpos=0.
	icntpos=0
	icntneg=0
	avneg=0.
	avpos=0.
	iavneg=0
	iavpos=0
	do ii=1,ndnr
	if (x2(ii).gt.xmd) then
	sumpos=sumpos+(x2(ii)-xmd)**2
	icntpos=icntpos+1
	endif
	if (x2(ii).lt.xmd) then
	sumneg=sumneg+(x2(ii)-xmd)**2
	icntneg=icntneg+1
	endif
	if (x2(ii).gt.aver) then
	avpos=avpos+(x2(ii)-aver)**2
	iavpos=iavpos+1
	endif
	if (x2(ii).le.aver) then
	avneg=avneg+(x2(ii)-aver)**2
	iavneg=iavneg+1
	endif
	enddo
	STD(M)=0.0
	STD1(M)=0.0
	STD2(M)=0.0
	if ((icntneg+icntpos).gt.0) then
	STD(M)=sqrt((sumneg+sumpos)/(icntneg+icntpos))
	else
	if((iavneg+iavpos).gt.0) STD(M)=sqrt((avneg+avpos)/
	&(iavneg+iavpos))
	endif
	if (icntneg.gt.0) then 
	if (iavneg.gt.0) STD1(M)=sqrt(sumneg/icntneg)
	else
	if (iavneg.gt.0) STD1(M)=sqrt(avneg/iavneg)
	endif
	if (icntpos.gt.0) then 
	STD2(M)=sqrt(sumpos/icntpos)
	else
	if (iavpos.gt.0) STD2(M)=sqrt(avpos/iavpos)
	endif
	if (std1(m).eq.0.0) STD1(M)=STD(M)
	if (std2(m).eq.0.0) STD2(M)=STD(M)
 1030 CONTINUE
C 
c
 1050 RETURN
      END