C+ C NAME: C Translate_4r.f C PURPOSE: C Read the externally-provided volumetric elements provided by an external 3D-MHD program, and translate these into C the volumetric elements used by the tomography program C CATEGORY: C I/O C CALLING SEQUENCE: C run translate_4r C C INPUTS: C 3-D data file from 3D-MHD program in a specified coordinate system C C OUTPUTS: C 3-D data file that can be read by the UCSD tomography program C C MODIFICATION HISTORY: C November-2015, Bernard Jackson (UCSD) C- program translate_4r parameter (NLmaxG = 30000, ! Max # G data points (13000) (400,000) & NGREST = 1, ! G temporal resolution factor better than 1 day & NGRESS = 1, ! Spatial G resolution factor better than 20 degrees & NF = 2, ! LOS resolution factor (Should be 1 if Spatial resolution is >= 2) & nTmaxG = 55*NGREST) ! Max # G times (50) (100) (150) parameter (nLng = 54*NGRESS+1, ! # Longitudes, resolution is 360/(nLng-1) (55) (109) (217) (325) & iLngE = 18+1, ! # Longitudes to use in input ENLIL map (19) (37) (73) (109) & iLngM = 72+1, ! # Longitudes to use in input MS-FLUKSS map (73) (37) (73) (109) & nLat = 9*NGRESS+1, ! # Latitudes , resolution is 180/(nLat-1) (10) (19) (37) (55) & iLat = 9+1, ! # Latitudes , resolution is 180/(nLat-1) (10) (19) (37) (55) & nMap = 121, ! # Radial heights. Range covered is [0,(nMap-1)*dRR) (31) (61) (121) & nLngLat = nLng*nLat, & nLngLatnTG =nLng*nLat*nTMaxG) real XCEA(nTmaxG), & DD1E(iLngE,iLat,nMap,nTmaxG), & TT1E(iLngE,iLat,nMap,nTmaxG), & VV3E(iLngE,iLat,nMap,nTmaxG,3), & BB3E(iLngE,iLat,nMap,nTmaxG,3), & DD1M(iLngM,iLat,nMap,nTmaxG), & TT1M(iLngM,iLat,nMap,nTmaxG), & VV3M(iLngM,iLat,nMap,nTmaxG,3), & BB3M(iLngM,iLat,nMap,nTmaxG,3), & DDD1(nLng,nLat,nMap,nTmaxG), & TTT1(nLng,nLat,nMap,nTmaxG), & VVV3(nLng,nLat,nMap,nTmaxG,3), & BBB3(nLng,nLat,nMap,nTmaxG,3), & RADS(nMap) real FALLOFFD /2.0/, & FALLOFFV /0.0/, & FALLOFFBR /2.0/, & FALLOFFBT /1.0/, & FALLOFFBN /1.34/ real XCMargin /.5/ parameter (xInc = 0.5) parameter (nCar = 2000) real*8 JDCar(nCar), & JD, JEpoch, JEpoch1, JEpoch2,JDdate, & dLngSun,dLatSun,dDisSun, & DDoy,XC, & MJDref,MJDref0,MJDrefB,MJDrefE, & MJDcntr real*8 XCintDG(nTmaxG) real*8 bad8 integer MJDfrst / 46100/, & MJDlast / 60000/ character cdate*8 character cfile*80 character ccfile*80 character cWild*80 character cheader*132 character cWild3DENLIL*80 character cWild3DMSFLUKSS*80 character cWild3DMSFLUKSSO*80 character cFmtF4*9 /'(55F10.4)'/ character cMon*3 character cStr*80 logical bForecast logical bENLIL, bMSFLUKSS, bH3DMHD external EARTH8 external EARTH include 'dirspec.h' include 'sun.h' bad8 = -1.0d0 print *, ' ' C call AskWhat('Use an external input from: ENLIL, MS-FLUKSS, H3D-MHD?$1$1',I) call AskWhat('Use an external input from: ENLIL, MS-FLUKSS, H3D-MHD?$1$2',I) bENLIL = I .eq. 1 bMSFLUKSS = I .eq. 2 bH3DMHD = I .eq. 3 if(bENLIL) MHDS = 1 if(bMSFLUKSS) MHDS = 2 if(bH3DMHD) MHDS = 3 C Set up program iOffUT = 7 NmidHRV = -9 call AskI4('Number of hours of tomography reconst. from UT midnight?',NmidHRV) Doy1985 = 1. iY1985 = 1985 call MAP_TZERO(EARTH,iY1985,Doy1985,.01,nCar,JDCar) NCoff = N_CARRINGTON(iY1985,Doy1985) call Julian(11,iYr1985,Doy1985,dble(MJDfrst),JEpoch) XClo = NCoff+XMAP_SC_POS(EARTH,iYr1985,Doy1985,nCar,JDCar) XClo = XClo-.5 C call Local2UT(iOffUT,iYr1985,Doy1985) call Local2UT(NmidHRV,iYr1985,Doy1985) ! Changed 5/2/11 BVJ XChi = max(NCoff+XMAP_SC_POS(EARTH,iYr1985,Doy1985,nCar,JDCar)-1.,XClo) XChi = XChi + 1.0 ! Changed 4/23/18 BVJ C------ C iOffUT is an offset time (in hours) between the local computer time and UT. This is used only if the program C is operated in forecast mode. Note that the local time is assumed to be the same as the C computer time. The default offset is 7 hours (the difference between PDT and GST). C call AskYN('Forecast mode?$no',bForecast) ! Regular or forecast mode if (bForecast) then call AskI4('Offset between local computer time and UT (hours)?',iOffUT) NLVmin = NLVmin/2 NLGmin = NLGmin/2 end if C*************** if (bForecast) then ! Forecast mode call Local2UT(iOffUT,iYr,Doy) ! calls the local computer time clock. With the offset calculates the Doy = universal time call Julian(10,iYr,Doy,JD,JEpoch) ! Current time in MJD write (*,'(/,10X,A,I6,A,I6,A,F11.3)') 'First MJD =',MJDfrst, & ' Last MJD =',MJDlast,' Current MJD =',JD MJDref = JD ! MJD at center of plot (Earth position) write (cStr,'(A,I5,A,F7.1,A)') '$',MJDfrst,'$',JD,'$0$' print *, ' ' call AskR8('Modified Julian Day (0 = STOP)'//cStr,MJDref) if (MJDref .le. 0) call Say('IPSD2020','I','StopPlay','StopPlay') C C------- Make sure iYr,Doy match MJDref C Get heliographic coordinates for center of plot (Earth position) C MJDref = MJDref + 7 C call Julian(11,iYr,Doy,MJDref,JEpoch) Idoy = int(Doy) HR = (Doy-float(Idoy))*24.0 IH = int(HR) IMN = nint((HR-float(IH))*60.0) call DATE_DOY(1,iYr,cMon,iM,iDay,iDoy) write(*,'(A,I2.2,A,I2.2,A,I2.2,A,A,I5,A,F8.3,A)') & ' (This MJD is ',IH,':',IMN,' UT ',iDay,' ',cMon,iYr,' -- DOY',Doy,')' XCEarth = XMAP_SC_POS(EARTH,iYr,Doy,nCar,JDCar) call AdjustJDCar(XCEarth, nCar, JDCar, NCoff) XCbeg = XCEarth-.5 XCtest0 = XCbeg-1.0 ! Beginning of three rotations XCtest1 = XCbeg-XCmargin ! Beginning of needed input files XCtest2 = XCbeg+NC+1.0 ! Ending of three rotations C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ else ! Normal mode write (cStr,'(A,F7.2,A,F7.2,A)') '$',XClo,'$',XChi,'$0$' XCbeg = max(XClo,NCoff+XCbeg) if(bENLIL) XCbeg = 2114. if(bMSFLUKSS) XCbeg = 2114. call AskR4('Carrington rotation (0 = STOP)?'//cStr,XCbeg) if (XCbeg .le. 0) call Say('IPSD2020','I','StopPlay','StopPlay') NintF = 1 C call AskI4('# Intermediate rotations interpolated',NintF) NC = 1 C call AskI4('# Rotations averaged',NC) NC = min(NC,1+int(XChi-XCbeg) ) ! Stay within available data range NCC = NC XCbeg = XCbeg-NCoff call AdjustJDCar(XCbeg, nCar, JDCar, NCoff) C------- Info only: start time of Carrington rotation call Julian(1,iYr,Doy,JDCar(int(XCbeg)),JEpoch) iDoy = Doy call DATE_DOY(1,iYr,cMon,iMon,iDay,iDoy) write (*,'(/,10X,A,I4,A,I2,3A,I4,A,I3,A,/)') 'Carrington rotation ',NCoff+int(XCbeg), & ' starts on date ',iDay,'-',cMon,'-',iYr,' ( Doy ',iDoy,' )' C------- Get modified Carrington variable for center of plot C Get heliographic coordinates for center of plot (Earth position) XCEarth = XCbeg+.5 ! Center of plot I = XCEarth ! Rotation containing plot center MJDcntr = JDCar(I)+(XCEarth-I)*(JDCar(I+1)-JDCar(I)) ! Julian day for plot center call Julian(1,iYr,Doy,MJDcntr,JEpoch) ! iYr,Doy for plot center MJDcntr = MJDcntr-SUN__MJDtoJD ! Modified Julian day for plot center XCtest0 = XCbeg-1.0 ! Beginning of three rotations XCtest1 = XCbeg-XCmargin ! Beginning of needed input files XCtest2 = XCbeg+NC+1.0 ! Ending of three rotations end if I = XCtest0 MJDref0 = JDCar(I)-SUN__MJDtoJD+(XCtest0-I)*(JDCar(I+1)-JDCar(I)) ! Start of MJD for tomographic analysis I = XCtest1 MJDrefB = JDCar(I)-SUN__MJDtoJD+(XCtest1-I)*(JDCar(I+1)-JDCar(I)) ! Start MJD for data search I = XCtest2 MJDrefE = JDCar(I)-SUN__MJDtoJD+(XCtest2-I)*(JDCar(I+1)-JDCar(I)) ! End MJD for data search if(bforecast) then print *, ' ' write(*,'(A,5F10.5)') 'Forecast Carrington variables', XCbeg, XCtest0, XCtest1, XCtest2, XCEarth else write(*,'(A,4F10.5)') 'Non-forecast Carrington variables', XCbeg, XCtest0, XCtest1, XCtest2 end if if(bENLIL) then cWild3DENLIL = './ENLIL_inputs/enlil_ips_xxxxxxxxxx.dat' print *, ' ' INN = 26 C Write(*,'(A,A50)'), 'ENLIL source ', cWild3DENLIL c print *, ' ' iY = 2011 Mo = 8 iD = 13 iH = 3 dRR = 0.1 idaynum = NGREST call ASKI4('Number of ENLIL values per day?', idaynum) print *, ' ' call Julian(11,iYr,Doy,MJDrefB,JEpoch) ! Where do the ENLIL data files begin in year and day of year? idoy = int(Doy+1.0) call DATE_DOY(1,iYr,cMon,iM,iDay,iDoy) iY = iYr call AskI4('Beginning ENLIL test start year?', iY) Mo = iM call AskI4('Beginning ENLIL test start month?',Mo) iD = iDay call AskI4('Beginning ENLIL test start day?', iD) call AskI4('Beginning test start hour?', iH) ifiles = 48 call ASKI4('Number of ENLIL files input?', ifiles) print *, 'Here 0 - iY,cMon,Mo,iD',iY,cMon,Mo,iD write(cWild3DENLIL(INN:INN+9),'(I4,3I2.2)') iY,mo,iD,iH nT = ifiles LLfst = 1 LLEnd = nT C write(*,'(A,A50)') 'ENLIL begining file ',cWild3DENLIL cMon = ' ' ! Needed because cMon is set above and if not blank provides this. call DATE_DOY(0,iY,cMon,Mo,iD,iDoy) do I=1,nT DDOY = dble(iDOY) + dble(iH + ((24.0d0/idaynum)*(I-1)))/24.0d0 XCintDG(I) = DDOY C write(*,'(A,I3,F15.8)') 'I, XCintDG(I)', I, XCintDG(I) end do call enlil_read(cWild3DENLIL,INN,iLngE,iLat,nMap-1,nT,nTmaxG,LLfst,LLEnd,iY,XCintDG, ! ENLIL has 1 too few distances & RRS,dRR,RADS,DD1E,TT1E,VV3E,BB3E,iNUM) C print *, ' ' do L=1,0 C do L=1,nT do I=1,iLngE do J=1,iLat do N=1,nMap-1 if(N.eq.10.and.J.eq.5.and.I.eq.12) write(*,'(A,4I3,F6.2,F7.1,5F7.3,F7.2)') & 'After ERead',I,J,N,L, DD1E(I,J,N,L),VV3E(I,J,N,L,1),VV3E(I,J,N,L,2),VV3E(I,J,N,L,3), & BB3E(I,J,N,L,1),BB3E(I,J,N,L,2),BB3E(I,J,N,L,3),TT1E(I,J,N,L) end do end do end do end do C stop C call interpolate_enlil_m(iPr,iLngE,iLat,nMap-1,NNT,nTmaxG,nLng,nLat,nMap,iYrBG,XCintDG(I),nCar,JDCar,EarthLng,XCtest0,XCtest2, C & DD1E,TT1E,VV3E,BB3E,DDD1,TTT1,VVV3,BBB3) RRS = RADS(1) iYrBG = iY EarthLngL = 0 do I=1,nT NNT = I call SunNewcomb(0,iYrBG,sngl(XCintDG(I)),dLngSun,dLatSun,dDisSun) C**** C**** New attampted fix iskip = 1 if(iskip.eq.1) go to 999 call Julian(0,iYrBG,sngl(XCintDG(I)),JDdate,JEpoch1) call Julian(0,2000.0,0.0,JDdate,JEpoch2) sunlng19 = sngl(dLngSun) sunlat19 = sngl(dLatSun) call Precession_rot(JEpoch1,JEpoch2,sunlng19,sunlat19) EarthLng = mod(sngl(sunLng19)+180.,360.) EarthLat = -sunLat19 go to 998 999 continue C End New attempted fix C**** C C**** C As in ipstd 3D C call SunNewcomb(0,iYr,Doy,dLngSun,dLatSun,dDisSun) C EarthLng = mod(sngl(dLngSun)+180.,360.) C EarthLat = -dLatSun C call ECLIPTIC_HELIOGRAPHIC(0,iYr,Doy,EarthLng,EarthLat) C**** EarthLng = mod(sngl(dLngSun)+180.,360.) EarthLat = sngl(-dLatSun) 998 continue C iPr = 1 iPr = 0 XXXX = Earth(iYrBG,sngl(XCintDG(I))) if(ipr.eq.1) write(*,'(A,I5,3F12.6)')'Earth Heliographic sub-Earth point ',iYrBG,XCintDG(I),XXXX AriesLng = 0.0 AriesLat = 0.0 call ECLIPTIC_HELIOGRAPHIC(0,iYrBG,sngl(XCintDG(I)),AriesLng,AriesLat) if(iPr.eq.1) write(*,'(A,I5,3F12.6)')'iYrBG, DOY, Eclip Aries Lon and Lat',iYrBG,XCintDG(I),AriesLng,AiresLat if(iPr.eq.1) write(*,'(A,I5,3F12.6)')'iYrBG, DOY, Eclip Earth Lon and Lat',iYrBG,XCintDG(I),EarthLng,EarthLat call ECLIPTIC_HELIOGRAPHIC(0,iYrBG,sngl(XCintDG(I)),EarthLng,EarthLat) if(iPr.eq.1) write(*,'(A,I5,3F12.6)')'iYrBG, DOY, Helio Earth Lon and Lat',iYrBG,XCintDG(I),EarthLng,EarthLat HGILng = EarthLng - AriesLng - 75.7570 ! Why 75.7570 constant, I don't know if(HGILng.lt.0.0) HGILng = HGILng + 360.0 if(iPr.eq.1) write(*,'(A,I5,3F12.6)')'iYrBG, DOY, HGI Earth Longitude ',iYrBG,XCintDG(I),HGILng C**** EarthLng = HGILng C XCEA(I) = XMAP_SC_POS8(EARTH8,iYrBG,XCintDG(I),nCAR,JDCar) ! Earth location in Carrington coords. at DOY XCintDG if(iPr.eq.1) write(*,'(A,I12,F12.5,I6)') 'I, XCEA(I), NCOFF ',I, XCEA(I), NCOFF if(iPr.eq.1) write(*,'(A,3F12.5)') 'XCtest0, XCtest2, EarthLng',XCtest0, XCtest2, EarthLng C if(EarthLngL.ne.0.0) then C COvsINTI = (EarthLng-EarthLngL) ! Amount of IHG ENLIL rotation per UCSD interval in deg. C COvsINTE = (XCEA(I) - XCEAL)*1080./float(nLng-1) ! Amount of EHG UCSD rotation per UCSD interval in deg. C COvsINT = COvsINTE/COvsINTI ! Degrees amount of ENG relative to IHG in a UCSD interval at this time. C else C COvsINTI = (360.0*24.25/365.25)/float(idaynum) ! Average amount of IHG ENLIL rotation per UCSD interval in deg. C COvsINTE = (360.0*27.232/365.25)/float(idaynum) ! Average amount of EHG UCSD rotation per UCSD interval in deg. C COvsINT = COvsINTE/COvsINTI ! Degrees amount of ENG relative to IHG in a UCSD interval at this time. C end if if(iPr.eq.1) print *, ' ' if(iPr.eq.1) write(*,'(I4,F12.4,5F12.5)') I, XCintDG(I), EarthLng, EarthLngL, COvsINT, COvsINTI, COvsINTE iPr = 0 C iPr = 1 C iPr = 2 C NTVAL = 48 C C call interpolate_enlil(iLngE,iLat,nMap-1,NNT,nTmaxG,nLng,nLat,nMap,EarthLng,COvsINT,XCEA(I),XCtest0,XCtest2, C & DD1E,TT1E,VV3E,BB3E,DDD1,TTT1,VVV3,BBB3) C C call interpolate_enlil_n(iPr,iLngE,iLat,nMap-1,NNT,nTmaxG,nLng,nLat,nMap,iYrBG,XCintDG(I),nCar,JDCar,EarthLng,XCtest0,XCtest2, C & DD1E,TT1E,VV3E,BB3E,DDD1,TTT1,VVV3,BBB3) C C call interpolate_enlil_m(iPr,iLngE,iLat,nMap-1,NNT,nTmaxG,nLng,nLat,nMap,iYrBG,XCintDG(I),nCar,JDCar,EarthLng,XCtest0,XCtest2, C & DD1E,TT1E,VV3E,BB3E,DDD1,TTT1,VVV3,BBB3) C C if(I.eq.1.or.I.eq.30.or.I.eq.47) then C iPr = 1 C print *, ' ' C write(*,'(A,I3)') '***** Print ENLIL output in interpolate_enlil_nm **** at I =', I C print *, ' ' C end if C call interpolate_enlil_nm(iPr,iLngE,iLat,nMap-1,NNT,nTmaxG,nLng,nLat,nMap,iYrBG,XCintDG(I),nCar,JDCar,EarthLng, & XCtest0,XCtest2,DD1E,TT1E,VV3E,BB3E,DDD1,TTT1,VVV3,BBB3) C C print *, 'after interpolate_enlil ENLIL L, nt', L, nt C if(NNT.eq.NTVAL) stop iPr = 0 EarthLngL = EarthLng XCEAL = XCEA(I) end do C stop end if C****************************************** if(bMSFLUKSS) then cWild3DMSFLUKSS = './MS-FLUKSS_inputs/msflukssxxxxxxxxxx' cWild3DMSFLUKSSO = './MS-FLUKSS_nv3o_inputs/nv3oxxxxxxxxx_00001' print *, ' ' INN = 28 INNO = 29 write(*,'(A,A50)'), 'MS-FLUKSS source ', cWild3DMSFLUKSS print *, ' ' iY = 2011 Mo = 8 iD = 13 iH = 3 RRS = 0.1 ! Beginning radial height dRR = 0.025 ! Resolution (AU) in radial direction (0.1) (0.05) (0.025) idaynum = NGREST call ASKI4('Number of MS-FLUKSS values per day?', idaynum) print *, ' ' call Julian(11,iYr,Doy,MJDrefB,JEpoch) ! Where do the MS-FLUKSS data files begin in year and day of year? idoy = int(Doy) call DATE_DOY(1,iYr,cMon,iM,iDay,iDoy) iY = iYr call AskI4('Beginning MS-FLUKSS test start year?', iY) Mo = iM call AskI4('Beginning MS-FLUKSS test start month?',Mo) iD = iDay if(iY.eq.2011.and.Mo.eq.8) iD = iD + 1 ! This sets the MS-FLUKSS analysis right for CR 2114.0 call AskI4('Beginning MS-FLUKSS test start day?', iD) call AskI4('Beginning test start hour?', iH) ifiles = 48 call ASKI4('Number of MS-FLUKSS files input?', ifiles) write(cWild3DMSFLUKSS(INN:INN+9),'(I4,3I2.2)') iY,mo,iD,iH nT = ifiles LLfst = 1 LLEnd = nT C write(*,'(A,A50)') 'MS-FLUKSS begining file ',cWild3DMSFLUKSS cMon = ' ' ! Needed because cMon is set above and if not blank provides this. call DATE_DOY(0,iY,cMon,Mo,iD,iDoy) iYrBG = iY iDoyy = iDoy do I=1,nT DDOY = dble(iDoyy) + dble(iH + ((24.0d0/idaynum)*(I-1)))/24.0d0 XCintDG(I) = DDOY XCEA(I) = XMAP_SC_POS8(EARTH8,iYrBG,XCintDG(I),nCAR,JDCar) ! Earth location in Carrington coords. at DOY XCintDG write(*,'(A,I3,2F15.8)') 'I, DO XCintDG(I) XCEA(I)', I, XCintDG(I), XCEA(I) if(I.eq.2) then print *, ' ' write(cwild3DMSFLUKSS(INN:INN+8),'(F9.4)') XCEA(I) + NCOFF if(I.eq.2) write(*,'(I4,A,A50)') I, ' First MS-FLUKSS file ',cWild3DMSFLUKSS iDoy = int(XCintDG(I)) call DATE_DOY(1,iY,cMon,Mo,iD,iDoy) write(*,'(A,5I6)') ' iYrBG, mo, iD, iH, iDoy:',iYrBG, mo, iD, iH, iDoy end if end do C !MS-FLUXSS inputs at UCSD Time locations call msflukss_readt(cWild3DMSFLUKSS,INN,iLngM,iLat,nMap-1,nT,nTmaxG,LLfst,LLEnd,iY,XCintDG, ! MS-FLUKSS has 1 too few distances & RRS,dRR,RADS,DD1M,TT1M,VV3M,BB3M,iNUM) iYrBG = iY iMap = nMap-1 nMapm = nMap-iMap C print *, iYrBG, nMapM, nMap, iMap C******************* print *, 'before interpolate nt', nt do L=1,nT do I=1,0 C do I=1,iLngM do J=1,iLat do N=1,nMap-1 if(N.eq.10.and.J.eq.5) & write(*,'(4I4,4F10.4)') L,I,J,N, DD1M(I,J,N,L),VV3M(I,J,N,L,1),BB3M(I,J,N,L,1),TT1M(I,J,N,L) end do end do end do NNT = L call SunNewcomb(0,iYrBG,sngl(XCintDG(L)),dLngSun,dLatSun,dDisSun) EarthLng = mod(sngl(dLngSun)+180.,360.) EarthLat = sngl(-dLatSun) AriesLng = 0.0 AriesLat = 0.0 call ECLIPTIC_HELIOGRAPHIC(0,iYrBG,sngl(XCintDG(L)),AriesLng,AriesLat) C write(*,'(A,I5,3F12.6)')'iYrBG, DOY, Eclip Aries Lon and Lat',iYrBG,XCintDG(L),AriesLng,AiresLat C write(*,'(A,I5,3F12.6)')'iYrBG, DOY, Eclip Earth Lon and Lat',iYrBG,XCintDG(L),EarthLng,EarthLat call ECLIPTIC_HELIOGRAPHIC(0,iYrBG,sngl(XCintDG(L)),EarthLng,EarthLat) C write(*,'(A,I5,3F12.6)')'iYrBG, DOY, Helio Earth Lon and Lat',iYrBG,XCintDG(L),EarthLng,EarthLat HGILng = EarthLng - AriesLng - 75.7570 ! Why 75.7570 constant, I don't know if(HGILng.lt.0.0) HGILng = HGILng + 360.0 C write(*,'(A,I5,3F12.6)')'iYrBG, DOY, HGI Earth Longitude ',iYrBG,XCintDG(L),HGILng EarthLng = HGILng C XCEA(I) = XMAP_SC_POS8(EARTH8,iYrBG,XCintDG(L),nCAR,JDCar) ! Earth location in Carrington coords. at DOY XCintDG C write(*,'(A,I12,F12.5,I6)') 'L, XCEA(L), NCOFF ',L, XCEA(L), NCOFF C write(*,'(A,3F12.5)') 'XCtest0, XCtest2, EarthLng',XCtest0, XCtest2, EarthLng C iPr = 0 C iPr = 1 C C call interpolate_enlil(iLngM,iLat,nMap-1,NNT,nTmaxG,nLng,nLat,nMap,EarthLng,XCEA(L),XCtest0,XCtest2, C & DD1M,TT1M,VV3M,BB3M,DDD1,TTT1,VVV3,BBB3) C C if(L.eq.1.or.L.eq.30.or.L.eq.47) then C iPr = 1 C print *, ' ' C write(*,'(A,I3)') '***** Print MS-FLUKSS output in interpolate_enlil_nm **** at L =', L C print *, ' ' C end if call interpolate_enlil_nm(iPr,iLngM,iLat,nMap-1,NNT,nTmaxG,nLng,nLat,nMap,iYrBG,XCintDG(L),nCar,JDCar,EarthLng, & XCtest0,XCtest2,DD1M,TT1M,VV3M,BB3M,DDD1,TTT1,VVV3,BBB3) C if(iPr.eq.1) print *, 'after interpolate_enlil_nm MS-FLUKSS L, nt', L, nt end do C******************* end if C****************************************** C print *, ' ' do L=1,0 C do L=1,nT do I=1,nLng do J=1,nLat do N=1,nMap if(N.eq.10.and.J.eq.5.and.I.eq.22) write(*,'(A,4I3,F6.2,F7.1,5F7.3,F7.2)') & 'Before ExWrite',I,J,N,L, DDD1(I,J,N,L),VVV3(I,J,N,L,1),VVV3(I,J,N,L,2),VVV3(I,J,N,L,3), & BBB3(I,J,N,L,1),BBB3(I,J,N,L,2),BBB3(I,J,N,L,3),TTT1(I,J,N,L) end do end do end do end do C XC = XCEA(1) + float(NCOFF) C write(*,'(A,2I6,F12.5)') 'Before External Write', iYrBG, NCOFF, XC C stop print *, ' ' if(iNUM.ge.1) then ! At least one file was read from the MHD inputs call ExternalWrite(MHDs,nLng,nLat,nMap,nT,nTmaxG,LLfst,LLEnd,JDCar,nCar,XCEA,NCOFF,iYrBG,XCintDG, & RRS,dRR,FALLOFFD,FALLOFFT,FALLOFFBR,FALLOFFBT,FALLOFFBN,DDD1,TTT1,VVV3,BBB3) end if stop end