C+
C NAME:
C       Write3Dinfotd1n
C PURPOSE:
C       Write 3D information over the nT (or at nTF, if different) times.
C CALLING SEQUENCE:
C       call Write3DInfotd1(Mo,Mn,cStrID,RR,dRR,XCbe,XCtbeg,XCtend,
C                       nLng,nLat,nMap,nT,XCint,VF,DF,V,D,
C                       XCshift,Vratio,Dratio,NCoff,Vtmp,Dtmp,CONRV,CONVT,PWV,PWR)
C INPUTS:
C       Mo              integer         0 Write out best velocities and densities
C                                       1 Write out above plus best filled velocities
C                                       2 Write out all above plus rest of info
C	Mn		integer		0 Do not produce interpolated files
C					# produce Mn iterated files between others
C       cStrID(*)       character       Carrington Rotation file name
C       RR              real            Height of deconvolution (source surface)
C       dRR             real            distance between heights
C       XCbe(2,nT)      real            starting and ending Carrington variable
C       XCtbeg          real            beginning Carrington time of arrays
C       XCtend          real            ending Carrington time of arrays
C       nLng            integer         # longitudes
C       nLat            integer         # latitudes
C       nMap            integer         # heights
C       nT              integer         # times
C       XCint(nT)       real            beginning and ending times of intervals
C       VF(nLng,nLat,nT)real            velocity map at height RR with holes
C       DF(nLng,nLat,nT)real            density map at height RR with holes
C       V(nLng,nLat,nT) real            filled velocity map at height RR
C       D(nLng,nLat,nT) real            filled density map at height RR
C       XCshift(nLng,nLat,nMap,nT) real shift amount from array
C       Vratio(nLng,nLat,nMap,nT)  real velocity difference at height
C       Dratio(nLng,nLat,nMap,nT)  real density difference at height
C       NCoff           integer         Carrington Rotation XC offset
C       Vtmp(nLng,nLat) real            Scratch array
C       Dtmp(nLng,nLat) real            Scratch array
C	CONRV		real		Velocity spatial filter
C	CONVT		real		Velocity temporal filter
C	PWV		real		G-level power to fit velocity
C	PWR		real		Radial G-level falloff to fit density
C OUTPUTS:
C FUNCTIONS/SUBROUTINES:
C       WR2DARR, Get3dTval
C-
	subroutine Write3DInfotd1n(Mo,Mn,cStrID,RR,dRR,XCbe,XCtbeg,XCtend,
     &			nLng,nLat,nMap,nT,XCint,VF,DF,V,D,
     &			XCshift,Vratio,Dratio,NCoff,Vtmp,Dtmp,CONRV,CONVT,PWV,PWR)
	character	cStrID*(*)
	real		VF	(nLng,nLat,nT),	! Maps with holes
     &			DF	(nLng,nLat,nT),
     &			V	(nLng,nLat,nT),	! Maps without holes
     &			D	(nLng,nLat,nT),
     &			XCbe	(2,nT),
     &			XCint	(nT),
     &			XCshift	(nLng,nLat,nMap,nT),
     &			Vratio	(nLng,nLat,nMap,nT),
     &			Dratio	(nLng,nLat,nMap,nT),
     &			Vtmp	(nLng,nLat),
     &			Dtmp	(nLng,nLat)

	character	cStr(10)*80,
     &			cFile*80,
     &			cSym*17		/'LIB__SIGNAL_QUIET'/,
     &			cQuiet*5,
     &			cFmt2D*7	/'(F10.4)'/

	integer		Str2Str, Flt2Str

	include		'MAPCOORDINATES.H'

	Bad = BadR4()

	iSym = iGetSymbol(cSym,cQuiet)
	I = iSetSymbol(cSym,'I',2)

       Mnp = Mn + 1
       AMnp = Mnp
       XCstr = XCtfull(1) - XCint(1)
