C+ C NAME: C interpolate_enlil.f C PURPOSE: C To interpolate from ENLIL maps in IHG coordinates to the UCSD rotating coordinating system at one given time. C CATEGORY: C Data processing C CALLING SEQUENCE: C call interpolate_enlil(iLng,iLat,iMap,NNT,nTmax,nLng,nLat,nMap,EarthLng,COvsINT,XCEarth,XCEA,XCtest0,XCtest2, C & DD1,TT1,VV3,BB3,DDD1,TTT1,VVV3,BBB3) C INPUTS: C iLng integer # ENLIL Longitude points C iLat integer # ENLIL Latitude points C iMap integer # ENLIL solar distances. For ENLIL, nMap - 1 C NNT integer specific time index of nTmax value C nTmax integer Maximum number of times allowed C nLng integer # UCSD tomography Longitude points C nLat integer # UCSD tomography Latitude points C nMap integer # UCSD tomography solar distances C EarthLng real Earth heliographic longitude at the time of the volume C COvsINT real Corotating vs Inertial C XCEA(I) real Earth Carrington variable value at the time of the volume C XCtest0 real Begining of the UCSD three Carrington rotation volume C XCtest2 real Ending of the UCSD three Carrington rotation volume C DD1(iLng,iLat,iMap+1,nTmax) real ENLIL 4-dimensional densities C TT1(iLng,iLat,iMap+1,nTmax) real ENLIL 4-dimensional temperatures C VV3(iLng,iLat,iMap+1,nTmax,3) real ENLIL 5-dimensional velocities C BB3(iLng,iLat,iMap+1,nTmax,3) real ENLIL 5-dimensional magnetic fields C OUTPUTS: C DDD1(nLng,nLat,nMap,nTmax) real Tomography 4-dimensional densities C TTT1(nLng,nLat,nMap,nTmax) real Tomography 4-dimensional temperatures C VVV3(nLng,nLat,nMap,nTmaxG,3) real Tomography 5-dimensional velocities C BBB3(nLng,nLat,nMap,nTmaxG,3) real Tomography 5-dimensional magnetic fields C FUNCTIONS/SUBROUTINES: C PROCEDURE: C Bad values (indicated by BadR4()) are processed C MODIFICATION HISTORY: C MAY, 1999 B. Jackson (STEL,UCSD), Modified 2008/03/03 BVJ C- subroutine interpolate_enlil(iLng,iLat,iMap,NNT,nTmax,nLng,nLat,nMap,EarthLng,COvsINT,XCEarth,XCtest0,XCtest2, & DD1,TT1,VV3,BB3,DDD1,TTT1,VVV3,BBB3) real DD1 (iLng,iLat,iMap+1,nTmax), & TT1 (iLng,iLat,iMap+1,nTmax), & VV3 (iLng,iLat,iMap+1,nTmax,3), & BB3 (iLng,iLat,iMap+1,nTmax,3), & DDD1(nLng,nLat,nMap, nTmax), & TTT1(nLng,nLat,nMap, nTmax), & VVV3(nLng,nLat,nMap, nTmax,3), & BBB3(nLng,nLat,nMap, nTmax,3) real XCbe(2,nTmax), ! Boundary values of time maps & XCint(nTmax), ! Start and end times of intervals & DMap(nLng,nLat,nTmax) ! Map C C Average the other non-Bad longitude values. C C print *, 'Into interpolate_enlil.f', NNT C Edeg = 360.0/float(iLng-1) ! Number of degrees per interval of ENLIL input volume C C Edeg Must be modified by the rotation difference between corotation and inertial. If you begin with inertial, then Edeg is slightly C smaller to fit per rotation by 24.25/27.232 modified by the location of the Earth in its orbit. C Tdeg = COvsINT*(1080.0/float(nLng-1)) ! Number of degrees per UCSD corotating volume ELOC = (XCTEST2-XCEarth)*360.0 ! Earth location in degrees from left end of the Tomography Carrington map C ELOCB = Edeg ! ENLIL beginning elongation from map left end if(ELOC. gt. 360.0) ELOC = ELOC-360.0 if(ELOC. lt. 360.0) ELOC = ELOC+360.0 C if(ELOC. gt. 0.0.and.ELOC.lt.360.0) ELOC=ELOC-360.0 write(*,'(A,F10.4,F8.3,F10.6,2F10.5)') 'Earth location, map, degs', ELOC-360., XCTEST2, XCEARTH, Edeg, Tdeg C return nMapm = nMap-iMap Icheck = 1 do I=1,nLng if(Icheck.ne.I) then write(*,'(A,I3,A,I3)') 'At NNT =',NNT,' in interpolate ENLIL no longitude value was found for I =',ICheck write(*,'(A,3I4,3F10.5)') 'I, II, IIm1, TLOC, ELOCB, ELOCBB', I-1, II, IIm1, TLOC, ELOCB, ELOCBB print *, ' ' end if ICheck = I TLOC = (I-1)*Tdeg - ELOC ! Goes from 0 to 1080 degrees - Earth's longitude if(TLOC .le. -720.0) TLOC = TLOC + 1080.0 if(TLOC .le. -360.0) TLOC = TLOC + 720.0 if(TLOC .le. 0.0) TLOC = TLOC + 360.0 if(TLOC .gt. 720.0) TLOC = TLOC - 720.0 if(TLOC .gt. 360.0) TLOC = TLOC - 360.0 if(TLOC.gt.360.0.or.TLOC.lt.0.0) then write(*,'(A,F12.5,A)') 'TLOC =', TLOC, ' and is out of bounds. STOP!' stop end if Istart = 0 9998 continue do II=1,iLng ELOCB = Edeg*float(II-1) - Earthlng ELOCBB = Edeg*float(II) - Earthlng if(ELOCB.lt.0.0.and.Istart.eq.0) then ! 1st time through ELOCB = ELOCB + 360.0 ELOCBB = ELOCBB + 360.0 end if C if(NNT.eq.14) write(*,'(A,I4,3F12.5)') 'I, TLOC, ELOCB, ELOCBB', I, TLOC, ELOCB, ELOCBB if(TLOC.ge.ELOCB.and.TLOC.lt.ELOCBB) then Icheck = Icheck + 1 C if(NNT.eq.14) write(*,'(A,I4,3F12.5)') 'Success: I, TLOC, ELOCB, ELOCBB', I, TLOC, ELOCB, ELOCBB C if(NNT.eq.14) print *, ' ' do J=1,nLat do N=1,iMap C C Now linear interpolation C IIm1 = II-1 if(IIm1.eq.0) then IIm1 = iLng - 1 if(Edeg.gt.Tdeg) then write(*,'(A,F12.6,A,F12.6,A)') 'Edeg =', Edeg, ' is greater than Tdeg =',Tdeg,' stop!' stop end if end if if(ELOCB.lt.0.0.and.ELOCBB.gt.0.0.and.Istart.eq.1) then ! Not successful 1st time through IIm1 = iLng - 1 end if DDD1(I,J,N,NNT) = (TLOC - ELOCB)*(DD1(II,J,N,NNT) - DD1(IIm1,J,N,NNT))/Edeg + DD1(IIm1,J,N,NNT) TTT1(I,J,N,NNT) = (TLOC - ELOCB)*(TT1(II,J,N,NNT) - TT1(IIm1,J,N,NNT))/Edeg + TT1(IIm1,J,N,NNT) VVV3(I,J,N,NNT,1) = (TLOC - ELOCB)*(VV3(II,J,N,NNT,1) - VV3(IIm1,J,N,NNT,1))/Edeg + VV3(IIm1,J,N,NNT,1) VVV3(I,J,N,NNT,2) = (TLOC - ELOCB)*(VV3(II,J,N,NNT,2) - VV3(IIm1,J,N,NNT,2))/Edeg + VV3(IIm1,J,N,NNT,2) VVV3(I,J,N,NNT,3) = (TLOC - ELOCB)*(VV3(II,J,N,NNT,3) - VV3(IIm1,J,N,NNT,3))/Edeg + VV3(IIm1,J,N,NNT,3) BBB3(I,J,N,NNT,1) = (TLOC - ELOCB)*(BB3(II,J,N,NNT,1) - BB3(IIm1,J,N,NNT,1))/Edeg + BB3(IIm1,J,N,NNT,1) BBB3(I,J,N,NNT,2) = (TLOC - ELOCB)*(BB3(II,J,N,NNT,2) - BB3(IIm1,J,N,NNT,2))/Edeg + BB3(IIm1,J,N,NNT,2) BBB3(I,J,N,NNT,3) = (TLOC - ELOCB)*(BB3(II,J,N,NNT,3) - BB3(IIm1,J,N,NNT,3))/Edeg + BB3(IIm1,J,N,NNT,3) C if(NNT.eq.21) write(*,'(3I4,6F10.4)'), I,J,N,DDD1(I,J,N,NNT), DD1(II,J,N,NNT), DD1(II,J,N,NNT), VVV3(I,J,N,NNT,1), C & VV3(II,J,N,NNT,1), VV3(II,J,N,NNT,1) end do end do go to 9999 end if end do if(Icheck.ne.(I+1).and.Istart.eq.0) then ! There was no success. A negative value of ELOCB was skipped C write(*,'(A,I3,A,I3)') 'At NNT =',NNT,' in interpolate ENLIL no longitude value was found for I =',I C write(*,'(A)') 'Set Istart to 1 and try again' Istart = 1 go to 9998 else write(*,'(A)') 'No solution at NNT =',NNT,' in interpolate ENLIL for longitude value at I =',I go to 9999 end if 9999 continue do J=1,nLat do N=1,nMapm DDD1(I,J,N+iMap,NNT) = DDD1(I,J,iMap,NNT) TTT1(I,J,N+iMap,NNT) = TTT1(I,J,iMap,NNT) VVV3(I,J,N+iMap,NNT,1) = VVV3(I,J,iMap,NNT,1) VVV3(I,J,N+iMap,NNT,2) = VVV3(I,J,iMap,NNT,2) VVV3(I,J,N+iMap,NNT,3) = VVV3(I,J,iMap,NNT,3) BBB3(I,J,N+iMap,NNT,1) = BBB3(I,J,iMap,NNT,1) BBB3(I,J,N+iMap,NNT,2) = BBB3(I,J,iMap,NNT,2) BBB3(I,J,N+iMap,NNT,3) = BBB3(I,J,iMap,NNT,3) end do end do end do C print *, 'End of interpolate_enlil' return end