C+
C NAME:
C	CopyDtoDVN
C PURPOSE:
C	To make density maps at the velocity map time cadence. Velocity and density maps of the same latitude and longitude 
C	grid spacing must be input. The Carrington rotation beginning and endings at each time are imput, but can change 
C	with time as long as their interval lengths remain constant.  Thus, either a corotating (tested) or an inertial
C	set of maps can be converted.
C	If Mode = 0, then only an interpolated value of the density map is transferred.
C	If Mode = 1, then an interpolated value of the density map is transferred plus an interpolated map out to many
C		sigma in time at that longitude.
C	If Mode = 2, then the interpolated value is cut off at the 1.5 times the velocity temporal interval - 
C		i.e. one time from the non-valid value in time.
C CATEGORY:
C	Data processing
C CALLING SEQUENCE:
C	call CopyDtoDVN(Mode,XCbe,XCtbegV,XCtendV,XCtbegG,XCtendG,nLng,nLat,nTG,nTV,ConsT,DD,DDV,DDVV,DDT,DDD)
C INPUTS:
C	Mode		   integer	0 - transfer good time values only
C					1 - transfer both good and average time values
C					2 - Same as 1, but cut off filter after 1 time away
C	XCbe		   real		Beginning and ending V rotation V values
C	XCtbegV		   real*8	Beginning V time
C	XCtendV		   real*8	Ending V time
C	XCtbegG		   real*8	Beginning G time
C	XCtendG		   real*8	Ending G time
C	nLng		   integer	# Longitude points
C	nLat		   integer	# Latitude points
C	nTG		   integer	# G-level times
C	nTV		   integer	# velocity times
C	ConsT		   real		Time filter constant (in terms of the velocity time interval)
C	DD(nLng,nLat,nTG)  real		Density map
C OUTPUTS:
C	DDV(nLng,nLat,nTV) real		Density map interpolated at the velocity time cadence
C SCRATCH ARRAYS:
C	DDVV(nLnf,nTV)	   real
C	DDT(nTG)	   real
C	DDD(nLng)	   real
C FUNCTIONS/SUBROUTINES:
C PROCEDURE:
C	Bad values (indicated by BadR4()) are not processed in Mode 0
C MODIFICATION HISTORY:
C	NOV, 1999 B. Jackson (STEL,UCSD)
C-
	subroutine CopyDtoDVN(Mode,XCbe,XCtbegV,XCtendV,XCtbegG,XCtendG,nLng,nLat,nTG,nTV,ConsT,DD,DDV,
     &			      DDVV,DDT,DDD)

	real		DD(nLng,nLat,nTG),	! Density map
     &			DDV(nLng,nLat,nTV)	! Velocity density map
	real		XCbe(2,nTV)		! Velocity longitude index

	real		DDVV(nLng,nTV)		! Scratch arrays
	real		DDT(nTG)
	real		DDD(nLng)
	
	real		dRR	/.1/
	
	real*8		XCtbeg,XCtend,XCtfull,X,XCtbegV,XCtendV,XCtbegG,XCtendG
	real*8		XCtfullV,XCnTV,XCnTG
C
	include		'MAPCOORDINATES.H'
C
	Bad = badR4()
	ExpPWC = 50.	
	do I=1,nLng					! Longitudes and 
	    do J=1,nLat					! Latitudes are at the same grid locations at the boundary
		do N=1,nTG
		    DDT(N) = DD(I,J,N)
		end do
		do N=1,nTV				! Interpolate for available density values at the velocity times
		    nT  = nTV				! Fixed 1/12/04
		    nT1 = nTV-1
		    XCtbeg = XCtbegV
		    XCtend = XCtendV
		    XCtfullV = XCtfull(N)
		    nT  = nTG				! Fixed 1/12/04
		    nT1 = nTG-1
		    XCtbeg = XCtbegG
		    XCtend = XCtendG    
		    TI  =  XCtvar(XCtfullV)
		    dv = FLINT(-1,nTG,DDT,TI,0.00002)
		    if(dv.ne.Bad) then
			DDV(I,J,N) = dv
		    else
			DDV(I,J,N) = Bad
		    end if