C	print *, 'In write3dinfotd1N - XCtfull(1), XCint(1), XCstr', XCtfull(1), XCint(1), XCstr

       do N=1,nT
	
	XCbeg = XCbe(1,N)
	XCend = XCbe(2,N)

	do L=1, Mnp

	if(N.lt.nT) AnewI = (XCint(N+1) - XCint(N))/(AMnp) 

C	XCint1 = XCint(N) + .5/SUN$SiderealP -.5 
C	XCint2 = XCint(N) + .5/SUN$SiderealP +.5 

	XCint1 = XCint(N) + AnewI*(L-1) -.5 
	XCint2 = XCint(N) + AnewI*(L-1) +.5

	if((XCint(N)+(AnewI*(L-1))).gt.(XCbeg+1).and.(XCint(N)+(AnewI*(L-1))).lt.(XCbeg+2)) then 

	J = 0

	J = J+1

	I = 0
	I = I+Str2Str('Parameters CONR?, CON?T, PW?, PWR for this run :',cStr(J)(I+1:))+1
	I = I+Flt2Str(CONRV,-1,cStr(J)(I+1:))+1
	I = I+Flt2Str(CONVT,-2,cStr(J)(I+1:))+1
	I = I+Flt2Str(PWV,-2,cStr(J)(I+1:))+1
	I = I+Flt2Str(PWR,-2,cStr(J)(I+1:))

	J = J+1

	I = 0
	I = I+Str2Str('This time in C.R.s :',cStr(J)(I+1:))+1
	I = I+Flt2Str( XCtfull(N),3,cStr(J)(I+1:))

	J = J+1

	I = 0
	I = I+Str2Str('Start and end of Carrington rotation:',cStr(J)(I+1:))+1
	I = I+Flt2Str(NCoff+XCint1,-3,cStr(J)(I+1:))+1
	I = I+Flt2Str(NCoff+XCint2,-3,cStr(J)(I+1:))+1
	I = I+Flt2Str(NCoff+XCbeg,5,cStr(J)(I+1:))+1
	I = I+Flt2Str(NCoff+XCend,5,cStr(J)(I+1:))

	J = J+1

	I = 0
	I = I+Str2Str('Radial reference and grid distances (AU) :',cStr(J)(I+1:))+1
	I = I+Flt2Str( RR,5,cStr(J)(I+1:))+1
	I = I+Flt2Str(dRR,5,cStr(J)(I+1:))

	J = J+1

	I = 0
	I = I+Str2Str('Dimensions (nLng,nLat,nMap,nT) :',cStr(J)(I+1:))+1
	I = I+Int2Str(nLng,cStr(J)(I+1:))+1
	I = I+Int2Str(nLat,cStr(J)(I+1:))+1
	I = I+Int2Str(nMap,cStr(J)(I+1:))+1
	I = I+Int2Str(nT,cStr(J)(I+1:))

	J = J+1

C	print *, 'File number, NCoff, XCINT(N), N, AnewI, L', NCoff+XCint(N)+(AnewI*(L-1)),NCoff, XCint(N), N, AnewI, L

        write (cStrID,'(A,F8.3)') '_',NCoff+XCint(N)+(AnewI*(L-1))
	cFile = 't3d'//cStrID

	cStr(J) = 'Velocity map at source surface'
	call WR2DARR(1,nLng,nLat,VF(1,1,N),cFile,cFmt2D,.FALSE.,J,cStr)

	cStr(1) = 'Normalized density map at source surface'
	call WR2DARR(0,nLng,nLat,DF(1,1,N),cFile,cFmt2D,.FALSE.,1,cStr(1))

	  do K=1,nMap
	    do J=1,nLat
	        do I=1,nLng

		 if(N.lt.NT) then
		  XCsh	= XCshift(I,J,K,N) + (L-1)*(XCshift(I,J,K,N+1)-XCshift(I,J,K,N))/AMnp
		  Drat	= Dratio(I,J,K,N) + (L-1)*(Dratio(I,J,K,N+1)-Dratio(I,J,K,N))/AMnp
		  Vrat	= Vratio(I,J,K,N) + (L-1)*(Vratio(I,J,K,N+1)-Vratio(I,J,K,N))/AMnp
		 else
		  XCsh	= XCshift(I,J,K,N)
		  Drat	= Dratio(I,J,K,N)
		  Vrat	= Vratio(I,J,K,N)
		 end if

		  XC	= XCfull(I) + XCsh
		  XC    = max(XCbeg,min(XC,XCend))
