subroutine readASCII(input,if_parseobs,num_stations,lat_min,lat_max,& lon_min,lon_max,vertical_min,vertical_max,time_min,time_max) ! ! Purpose: To read a GTSPP ASCII file and return the number of stations, ! the minimum and maximum values of latitude, longitude, vertical depth, ! observation time (in Julian days). ! ! Warning: This version (1.0) will not return any observational variables ! even the variable, 'if_parseobs', is set to '.TRUE.' ! implicit none character*(*), parameter :: program_name='readASCII.f90', & VERSION='1.0', VERSION_DATE='2012-05-21' character*(*), intent ( in) :: input logical, intent(in) :: if_parseobs integer, intent (out) :: num_stations real, intent(out) :: lat_min,lat_max,lon_min,lon_max real, intent(out) :: vertical_min,vertical_max real*8, intent(out) :: time_min,time_max character :: Cruise_ID*10,Data_Type*2,Iumsgno*12,Stream_Source,& Uflag,MEDS_Sta*8,Up_Date*8,Bul_Time*12,Bul_Header*6,Source_ID*4,& Stream_Ident*4,Data_Avail,QC_Version*4,Profile_Seg*2,Profile_Type*4 character :: Obs_Year*4,Obs_Month*2,Obs_Day*2,Obs_Hour*2,Obs_Minute*2 integer :: MKey,One_Deg_sq,ObsYear,ObsMonth,ObsDay,ObsHour,& ObsMinute,No_Prof,Nparms,Nsurfc,Num_Hists real :: Latitude,Longitude integer :: Q_Pos,Q_Date_Time,Q_Record,b_record character Prof_Inf*14,SPGp*15,SCGp*15,HGp*42 real, dimension(:), allocatable :: Deep_Depth integer, dimension(:), allocatable :: No_Seg character, dimension(:), allocatable :: Prof_Type*4,Dup_Flag,D_P_Code character, dimension(:), allocatable :: Pcode*4,Parm*10,Q_Parm character, dimension(:), allocatable :: SRFC_Code*4,SRFC_Parm*10,SRFC_Q_Parm character*2, dimension(:), allocatable :: Ident_Code character*4, dimension(:), allocatable :: PRC_Code character*4, dimension(:), allocatable :: Hist_Version character*8, dimension(:), allocatable :: PRC_Date character*2, dimension(:), allocatable :: Act_Code character*4, dimension(:), allocatable :: Act_Parm character*8, dimension(:), allocatable :: Aux_ID character*10,dimension(:), allocatable :: Previous_Val character, dimension(:), allocatable :: Digit_Code,Standard real, dimension(:,:,:), allocatable :: DepPres,Pval integer, dimension(:,:,:), allocatable :: Pdpq,Ppvq real, dimension(:,:), allocatable :: DepPres_total,Pval_total integer, dimension(:,:), allocatable :: Pdpq_total,Ppvq_total integer No_Depths real*8, parameter :: fildoubl = 9.9692099683868690e+36 real, parameter :: filfloat = 9.9692099683868690e+36 integer, parameter :: IN=10,max_num=1500,max_no_seg=2,maxStrLen=48000 character (len=maxStrLen) string integer :: i,j,k,n,B real*8 julday,day real obs_depth time_min = fildoubl time_max = -1.0D0 * fildoubl lat_min = filfloat lat_max = -1.0 * filfloat lon_min = filfloat lon_max = -1.0 * filfloat vertical_min = filfloat vertical_max = -1.0 * filfloat num_stations = 0 open(IN,file=input,STATUS='OLD') DO WHILE(.TRUE.) read(IN,'(a)',END=900) string call check_string_length(string,maxStrLen) num_stations = num_stations + 1 read(string(1: 8),*) MKey read(string(9:16),*) One_Deg_sq read(string(27:30),*) ObsYear read(string(31:32),*) ObsMonth read(string(33:34),*) ObsDay read(string(35:36),*) ObsHour read(string(37:38),*) ObsMinute day = real(ObsDay) call greg2jd(ObsMonth,day,ObsYear,julday) julday = julday + real(ObsHour)/24.0 + real(ObsMinute)/1440.0 time_min = min(julday,time_min) time_max = max(julday,time_max) read(string(63:70),*) Latitude read(string(71:79),*) Longitude lat_min = min(Latitude,lat_min) lat_max = max(Latitude,lat_max) lon_min = min(Longitude,lon_min) lon_max = max(Longitude,lon_max) if (if_parseobs) then Cruise_ID = string(17:26) Data_Type = string(39:40) Iumsgno = string(41:52) Stream_Source = string(53:53) Uflag = string(54:54) MEDS_Sta = string(55:52) read(string(63:70),*) Latitude read(string(71:79),*) Longitude read(string(80:80),*) Q_Pos read(string(81:81),*) Q_Date_Time read(string(82:82),*) Q_Record Up_Date = string(83:90) Bul_Time = string(91:102) Bul_Header = string(103:108) Source_ID = string(109:112) Stream_Ident = string(113:116) QC_Version = string(117:120) Data_Avail = string(121:121) read(string(122:123),*) No_Prof read(string(124:125),*) Nparms read(string(126:127),*) Nsurfc read(string(128:130),*) Num_Hists allocate(No_Seg(No_Prof)) allocate(Prof_Type(No_Prof)) allocate(D_P_Code(No_Prof)) allocate(Dup_Flag(No_Prof)) allocate(Digit_Code(No_Prof)) allocate(Standard(No_Prof)) allocate(Deep_Depth(No_Prof)) B = 131 do n = 1, No_Prof Prof_Inf = string(B:B+13); read(Prof_Inf(1:2),*) No_Seg(n) Prof_Type(n) = Prof_Inf(3:6) Dup_flag(n) = Prof_Inf(7:7) Digit_Code(n) = Prof_Inf(8:8) Standard(n) = Prof_Inf(9:9) read(Prof_Inf(10:14),*) Deep_Depth(n) B = B + 14 end do if (Nparms .ne. 0) then allocate(Pcode(Nparms)) allocate(Parm(Nparms)) allocate(Q_Parm(Nparms)) do n = 1, Nparms SPGp = string(B:B+14) Pcode(n) = SPGp(1:4) Parm(n) = SPGp(5:9) Q_Parm(n) = SPGp(15:15) B = B + 15 end do deallocate(Pcode) deallocate(Parm) deallocate(Q_Parm) endif ! if (Nsurfc .NE. 0) then allocate(SRFC_Code(Nsurfc)) allocate(SRFC_Parm(Nsurfc)) allocate(SRFC_Q_Parm(Nsurfc)) do n = 1, Nsurfc SCGp = string(B:B+14) SRFC_Code(n) = SCGp(1:4) SRFC_Parm(n) = SCGp(5:14) SRFC_Q_Parm(n) = SCGp(15:15) B = B + 15 end do deallocate(SRFC_Code) deallocate(SRFC_Parm) deallocate(SRFC_Q_Parm) endif ! if (Num_Hists .NE. 0) then allocate(Ident_Code(Num_Hists)) allocate(PRC_Code(Num_Hists)) allocate(Hist_Version(Num_Hists)) allocate(PRC_Date(Num_Hists)) allocate(Act_Code(Num_Hists)) allocate(Act_Parm(Num_Hists)) allocate(Aux_ID(Num_Hists)) allocate(Previous_Val(Num_Hists)) do n = 1, Num_Hists HGp = string(B:B+41) B = B + 42 Ident_Code(n) = HGp(1:2) PRC_Code(n) = HGp(3:6) Hist_Version(n) = HGp(7:10) PRC_Date(n) = HGp(11:18) Act_Code(n) = HGp(19:20) Act_Parm(n) = HGp(21:24) Aux_ID(n) = HGp(25:32) Previous_Val(n) = HGp(33:42) end do deallocate(Ident_Code) deallocate(PRC_Code) deallocate(Hist_Version) deallocate(PRC_Date) deallocate(Act_Code) deallocate(Act_Parm) deallocate(Aux_ID) deallocate(Previous_Val) endif allocate(DepPres(No_Prof,max_no_seg,max_num)) allocate(Pval(No_Prof,max_no_seg,max_num)) allocate(Pdpq(No_Prof,max_no_seg,max_num)) allocate(Ppvq(No_Prof,max_no_seg,max_num)) else read(string(122:123),*) No_Prof read(string(124:125),*) Nparms read(string(126:127),*) Nsurfc read(string(128:130),*) Num_Hists allocate(No_Seg(No_Prof)) allocate(DepPres(No_Prof,max_no_seg,max_num)) B = 131 do n = 1, No_Prof Prof_Inf = string(B:B+13) B = B + 14 read(Prof_Inf(1:2),*) No_Seg(n) end do end if do i = 1, No_Prof do j = 1, No_Seg(i) read(IN,'(a)',END=900) string call check_string_length(string,maxStrLen) if (if_parseobs) then read(string(1: 8),*) MKey read(string(9:16),*) One_Deg_sq Cruise_ID = string(17:26) Obs_Year = string(27:30) Obs_Month = string(31:32) Obs_Day = string(33:34) Obs_Hour = string(35:36) Obs_Minute = string(37:38) Data_Type = string(39:40) Iumsgno = string(41:52) Profile_Type = string(53:56) Profile_Seg = string(57:58) read(string(59:62),*) No_Depths D_P_Code = string(63:63) ! B = 64 do k = 1, No_Depths read(string(B:B+5),*) DepPres(i,j,k) vertical_min = min(DepPres(i,j,k),vertical_min) vertical_max = max(DepPres(i,j,k),vertical_max) B = B + 6 read(string(B:B),*) Pdpq(i,j,k) B = B + 1 read(string(B:B+8),*) Pval(i,j,k) B = B + 9 read(string(B:B),*) Ppvq(i,j,k) B = B + 1 end do else read(string(59:62),*) No_Depths B = 64 do k = 1, No_Depths read(string(B:B+5),*) DepPres(i,j,k) vertical_min = min(DepPres(i,j,k),vertical_min) vertical_max = max(DepPres(i,j,k),vertical_max) B = B + 17 end do endif end do end do deallocate(No_Seg) deallocate(DepPres) if (if_parseobs) then deallocate(Pval) deallocate(Pdpq) deallocate(Ppvq) deallocate(Prof_Type) deallocate(Deep_Depth) deallocate(D_P_Code) deallocate(Dup_Flag) deallocate(Digit_Code) deallocate(Standard) end if END DO 900 close(IN) return end subroutine readASCII subroutine check_string_length(string,maxStrLen) character*(*), intent ( in) :: string integer, intent ( in) :: maxStrLen if (len(trim(string)) .EQ. maxStrLen) then write(*,*)'WARNING: The length of input string may exceed the maximum & length allowed in the subroutine "readASCII".' write(*,*)'The length of input string=',len(trim(string)) write(*,*)'The maximum length of string allowed=',maxStrLen write(*,*)'SUGGESTION: Increase the value of the variable, maxStrLen, & recompile the program, and run it again.' end if return end subroutine check_string_length FUNCTION GREGORIAN (YEAR, MONTH, DAY, GREG_START) RESULT (GREG_FLAG) IMPLICIT NONE TYPE :: DATE_TYPE INTEGER :: YEAR_J ! year of end of Julian calendar INTEGER :: MONTH_J ! month of end of Julian calendar INTEGER :: DAY_J ! day of end of Julian calendar INTEGER :: YEAR_G ! year of start of Gregorian calendar INTEGER :: MONTH_G ! month of start of Gregorian calendar INTEGER :: DAY_G ! day of start of Gregorian calendar END TYPE DATE_TYPE LOGICAL :: GREG_FLAG ! result flag (.TRUE. for Gregorian) INTEGER, INTENT(IN) :: YEAR ! input year INTEGER, INTENT(IN) :: MONTH ! input month INTEGER, INTENT(IN) :: DAY ! input day of month TYPE (DATE_TYPE), INTENT(IN) :: GREG_START ! contains Julian stop/Gregorian start dates INTEGER :: CALTYPE = 0 ! 0=unknown, 1=Julian, 2=Gregorian IF (YEAR .LT. GREG_START%YEAR_J) THEN ! if year before end of Julian calendar.. CALTYPE = 1 ! ..then this is a Julian date ELSE IF (YEAR .EQ. GREG_START%YEAR_J) THEN ! if this is the last year of the Julian cal.. IF (MONTH .LT. GREG_START%MONTH_J) THEN ! ..then if this is before the ending month.. CALTYPE = 1 ! ..then this is a Julian date ELSE IF (MONTH .EQ. GREG_START%MONTH_J) THEN ! if this is the ending month.. IF (DAY .LE. GREG_START%DAY_J) THEN ! ..then if this is before/at the ending date.. CALTYPE = 1 ! ..then this is a Julian date END IF END IF END IF IF (YEAR .GT. GREG_START%YEAR_G) THEN ! if year after start of Gregorian calendar.. CALTYPE = 2 ! ..then this is a Gregorian date ELSE IF (YEAR .EQ. GREG_START%YEAR_G) THEN ! if this is the first year of the Greg. cal.. IF (MONTH .GT. GREG_START%MONTH_G) THEN ! ..then if this is after the starting month.. CALTYPE = 2 ! ..then this is a Gregorian date ELSE IF (MONTH .EQ. GREG_START%MONTH_G) THEN ! if this is the starting month.. IF (DAY .GE. GREG_START%DAY_G) THEN ! ..then if this is at/after the starting date.. CALTYPE = 2 ! ..then this is a Gregorian date END IF END IF END IF SELECT CASE (CALTYPE) ! check calendar type CASE (0) ! if unknown, we have an invalid date WRITE (UNIT=*, FMT='(A)') ' No such date.' ! print error message STOP ! stop program CASE (1) ! if Julian date.. GREG_FLAG = .FALSE. ! ..set return value to .false. CASE (2) ! if Gregorian date.. GREG_FLAG = .TRUE. ! ..set return value to .true. END SELECT END FUNCTION GREGORIAN SUBROUTINE GREG2JD(M,D,Y,JD) TYPE :: DATE_TYPE INTEGER :: YEAR_J ! year of end of Julian calendar INTEGER :: MONTH_J ! month of end of Julian calendar INTEGER :: DAY_J ! day of end of Julian calendar INTEGER :: YEAR_G ! year of start of Gregorian calendar INTEGER :: MONTH_G ! month of start of Gregorian calendar INTEGER :: DAY_G ! day of start of Gregorian calendar END TYPE DATE_TYPE INTEGER :: A, B ! intermediate variables DOUBLE PRECISION :: D ! day of month (+ fraction) DOUBLE PRECISION :: JD ! Julian day INTEGER :: M ! month (1-12) INTEGER :: Y ! year LOGICAL :: GREGORIAN_FLAG ! .TRUE. for Gregorian date, .FALSE. for Julian TYPE (DATE_TYPE), DIMENSION (3) :: GREGORIAN_START = & (/ DATE_TYPE (1582, 10, 4, 1582, 10, 15), & ! 1: Decree by Pope Gregory XIII DATE_TYPE (1752, 9, 2, 1752, 9, 14), & ! 2: Great Britain DATE_TYPE (1918, 1, 31, 1918, 2, 14) /) ! 3: Russia INTEGER, PARAMETER :: GREGORIAN_CHOICE = 1 ! set to 1 for 1582 date, 2 for 1752 date, etc. LOGICAL :: GREGORIAN GREGORIAN_FLAG = GREGORIAN(Y, M, INT(D), GREGORIAN_START(GREGORIAN_CHOICE)) ! test for Gregorian calendar IF (M .LE. 2) THEN Y = Y - 1 M = M + 12 END IF IF (GREGORIAN_FLAG) THEN ! Gregorian calendar A = Y/100 B = 2 - A + A/4 ELSE ! Julian calendar B = 0 END IF JD = INT(365.25D0*(Y+4716)) + INT(30.6001D0*(M+1)) + D + B - 1524.5D0 IF (.NOT. GREGORIAN_FLAG) THEN ! print msg if Julian calendar in effect WRITE (UNIT=*, FMT='(/,A)') ' Julian calendar.' END IF END SUBROUTINE GREG2JD