C	if(I.eq.52.and.J.eq.9.and.N.eq.39) print *, I,J,N,N*2-2,nT,nTG,TI,DDT(N*2-2),DDV(I,J,N)
C	if(I.eq.52.and.J.eq.9.and.N.eq.39) print *, I,J,N,N*2-1,nT,nTG,TI,DDT(N*2-1),DDV(I,J,N)
		end do
	    end do
	end do
C
C  If no density values are present at the given velocity time, then fill in these values if possible using interpolated 
C  non-Bad time values at the same longitude and latitude, but at different density times.
C
	if(Mode.ge.1) then
	    nTV1 = nTV - 1
	    do J=1,nLat
		do N=1,nTV
		    XCbeg = XCbe(1,N)
		    XCend = XCbe(2,N)			! This will allow for an inertial reference frame if it is ever used.
		    XCtbeg = XCtbegV
		    XCtend = XCtendV
		    nT  = nTV
		    nT1 = nTV-1
		    XCnTV = XCtfull(N)			! This is the time location of the velocity time
		    XCtbeg = XCtbegG
		    XCtend = XCtendG
		    nT  = nTG
		    nT1 = nTG-1
		    do I=1,nLng

C	if(I.eq.1.and.J.eq.1.and.N.eq.1) print *, I,J,N,DDV(I,J,N)
C	if(I.eq.40.and.J.eq.9.and.N.eq.1) print *, I,J,N,DDV(I,J,N)

		      if(DDV(I,J,N).eq.Bad) then
			DDNT   = 0.
			ANDDNT = 0.
			X = (XCfull(I))
			XC  = XCvar(X)		! This is the density and/or velocity longitude location at this time
			do NN=1,nTG
			    do II=1,nLng
				DDD(II) = DD(II,J,NN)
			    end do
			    DDF = FLINT(-1,nLng,DDD,XC,0.00002)
			    if(DDF.ne.Bad) then
				XCnTG = XCtfull(NN)	! This is the time location of the density for a non-bad density value
			        XCMXC = sngl(nTV1*(XCnTV-XCnTG)/(XCtendV-XCtbegV))
				ExpPW = (XCMXC*ConsT)**2	! This is in interval units
				if(Mode.eq.2.and.abs(XCMXC).lt.1.5.and.ExpPW.lt.ExpPWC) then
				    Expot  = Exp(-ExpPW)
				    DDNT   = DDNT   + DDF*Expot
			  	    ANDDNT = ANDDNT + Expot
				else if(ExpPW.lt.ExpPWC) then
				    Expot  = Exp(-ExpPW)
				    DDNT   = DDNT   + DDF*Expot
			  	    ANDDNT = ANDDNT + Expot				
				end if
			    end if

C	if(I.eq.40.and.J.eq.9) print *, I,J,NN,XCnTV,XCnTG,XCMXC,ExpPW,DDF,DDNT,ANDDNT
C	if(I.eq.40.and.J.eq.9) print *, DDD

			end do
			if(ANDDNT.ne.0.) then
			    DDVV(I,N) = DDNT/ANDDNT
			else
			    DDVV(I,N) = Bad
			end if
		      end if

C	if(I.eq.1.and.J.eq.1) print *, I,J,N,DDV(I,J,N)
C	if(I.eq.40.and.J.eq.9) print *, I,J,N,DDV(I,J,N)

		    end do
		end do
	        do N=1,nTV
		    do I=1,nLng
		        if(DDV(I,J,N).eq.Bad) DDV(I,J,N) = DDVV(I,N)	! Replace the bad value with the interpolated value
		    end do
		end do	
	    end do
	end if

	return
	end