C		  XCt	= XCtfull(N) + XCsh + AnewI*(L-1)
C	if(N.eq.30.and.J.eq.9.and.K.eq.1.and.I.eq.50) print *, N,I,XCt 
C	if(N.eq.31.and.J.eq.9.and.K.eq.1.and.I.eq.50) print *, N,I,XCt 
		  XCt	= XCint(N) + XCstr + XCsh + AnewI*(L)
C	if(N.eq.30.and.J.eq.9.and.K.eq.1.and.I.eq.50) 
C     &         print *, N,I,XCt,XCtfull(N),XCsh,L,AnewI*(L) 
C	if(N.eq.31.and.J.eq.9.and.K.eq.1.and.I.eq.50) 
C     &         print *, N,I,XCt,XCtfull(N),XCsh,L,AnewI*(L) 
		  XCt   = max(XCtbeg,min(XCt,XCtend))
		  call Get3DTval(XCbe,XCtbeg,XCtend,nLng,nLat,nT,DF,1,
     &		    XC,XLdeg(J),XCt,DEN)
		  call Get3DTval(XCbe,XCtbeg,XCtend,nLng,nLat,nT,VF,1,
     &		    XC,XLdeg(J),XCt,VEL)
C	if(N.eq.30.and.J.eq.9.and.K.eq.1.and.I.eq.50) print *, N,I,L,DEN,VEL
C	if(N.eq.31.and.J.eq.9.and.K.eq.1.and.I.eq.50) print *, N,I,L,DEN,VEL
		  Dtmp(I,J) = Bad
		  Vtmp(I,J) = Bad
C	if(K.eq.1.and.J.eq.1.and.I.eq.1) 
C     &	print *, N,K,J,I,XC,XLdeg(J),XCt,DEN,VEL,Drat,Vrat
		  if(DEN.ne.Bad) Dtmp(I,J) = Drat*DEN
C
C  Remove an additional 1/R density beyond 1 AU with the next 3 lines
C
C			if(K.gt.11) then
C				if(DEN.ne.Bad) Dtmp(I,J) = Drat*DEN*(10.)/(K-1.)
C			end if

		  if(VEL.ne.Bad) Vtmp(I,J) = Vrat*VEL
		end do
	    end do

	    I = 0
	    I = I+Str2Str('Velocity at radial distance:',cStr(1)(I+1:))+1
	    I = I+Flt2Str(RRhght(K),5,cStr(1)(I+1:))
	    I = I+Str2Str(' AU',cStr(1)(I+1:))+1
	    I = I+Flt2Str(XCtfull(N),3,cStr(1)(I+1:))
	    call WR2DARR(0,nLng,nLat,Vtmp,cFile,cFmt2D,.FALSE.,1,cStr(1))

	    I = 0
	    I = I+Str2Str('Normalized density at radial distance:' ,cStr(1)(I+1:))+1
	    I = I+Flt2Str(RRhght(K),5,cStr(1)(I+1:))
	    I = I+Str2Str(' AU',cStr(1)(I+1:))+1
	    I = I+Flt2Str(XCtfull(N),3,cStr(1)(I+1:))
	    call WR2DARR(0,nLng,nLat,Dtmp,cFile,cFmt2D,.FALSE.,1,cStr(1))
	end do

	if(Mo.ge.1) then
	  do K=1,nMap
	    do J=1,nLat
		do I=1,nLng

		 if(N.lt.NT) then
		  XCsh	= XCshift(I,J,K,N) + (L-1)*(XCshift(I,J,K,N+1)-XCshift(I,J,K,N))/AMnp
		  Drat	= Dratio(I,J,K,N) + (L-1)*(Dratio(I,J,K,N+1)-Dratio(I,J,K,N))/AMnp
		  Vrat	= Vratio(I,J,K,N) + (L-1)*(Vratio(I,J,K,N+1)-Vratio(I,J,K,N))/AMnp
		 else
		  XCsh	= XCshift(I,J,K,N)
		  Drat	= Dratio(I,J,K,N)
		  Vrat	= Vratio(I,J,K,N)
		 end if

		  XC	= XCfull(I) + XCsh
		  XC    = max(XCbeg,min(XC,XCend))
