C+ C NAME: C MkPostd_ins C PURPOSE: C To determine line of sight positions and their projections at a given C height above the Sun. C CATEGORY: C Data processing C CALLING SEQUENCE: C call MkPostd_ins(XCbe,XCtbeg,XCtend,Idaybeg,RRS,dRR,DLOS,nLng,nLat,nMap,nT, C & XCshift,NL,NS1,NI,NIN,DIST,XLS,XCE,XP,XD,IYRS,DOYS8,NBS, C & XCproj,XClos,XLproj,XLlos,RRlos,XCtpr,XCtim,LON,LAT,ITIM, C & N_IN,XCproj_in,XClos_in,XLproj_in,XLlos_in,RRlos_in,XCtpr_in,XCtim_in,LON_in,LAT_in,ITIM_in) C INPUTS: C XCbe(2,nT) real Time interval map beginnings and endings C XCtbeg real Time interval beginning C XCtend real Time interval ending C Idaybeg integer Total integer number of days in beginning year associated with XCtbeg C RRS real Distance above Sun for deconvolution maps C dRR real Heights between individual maps (in AU) C DLOS real resolution along LOS (in units of DIST) C nLng integer # longitudes C nLat integer # latitudes C nMap integer # TS MAPS (height begins at 0 AU) C nT integer # time intervals C XCshift(nLng,nLat,nMap,nT,3) real Shifts to map to height RR (set by MkXCShifts) C NL integer # IPS data points C NS1 integer # segments of length DLOS along LOS + 1 C (the extra element is used to store Point-P info) C NI(NIN) integer Number of data values for this data set C NIN integer Total # of data sets used - 1 for IPS 1, 2 for IPS2, etc., 3 in-situ #1, 4 for in-situ #2, etc C DIST(NL) real Observer-Sun distance at source time C XLS(NL) real Ecliptic longitude of Sun C XCE(NL) real Carrington variable of Earth C XP(NL) real Longitude diff. relative to Earth-Sun C XD(NL) real Declination of Source C IYRS(NL) integer Year of source observation C DOYS8(NL) real*8 Day of year (and fraction) of source observation C NBS(NL) integer +N_IN C OUTPUTS: C XCproj(NS1,NL) real Modified Carrington variables for points on LOS, C after projection to reference surface RR C XClos(NS1,NL) real Modified Carrington variables for points on LOS C XLproj(NS1,NL) real Heliographic latitude for L.O.S. points on the reference surface C XLlos(NS1,NL) real Heliographic latitude for points on LOS C RRlos(NS1,NL) real Height above Sun for points on L.O.S. C XCtpr(NS1,NL) real*8 Time of observation on reference surface C XCtim(NS1,NL) real*8 Time of observation on L.O.S. C LON (NS1,NL) integer Int value of the heliographic longitude C LAT (NS1,NL) integer Int value of the heliographic latitude C ITIM (NS1,NL) integer Int value of time on reference surface C N_IN integer Total # of insitu data points C XCproj_in(N_IN) real Modified Carrington variables for points on LOS, C after projection to reference surface RR C XClos_in (N_IN) real Modified Carrington variables for points on LOS C XLproj_in(N_IN) real Heliographic latitude for L.O.S. points on the reference surface C XLlos_in (N_IN) real Heliographic latitude for points on LOS C RRlos_in (N_IN) real Height above Sun for points on L.O.S. C XCtpr_in (N_IN) real*8 Time of observation on reference surface C XCtim_in (N_IN) real*8 Time of observation on L.O.S. C LON_in (N_IN) integer Int value of the heliographic longitude C LAT_in (N_IN) integer Int value of the heliographic latitude C ITIM_in (N_IN) integer Int value of time on reference surface C FUNCTIONS/SUBROUTINES: C POINT_ON_LOS, ECLIPTIC_HELIOGRAPHIC, XMAP_OBS_POS C PROCEDURE: C > XClos and XCproj are on the same scale as XCE. C MODIFICATION HISTORY: C NOV, 1995 B. Jackson (STEL,UCSD) C- subroutine MkPostd_ins(XCbe,XCtbeg,XCtend,Idaybeg,RRS,dRR,DLOS,nLng,nLat,nMap,nT, & XCshift,NL,NS1,NI,NIN,DIST,XLS,XCE,XP,XD,IYRS,DOYS8,NBS, & XCproj,XClos,XLproj,XLlos,RRlos,XCtpr,XCtim,LON,LAT,ITIM, & N_IN,XCproj_in,XClos_in,XLproj_in,XLlos_in,RRlos_in,XCtpr_in,XCtim_in,LON_in,LAT_in,ITIM_in) real XCshift(nLng,nLat,nMap,nT,3),! Map of TS shifts to RR map height & XCbe(2,nt), ! Beginning and ending times & DIST(NL+N_IN), ! Distance of Sun to Earth & XLS(NL+N_IN), ! Ecliptic longitude of Sun (from Earth or Spacecraft) & XCE(NL+N_IN), ! Modified Carrington variable of Earth or Spacecraft & XP(NL+N_IN), ! Ecliptic Longitude diff. relative to Earth-Sun (Spacecraft - Sun) & XD(NL+N_IN), ! Ecliptic Latitude of Source & XCproj(NS1,NL), ! Projected Carr. rot. value of point on L.O.S. & XClos (NS1,NL), ! Point on L.O.S. Carr. rot. value & XLproj(NS1,NL), ! Latitude on the reference surface & XLlos (NS1,NL), ! Heliographic lat. (in deg.) point on L.O.S. & RRlos (NS1,NL), ! Height above Sun of point on L.O.S. C & XCproj_in(N_IN), ! Projected Carr. rot. value of point on L.O.S. & XClos_in (N_IN), ! Point on L.O.S. Carr. rot. value & XLproj_in(N_IN), ! Latitude on the reference surface & XLlos_in (N_IN), ! Heliographic lat. (in deg.) point on L.O.S. & RRlos_in (N_IN), ! Height above Sun of point on L.O.S. & XCbeg,XCend, ! Beginning and ending Carrington variable values & XC(3) ! Internal array real*8 DOYS8(NL+N_IN) ! Day of year (1.0 at beginning) (and fraction) of L.O.S. obs. (for SMEI, double precision) real*8 Doy8, ! Day of Year (1.0 at beginning) to use for ECLIPTIC_HELIOGRAPHIC8 & XCtstrt, & Adaybeg real*8 XCtbeg,XCtend,XCtfull,X,tXC,XCtinc,T1 real*8 XCtpr (NS1,NL), ! Projected time on reference surface & XCtim (NS1,NL) ! Time in 3D space along LOS real*8 XCtpr_in(N_IN), ! Projected time on reference surface & XCtim_in(N_IN) ! Time in 3D space along LOS integer IYRS(NL+N_IN), ! Year of source observation & LON(NS1,NL), ! Nearest longitude index & LAT(NS1,NL), ! Nearest latitude index & ITIM(NS1,NL), ! Nearest time index & LON_in(N_IN), ! Nearest longitude index & LAT_in(N_IN), ! Nearest latitude index & ITIM_in(N_IN) ! Nearest time index integer NI(NIN), ! Number of data values for this data set & NBS(NL+N_IN) ! Bad source flag - 1 = good, 0 = Bad parameter (DLOS25 = 0.25) include 'MAPCOORDINATES.H' XCinc = (XCbe(2,1)-XCbe(1,1))/(2.*nLng1) ! Increment beyond beginning and end longitude XCtinc = (XCtend-XCtbeg)/(2.D0*dble(nT1)) ! Increment beyond end of time XCbe1 = XCbe(1,1) ! Beginning longitude limit XCtstrt = XCtbeg ! Beginning time in DOY XCfin = XCbe(2,nT) ! Ending longitude limit C print *, 'XCbe1, XCfin, XCtstrt', XCbe1,XCfin,XCtstrt NB = 1 do II=1,NIN ! Cycle through the number of data sets used C C print *, 'in mkpostd_ins', II,NIN,NB,NI(II) C if(II .eq. 3) print *, ' Change to STEREO-A' C if(II .eq. 4) print *, ' Change to STEREO-B' C if(NI(II) .ne. 0) then ! Only go through the loop below if the value of the last element of the loop is not zero do I=NB,NB+NI(II) - 1 ! Number of data values plus number of insitu values N_IN iYr = IYRS(I) Doy8 = DOYS8(I) ! Double precision DOY value of line-of-sight observation ObsR = DIST(I) ! Heliocentric distance observer ObsLng = XLS(I)-180. ! Heliocentric ecliptic longitude Earth or Spacecraft (Change from Topocentric) ObsLat = -XD(I) ! Heliocentric ecliptic latitude Earth or Spacecraft (put in place BVJ Mar 2014) (From Topocentric) ObsXC = XCE(I) ! Modified Carrington variable Earth or Spacecraft C If(I.lt.10) print *, 'Inside mkpostd_ins',IYRS(I),DOYS8(I),DIST(I),ObsLng,OBsxc C tXC = DOYS(I)-XCtstrt ! Time of source observation (days since start) in single precision tXC = DOYS8(I) ! Time of source observation (days since year beginning) in double precision Adaybeg = dble(Idaybeg) ! This is the number of days in the preceeding year if(XCtend.gt.Adaybeg.and.DOYS8(I).lt.dble(180.)) then ! Handle case where one year to another is crossed C tXC = DOYS8(I) + Adaybeg - 1.0d0 ! Add the number of days in the previous year minus 1, since DOY begins at 1.0 C Doy8 = Doy8 + Adaybeg - 1.0d0 ! Add the number of days in the previous year minus 1, since DOY begins at 1.0 tXC = DOYS8(I) + Adaybeg ! Add the number of days in the previous year. Fixed BVJ Jan 2014 Doy8 = Doy8 + Adaybeg ! Add the number of days in the previous year. Fixed BVJ Jan 2014 iYr = iYr - 1 ! The source year time should be one greater, so subtract one year end if C If(I.lt.3) then C print *, 'Inside mkpostd_ins',IYRS(I),iYr,tXC,DOYS8(I),DIST(I),ObsLng,OBsXC,Adaybeg C print *,XCinc,XCbe(2,1),XCbe(1,1),XCtinc,XCtend,Xctbeg,XCfin,Idaybeg C end if if(I.le.NL) then ! These are the IPS data points do J=1,NS1 ! Loop over all segments along LOS ! Define point P along LOS: P1 = XP(I) ! Topocentric ecliptic longitude diff. with Spacecraft or Earth-Sun line (0.0 for Earth) P2 = XD(I) ! Topocentric ecliptic latitude of Spacecraft or Earth (0.0 for Earth) P3 = (J-.5)*DLOS ! Distance along LOS (units of ObsR) ! Point-P distance along LOS (units of ObsR) if (J .eq. NS1) P3 = max(DLOS25,cosd(P2)*cosd(P1)) ! Convert to heliocentric coordinates T1 = tXC-dble(P3*ObsR*TCD) ! Time of L.O.S. segment obs. T1 = min(XCtend+XCtinc,max(XCtbeg-XCtinc,T1)) ! No T1 outside limits XCtim(J,I) = T1 call POINT_ON_LOS(P1,P2,P3,ELO,ELOS,iEorW) P3 = P3*ObsR ! Convert to AU P1 = ObsLng+P1 ! Heliocentric longitude of P call ECLIPTIC_HELIOGRAPHIC8(0,iYr,Doy8,P1,P2) ! Convert to heliographic coordinates P1 = XMAP_OBS_POS(ObsXC,P1) ! Modified Carrington variable of P P1 = min(XCfin+XCinc,max(XCbe1-XCinc,P1)) ! No P1 outside longitude limits C if(J.eq.1.and.I.lt.150) write(*,'(3I6,11F10.4)') I, iYr, NBS(I), DOY8, ObsR, ObsXC, XP(I), ObsLng, XD(I),ObsLat,P1,P2,P3,T1 C print *, ' ' C print *, 'before Get4Dval - I J T1 P1 P2 P3', I,J, T1,P1,P2,P3 call Get4Dval(3,XCbe,XCtbeg,XCtend,RRS,dRR,nLng,nLat,nMap,nT, & XCshift,1,T1,P1,P2,P3,XC) C if(J.eq.20.and.I.eq.1000)print *, 'after Get4Dval - XC =', XC(1), XC(2), XC(3) XClos(J,I) = P1 ! Modified Carrington variable of P XLlos(J,I) = P2 ! Heliographic latitude of P RRlos(J,I) = P3 ! Heliocentric distance of P (AU) XCproj(J,I)= P1+XC(1) ! Projected modified Carrington variable XCproj(J,I)=min(XCfin+XCinc,max(XCbe1-XCinc,XCproj(J,I))) ! No XCproj outside longitude limits XLproj(J,I)= P2+XC(2) ! Projected heliographic latitude of P XLproj(J,I)=max(-90.,min(90.,XLproj(J,I))) ! No XLproj outside latitude limits XCtpr(J,I) = T1+dble(XC(3)) ! Proj. time from beginning CL value XCtpr(J,I) = min(XCtend+XCtinc,max(XCtbeg-XCtinc,XCtpr(J,I))) ! No XCtpr outside limits ITIM(J,I)= nint(XCtvar(XCtpr(J,I))) ! Nearest interval time index C if(ITIM(J,I).eq.3.and.J.eq.1) print *, ' ' C if(ITIM(J,I).eq.3.and.J.eq.1) print *, ' ITIM', ITIM(J,I), J C if(ITIM(J,I).eq.3.and.J.eq.NS1) print *, ' ITIM', ITIM(J,I), J C print *, ITIM(J,I),XCtpr(J,I),T1,dble(XC(3)),T1+dble(XC(3)) ITIM(J,I)= max(1,min(nT,ITIM(J,I))) ! No IT outside limits LAT(J,I) = nint(XLindx(XLproj (J,I))) ! Nearest latitude index LAT(J,I) = max(1,min(nLat,LAT(j,I))) ! No LAT outside limits XCbeg = XCbe(1,ITIM(J,I)) XCend = XCbe(2,ITIM(J,I)) LON(J,I) = nint(XCvar(XCproj(J,I))) ! Nearest longitude index LON(J,I) = max(1,min(nLng,LON(J,I))) ! No LON outside limits C if(ITIM(J,I).eq.3.and.J.eq.1) then C if(I.eq.10.and.J.eq.1) then C print *, 'At end of Mkpos.', J,I,XCproj(J,I),XClos(J,I),XLproj(J,I),XLlos(J,I), C & RRlos(J,I),XCtpr(J,I),XCtim(J,I),LON(J,I),LAT(J,I),ITIM(J,I) C print *, 'At end of Mkpos.', XCbeg, XCend, XCtpr(J,I), C & XCtvar(XCtpr(J,I)), ITIM(J,I),XCproj(J,I),LON(J,I), LAT(J,I), J C end if C if(ITIM(J,I).eq.3.and.J.eq.NS1) then C if(I.eq.10.and.J.eq.NS1) then C print *, 'At end of Mkpos.', XCbeg, XCend, XCtpr(J,I), C & XCtvar(XCtpr(J,I)), ITIM(J,I),XCproj(J,I),LON(J,I), LAT(J,I), J C end if end do else ! These are the in-situ data points I >= NL P1 = XP(I) ! Topocentric ecliptic longitude diff. with Earth-Sun line (0.0 for Earth) P2 = XD(I) ! Topocentric ecliptic latitude of spacecraft (0.0 for Earth) P3 = ObsR ! Distance of the object (Earth or spacecraft from the Sun in AU C P1 = ObsLng+P1 ! Heliocentric longitude of the in-situ point (mod. by P1 not needed) P1 = ObsLng ! Heliocentric ecliptic longitude of the Earth or Spacecraft P2 = ObsLat ! Heliocentric ecliptic latitude of the Earth or Spacecraft call ECLIPTIC_HELIOGRAPHIC8(0,iYr,Doy8,P1,P2) ! Convert to heliographic coordinates (needed to determine P2) if(P1.gt.360.0) P1 = P1 - 360.0 if(P1.lt.360.0) P1 = P1 + 360.0 P1 = XMAP_OBS_POS(ObsXC,P1) ! Modified Carrington variable of the ObsLng (same as ObsXC) if(abs(P1-ObsXC).gt.0.00001) then write(*,'(A,F10.5,A,F10.5)') 'The value of P1 =',P1,' is different from ObsXC =',ObsXC write(*,'(A)') 'We need to stop and reprogram the spacecraft in-situ input' stop end if P1 = min(XCfin+XCinc,max(XCbe1-XCinc,P1)) ! No P1 outside longitude limits T1 = tXC T1 = min(XCtend+XCtinc,max(XCtbeg-XCtinc,T1)) ! No T1 outside limits if(T1 .ne. tXC) NBS(I) = 0 ! set source flag to zero if outside time limit. C inext44 = 0 if(inext44.eq.0) go to 44 if(I .gt. NL .and. I .lt. NL+2) write(*,'(A)') 'inside mkpostd_ins Earth ' if(I .lt. NL+5) & write(*,'(3I6,11F10.4)') I, iYr, NBS(I), DOY8, ObsR, ObsXC, XP(I), ObsLng, XD(I),ObsLat,P1,P2,P3,T1 C if(I .gt. (NL+1220) .and. I .lt. (NL+1222)) print *, ' ' if(I. gt. (NL+1220) .and. I .lt. (NL+1222)) write(*,'(A)') 'inside mkpostd_ins Earth-STEREO-A ' if(I .gt. (NL+1238) .and. I .lt. (NL+1247)) & write(*,'(3I6,11F10.4)') I, iYr, NBS(I), DOY8, ObsR, ObsXC, XP(I), ObsLng, XD(I),ObsLat,P1,P2,P3,T1 C if(I .gt. (NL+1242+1190) .and. I .lt. (NL+1242+1192)) print *, ' ' if(I .gt. (NL+1242+1190) .and. I .lt. (NL+1242+1192)) write(*,'(A,2I5,10F11.4)') 'inside mkpostd_ins STEREO-A-B ' if(I .gt. (NL+1242+1198) .and. I .lt. (NL+1242+1206)) & write(*,'(3I6,11F10.4)') I, iYr, NBS(I), DOY8, ObsR, ObsXC, XP(I), ObsLng, XD(I),ObsLat,P1,P2,P3,T1 C if(I .gt. (NL+1242+1202+1110) .and. I .lt. (NL+1242+1202+1112)) print *, ' ' if(I .gt. (NL+1242+1202+1110) .and. I .lt. (NL+1242+1202+1112)) write(*,'(A)') 'inside mkpostd_ins at STEREO-B End' if(I .gt. (NL+1242+1202+1112)) & write(*,'(3I6,11F10.4)') I, iYr, NBS(I), DOY8, ObsR, ObsXC, XP(I), ObsLng, XD(I),ObsLat,P1,P2,P3,T1 44 continue C XCtim_in(I-NL) = T1 call Get4Dval(3,XCbe,XCtbeg,XCtend,RRS,dRR,nLng,nLat,nMap,nT, & XCshift,1,T1,P1,P2,P3,XC) XClos_in(I-NL) = P1 ! Modified Carrington variable of in-situ point XLlos_in(I-NL) = P2 ! Heliographic latitude of in-situ point RRlos_in(I-NL) = P3 ! Heliocentric distance of in-situ point (AU) XCproj_in(I-NL)= P1+XC(1) ! Projected modified Carrington variable of in-situ point XCproj_in(I-NL)= min(XCfin+XCinc,max(XCbe1-XCinc,XCproj_in(I-NL))) ! No XCproj outside longitude limits XLproj_in(I-NL)= P2+XC(2) ! Projected heliographic latitude of in-situ point XLproj_in(I-NL)= max(-90.,min(90.,XLproj_in(I-NL))) ! No XLproj outside latitude limits XCtpr_in(I-NL) = T1+dble(XC(3)) ! Proj. time from beginning CL value XCtpr_in(I-NL) = min(XCtend+XCtinc,max(XCtbeg-XCtinc,XCtpr_in(I-NL))) ! No XCtpr outside limits ITIM_in(I-NL) = nint(XCtvar(XCtpr_in(I-NL))) ! Nearest interval time index ITIM_in(I-NL) = max(1,min(nT,ITIM_in(I-NL))) ! No IT outside limits LAT_in(I-NL) = nint(XLindx(XLproj_in(I-NL))) ! Nearest latitude index LAT_in(I-NL) = max(1,min(nLat,LAT_in(I-NL))) ! No LAT outside limits XCbeg = XCbe(1,ITIM_in(I-NL)) XCend = XCbe(2,ITIM_in(I-NL)) LON_in(I-NL) = nint(XCvar(XCproj_in(I-NL))) ! Nearest longitude index LON_in(I-NL) = max(1,min(nLng,LON_in(I-NL))) ! No LON outside limits end if end do end if NB = NB + NI(II) end do C stop return end