C+
C NAME:
C	CopyVtoVDN
C PURPOSE:
C	To make velocity maps at the density map time cadence. Velocity and density maps of the same latitude and longitude 
C	grid spacing must be input. The Carrington rotation beginning and endings at each time are imput, but can change 
C	with time as long as their interval lengths remain constant.  Thus, either a corotating (tested) or an inertial
C	set of maps can be converted.
C	If Mode = 0, then only an interpolated value of the velocity map is transferred.
C	If Mode = 1, then an interpolated value of the velocity map is transferred plus an interpolated map out to many
C		sigma in time at that longitude.
C	If Mode = 2, then the interpolated value is cut off at the 1.5 times the density temporal interval - 
C		i.e. one time from the non-valid value in time.
C CATEGORY:
C	Data processing
C CALLING SEQUENCE:
C	call CopyVtoVDN(Mode,XCbe,XCtbegG,XCtendG,XCtbegV,XCtendV,nLng,nLat,nTV,nTG,ConsT,VV,VVD,VVDD,VVT,VVV)
C INPUTS:
C	Mode		   integer	0 - transfer good time values only
C					1 - transfer both good and average time values
C					2 - Same as 1, but cut off filter after 1 time away
C	XCbe		   real		Beginning and ending G rotation values
C	XCtbegG		   real		Beginning G time
C	XCtendG		   real		Ending G time
C	XCtbegV		   real		Beginning V time
C	XCtendV		   real		Ending V time
C	nLng		   integer	# Longitude points
C	nLat		   integer	# Latitude points
C	nTV		   integer	# velocity times
C	nTG		   integer	# G-level times
C	ConsT		   real		Time filter constant (in terms of the density time interval)
C	VV(nLng,nLat,nTV)  real		Velocity map
C OUTPUTS:
C	VVD(nLng,nLat,nTG) real		Velocity map interpolated at the density time cadence
C SCRATCH ARRAYS:
C	VVDD(nLng,nTG)	   real
C	VVT(nTV)	   real
C	VVV(nLng)	   real
C FUNCTIONS/SUBROUTINES:
C PROCEDURE:
C	Bad values (indicated by BadR4()) are not processed in Mode 0
C MODIFICATION HISTORY:
C	NOV, 1999 B. Jackson (STEL,UCSD) Found error 11/14/00 in determining first part BVJ
C-
	subroutine CopyVtoVDN(Mode,XCbe,XCtbegG,XCtendG,XCtbegV,XCtendV,nLng,nLat,nTV,nTG,ConsT,VV,VVD,
     &			      VVDD,VVT,VVV)
	real		VV(nLng,nLat,nTV),	! Velocity map
     &			VVD(nLng,nLat,nTG)	! Density velocity map
	real		XCbe(2,nTG)		! Density longitude index

	real		VVDD(nLng,nTG)		! Scratch arrays
	real		VVT(nTV)
	real		VVV(nLng)
	real		dRR	/.1/
	
	real*8		XCtbeg,XCtend,XCtfull,X,XCtbegV,XCtendV,XCtbegG,XCtendG
	real*8		XCtfullG,XCnTG,XCnTV
C
	include 	'MAPCOORDINATES.H'
C
	Bad = badR4()
C	print *, 'Inputs Inside CopyVtoVDN', Mode,XCtbegG,XCtendG,XCtbegV,XCtendV,nLng,nLat,nTV,nTG,ConsT
	ExpPWC = 50.
	do I=1,nLng					! Longitudes and 
	    do J=1,nLat					! Latitudes are at the same grid locations at the boundary
		do N=1,nTV
		    VVT(N) = VV(I,J,N)
		end do
		do N=1,nTG				! Interpolate for available velocity values at the density times
		    nT  = nTG				! Fixed 1/12/04
		    nT1 = nTG-1
		    XCtbeg = XCtbegG
		    XCtend = XCtendG
		    XCtfullG = XCtfull(N)	
		    nT  = nTV				! Fixed 1/12/04
		    nT1 = nTV-1
		    XCtbeg = XCtbegV
		    XCtend = XCtendV    
		    TI  =  XCtvar(XCtfullG)
		    vd = FLINT(-1,nTV,VVT,TI,0.00002)
		    if(vd.ne.Bad) then
			VVD(I,J,N) = vd
		    else
			VVD(I,J,N) = Bad
		    end if
		    
C		    if(I.eq.50.and.J.eq.9)
C     &		    print *, I,J,N,(N+1)/2,(N+1)/2+1,nT,nTG,XCtbeg,XCtend,XCtfullG,TI,VVT((N+1)/2),VVT((N+1)/2+1),VVD(I,J,N)

		end do
	    end do
	end do