C		  XCt	= XCtfull(N) + XCsh + AnewI*(L-1)
		  XCt	= XCint(N) + XCstr + XCsh + AnewI*(L)
		  XCt   = max(XCtbeg,min(XCt,XCtend))
		  call Get3DTval(XCbe,XCtbeg,XCtend,nLng,nLat,nT,D,1,
     &		    XC,XLdeg(J),XCt,DEN)
		  call Get3DTval(XCbe,XCtbeg,XCtend,nLng,nLat,nT,V,1,
     &		    XC,XLdeg(J),XCt,VEL)
		  Dtmp(I,J) = Bad
		  Vtmp(I,J) = Bad
C	if(K.eq.1.and.J.eq.1.and.I.eq.1) 
C     &	print *, N,K,J,I,XC,XLdeg(J),XCt,DEN,VEL,Drat,Vrat

		  if(DEN.ne.Bad) Dtmp(I,J) = Drat*DEN
		  if(VEL.ne.Bad) Vtmp(I,J) = Vrat*VEL
		end do
	    end do

	    I = 0
	    I = I+Str2Str('Filled velocity at radial distance:',cStr(1)(I+1:))+1
	    I = I+Flt2Str(RRhght(K),5,cStr(1)(I+1:))
	    I = I+Str2Str(' AU',cStr(1)(I+1:))+1
	    I = I+Flt2Str(XCtfull(N),3,cStr(1)(I+1:))
	    call WR2DARR(0,nLng,nLat,Vtmp,cFile,cFmt2D,.FALSE.,1,cStr(1))

	    I = 0
	    I = I+Str2Str('Filled normalized density at radial distance:' ,cStr(1)(I+1:))+1
	    I = I+Flt2Str(RRhght(K),5,cStr(1)(I+1:))
	    I = I+Str2Str(' AU',cStr(1)(I+1:))+1
	    I = I+Flt2Str(XCtfull(N),3,cStr(1)(I+1:))
	    call WR2DARR(0,nLng,nLat,Dtmp,cFile,cFmt2D,.FALSE.,1,cStr(1))
	  end do
	end if

	if(Mo.ge.2) then

	 cStr(1) = 'Shift in Carrington variable (add to go back to reference surface)'
	 call WR2DARR(0,nLng,nLat*nMap,XCshift(1,1,1,N),cFile,cFmt2D,.FALSE.,1,cStr(1))

	 cStr(1) = 'Velocity ratios'
	 call WR2DARR(0,nLng,nLat*nMap,Vratio(1,1,1,N), cFile,cFmt2D,.FALSE.,1,cStr(1))

	 cStr(1) = 'Normalized density ratios'
	 call WR2DARR(0,nLng,nLat*nMap,Dratio(1,1,1,N), cFile,cFmt2D,.FALSE.,1,cStr(1))
	end if

	end if

       end do
       end do

	if (iSym .ne. 0) I = iSetSymbol   (cSym,cQuiet,2)
	if (iSym .eq. 0) I = iDeleteSymbol(cSym,2)

      return
      end