C
C  If no velocity values are present at the given density time, then fill in these values if possible using interpolated 
C  non-Bad velocity values at the same longitude and latitude, but at different times.
C
	if(Mode.ge.1) then
	    nTG1 = nTG-1
	    do J=1,nLat
	        do N=1,nTG
		    XCtbeg = XCtbegG
		    XCtend = XCtendG
		    nT  = nTG
		    nT1 = nTG - 1
		    XCnTG = XCtfull(N)			! This is the time location of the density time 
		    nT  = nTV
		    nT1 = nTV - 1
		    XCtbeg = XCtbegV
		    XCtend = XCtendV
		    XCbeg = XCbe(1,N)			! This will allow for an inertial reference frame if it is ever used. 
		    XCend = XCbe(2,N)
		    do I=1,nLng

C	if(I.eq.1.and.J.eq.1) print *, I,J,N,VVD(I,J,N)
C	if(I.eq.40.and.J.eq.4.and.N.eq.23) print *, I,J,N,VVD(I,J,N)
C	if(I.eq.40.and.J.eq.4.and.N.eq.24) print *, I,J,N,VVD(I,J,N)
C	if(I.eq.40.and.J.eq.9) print *, I,J,N,VVD(I,J,N)

		      if(VVD(I,J,N).eq.Bad) then
			VVNT   = 0.
			ANVVNT = 0.
			X = (XCfull(I))
			XC = XCvar(X)		! This is the density and/or velocity longitude location at this time
			do NN=1,nTV
			    do II=1,nLng
				VVV(II) = VV(II,J,NN)
			    end do
			    VVF = FLINT(-1,nLng,VVV,XC,0.00002)
C		if(I.eq.40.and.J.eq.4.and.N.eq.23) print *, I,J,N,NN,XC,VVF
			    if(VVF.ne.Bad) then
				XCnTV = XCtfull(NN)	! This is the time location of the velocity for a non-bad velocity value
				XCMXC = sngl(nTG1*(XCnTV-XCnTG)/(XCtendG-XCtbegG))
				ExpPW = (XCMXC*ConsT)**2	! This is in interval units
				if(Mode.eq.2.and.abs(XCMXC).lt.1.5.and.ExpPW.lt.ExpPWC) then
				    Expot  = Exp(-ExpPW)
C		if(I.eq.30.and.J.eq.4.and.N.eq.40) print *, I,J,N,NN,XCMXC,CONST,ExpPW,VVF,VVNT,Expot,ANVVNT
C		if(I.eq.40.and.J.eq.4.and.N.eq.23) print *, I,J,N,NN,XCMXC,CONST,ExpPW,VVF,VVNT,Expot,ANVVNT
C		if(I.eq.40.and.J.eq.4.and.N.eq.24) print *, I,J,N,NN,XCMXC,CONST,ExpPW,VVF,VVNT,Expot,ANVVNT
				    VVNT   = VVNT   + VVF*Expot
				    ANVVNT = ANVVNT + Expot
				else if(ExpPW.lt.ExpPWC) then
				    Expot  = Exp(-ExpPW)
				    VVNT   = VVNT   + VVF*Expot
				    ANVVNT = ANVVNT + Expot
				end if
			    end if

C	if(I.eq.40.and.J.eq.9) print *, I,J,NN,XCnTV,XCnTG,XCMXC,EXpPW,VVF,VVNT,ANVVNT
C	if(I.eq.40.and.J.eq.9) print *, VVV

			end do
			if(ANVVNT.ne.0.) then
			    VVDD(I,N) = VVNT/ANVVNT
			else
			    VVDD(I,N) = BAD
			end if
		      end if

C	if(I.eq.1.and.J.eq.1) print *, I,J,N,VVD(I,J,N)
C	if(I.eq.40.and.J.eq.4.and.N.eq.23) print *, 'FINAL VALUE',I,J,N,VVD(I,J,N)
C	if(I.eq.40.and.J.eq.4.and.N.eq.24) print *, 'FINAL VALUE',I,J,N,VVD(I,J,N)
C	if(I.eq.40.and.J.eq.9) print *, I,J,N,VVD(I,J,N)

		    end do
		end do
	        do N=1,nTG
		    do I=1,nLng
		        if(VVD(I,J,N).eq.Bad) VVD(I,J,N) = VVDD(I,N)	! Replace the bad value with the interpolated value
		    end do
		end do				
	    end do
	end if

	return
	end