       SUBROUTINE PCMIN3(SST_IN,PSL_IN,P,T_IN,R_IN,IMAX,JMAX,NA,PMIN,VMAX,IFL)
       implicit none
!
!   Revised on 9/24/2005 to fix convergence problems at high pressure
!
!   ***   This subroutine calculates the maximum wind speed        ***
!   ***             and mimimum central pressure                   ***
!   ***    achievable in tropical cyclones, given a sounding       ***
!   ***             and a sea surface temperature.                 ***
!
!  INPUT:   SST: Sea surface temperature in K
!
!           PSL: Sea level pressure (Pa)
!
!           P,T,R: One-dimensional arrays of dimension NA
!             containing pressure (mb), temperature (K),
!             and mixing ratio (kg/kg). The arrays MUST be
!             arranged so that the lowest index corresponds
!             to the lowest model level, with increasing index
!             corresponding to decreasing pressure. The temperature
!             sounding should extend to at least the tropopause and
!             preferably to the lower stratosphere, however the
!             mixing ratios are not important above the boundary
!             layer. Missing mixing ratios can be replaced by zeros.
!
!           NA: The dimension of P,T and R
!
!           N:  The actual number of points in the sounding
!                (N is less than or equal to NA)
!
!  OUTPUT:  PMIN is the minimum central pressure, in mb
!
!           VMAX is the maximum surface wind speed, in m/s
!                  (reduced to reflect surface drag)
!
!
!           IFL is a flag: A value of 1 means OK; a value of 0
!              indicates no convergence (hypercane); a value of 2
!              means that the CAPE routine failed.
!
!-----------------------------------------------------------------------------
       integer, intent(in)  :: NA, IMAX, JMAX
       real(4), intent(in)  :: SST_IN(IMAX,JMAX),PSL_IN(IMAX,JMAX),P(NA),T_IN(NA,IMAX,JMAX),R_IN(NA,IMAX,JMAX)
       real(4), intent(out) :: PMIN(IMAX,JMAX),VMAX(IMAX,JMAX)
       integer, intent(out) :: IFL
       real(4) :: SST, PSL, T(NA), R(NA)
       integer :: i,j,k, N
       N=NA
       do i=1,IMAX
        do j=1,JMAX
          PSL=PSL_IN(i,j)/100.0 !Pa=>mb
          SST=SST_IN(i,j)-273.15
          do k=1,NA
             T(k)=T_IN(k,i,j)-273.15
             R(k)=R_IN(k,i,j)*1000.0
          enddo
          call PCMIN(SST,PSL,P,T,R,NA,N,PMIN(i,j),VMAX(i,j),IFL)
        enddo
       enddo
       RETURN
       END
       SUBROUTINE PCMIN3_KFLAG(SST_IN,PSL_IN,P,T_IN,R_IN,IMAX,JMAX,NA,PMIN,VMAX,IFL)
       implicit none
!
!   Revised on 9/24/2005 to fix convergence problems at high pressure
!
!   ***   This subroutine calculates the maximum wind speed        ***
!   ***             and mimimum central pressure                   ***
!   ***    achievable in tropical cyclones, given a sounding       ***
!   ***             and a sea surface temperature.                 ***
!
!  INPUT:   SST: Sea surface temperature in K
!
!           PSL: Sea level pressure (Pa)
!
!           P,T,R: One-dimensional arrays of dimension NA
!             containing pressure (mb), temperature (K),
!             and mixing ratio (kg/kg). The arrays MUST be
!             arranged so that the lowest index corresponds
!             to the lowest model level, with increasing index
!             corresponding to decreasing pressure. The temperature
!             sounding should extend to at least the tropopause and
!             preferably to the lower stratosphere, however the
!             mixing ratios are not important above the boundary
!             layer. Missing mixing ratios can be replaced by zeros.
!
!           NA: The dimension of P,T and R
!
!           N:  The actual number of points in the sounding
!                (N is less than or equal to NA)
!
!  OUTPUT:  PMIN is the minimum central pressure, in mb
!
!           VMAX is the maximum surface wind speed, in m/s
!                  (reduced to reflect surface drag)
!
!
!           IFL is a flag: A value of 1 means OK; a value of 0
!              indicates no convergence (hypercane); a value of 2
!              means that the CAPE routine failed.
!
!-----------------------------------------------------------------------------
       integer, intent(in)  :: NA, IMAX, JMAX
       real(4), intent(in)  :: SST_IN(IMAX,JMAX),PSL_IN(IMAX,JMAX),P(NA),T_IN(NA,IMAX,JMAX),R_IN(NA,IMAX,JMAX)
       real(4), intent(out) :: PMIN(IMAX,JMAX),VMAX(IMAX,JMAX)
       integer, intent(out) :: IFL
       real(4) :: SST, PSL, T(NA), R(NA)
       real(4), parameter :: UNDEF=-9.99E+33
       integer :: i,j,k, N, NS
       N=NA
       do i=1,IMAX
        do j=1,JMAX
          PSL=PSL_IN(i,j)/100.0 !Pa=>mb
          SST=SST_IN(i,j)-273.15
          NS=NA
          do k=1,NA
             if (T_IN(k,i,j) .ne. UNDEF .and. R_IN(k,i,j) .ne. UNDEF) then
               NS=min(NS,K)
               T(k)=T_IN(k,i,j)-273.15
               R(k)=R_IN(k,i,j)*1000.0
             else
               T(k)=UNDEF
               R(k)=UNDEF
             endif
          enddo
         !print *, "I,J=",i,j, "T=",T_IN(:,i,j), "R=",R_IN(:,i,j), "PSL=",PSL, "SST=",SST
         !call PCMIN(SST,PSL,P,T,R,NA,N,PMIN(i,j),VMAX(i,j),IFL)
          call PCMIN(SST,PSL,P(NS:NA),T(NS:NA),R(NS:NA),NA-NS+1,NA-NS+1,PMIN(i,j),VMAX(i,j),IFL)
        enddo
       enddo
       RETURN
       END
       SUBROUTINE PCMIN2(SST_IN,PSL_IN,P,T_IN,R_IN,NA,N,PMIN,VMAX,IFL)
       implicit none
!
!   Revised on 9/24/2005 to fix convergence problems at high pressure
!
!   ***   This subroutine calculates the maximum wind speed        ***
!   ***             and mimimum central pressure                   ***
!   ***    achievable in tropical cyclones, given a sounding       ***
!   ***             and a sea surface temperature.                 ***
!
!  INPUT:   SST: Sea surface temperature in K
!
!           PSL: Sea level pressure (Pa)
!
!           P,T,R: One-dimensional arrays of dimension NA
!             containing pressure (mb), temperature (K),
!             and mixing ratio (kg/kg). The arrays MUST be
!             arranged so that the lowest index corresponds
!             to the lowest model level, with increasing index
!             corresponding to decreasing pressure. The temperature
!             sounding should extend to at least the tropopause and
!             preferably to the lower stratosphere, however the
!             mixing ratios are not important above the boundary
!             layer. Missing mixing ratios can be replaced by zeros.
!
!           NA: The dimension of P,T and R
!
!           N:  The actual number of points in the sounding
!                (N is less than or equal to NA)
!
!  OUTPUT:  PMIN is the minimum central pressure, in mb
!
!           VMAX is the maximum surface wind speed, in m/s
!                  (reduced to reflect surface drag)
!
!
!           IFL is a flag: A value of 1 means OK; a value of 0
!              indicates no convergence (hypercane); a value of 2
!              means that the CAPE routine failed.
!
!-----------------------------------------------------------------------------
       integer, intent(in)  :: NA,N
       real(4), intent(in)  :: SST_IN,PSL_IN,P(NA),T_IN(NA),R_IN(NA)
       real(4), intent(out) :: PMIN,VMAX
       integer, intent(out) :: IFL
       real(4) :: SST, PSL, T(NA), R(NA)
       integer :: k
       PSL=PSL_IN/100.0 !Pa=>mb
       SST=SST_IN-273.15
       do k=1,NA
          T(k)=T_IN(k)-273.15
          R(k)=R_IN(k)*1000.0
       enddo
       call PCMIN(SST,PSL,P,T,R,NA,N,PMIN,VMAX,IFL)
       RETURN
       END
!
       SUBROUTINE PCMIN(SST,PSL,P,T_IN,R_IN,NA,N,PMIN,VMAX,IFL)
       implicit none
!
!   Revised on 9/24/2005 to fix convergence problems at high pressure
!
!   ***   This subroutine calculates the maximum wind speed        ***
!   ***             and mimimum central pressure                   ***
!   ***    achievable in tropical cyclones, given a sounding       ***
!   ***             and a sea surface temperature.                 ***
!
!  INPUT:   SST: Sea surface temperature in C
!
!           PSL: Sea level pressure (mb)
!
!           P,T,R: One-dimensional arrays of dimension NA
!             containing pressure (mb), temperature (C),
!             and mixing ratio (g/kg). The arrays MUST be
!             arranged so that the lowest index corresponds
!             to the lowest model level, with increasing index
!             corresponding to decreasing pressure. The temperature
!             sounding should extend to at least the tropopause and
!             preferably to the lower stratosphere, however the
!             mixing ratios are not important above the boundary
!             layer. Missing mixing ratios can be replaced by zeros.
!
!           NA: The dimension of P,T and R
!
!           N:  The actual number of points in the sounding
!                (N is less than or equal to NA)
!
!  OUTPUT:  PMIN is the minimum central pressure, in mb
!
!           VMAX is the maximum surface wind speed, in m/s
!                  (reduced to reflect surface drag)
!
!
!           IFL is a flag: A value of 1 means OK; a value of 0
!              indicates no convergence (hypercane); a value of 2
!              means that the CAPE routine failed.
!
!-----------------------------------------------------------------------------
        integer, intent(in)  :: NA,N
        real(4), intent(in)  :: SST,PSL,P(NA),T_IN(NA),R_IN(NA)
        real(4), intent(out) :: PMIN,VMAX
        integer, intent(out) :: IFL
        real,parameter :: UNDEF=-9.99E33
!
        real(4) :: CAPEA
        real(4) :: CKCD,SIG,b,VREDUC,SSTK,ES0,T(NA),R(NA),PM,TP,RP,PP,TOA,CAPEM,CAPEMS,TOM,TOMS
        real(4) :: RAT,RS0,TV1,TVAV,CAT,PNEW,CATFAC,FAC,S
        integer :: IDISS,NK,I,NP,IFLAG
!
!   *** Initialize
        PMIN=UNDEF
        VMAX=UNDEF
        CAPEA=UNDEF
!
!   ***   Adjustable constant: Ratio of C_k to C_D    ***
!
	CKCD=0.9
!
!   ***   Adjustable constant for buoyancy of displaced parcels:  ***
!   ***    0=Reversible ascent;  1=Pseudo-adiabatic ascent        ***
!
      SIG=0.0
!
!   ***  Adjustable switch: if IDISS = 0, no dissipative heating is   ***
!   ***     allowed; otherwise, it is                                 ***
!
	IDISS=1
!
!   ***  Exponent, b, in assumed profile of azimuthal velocity in eye,   ***
!   ***   V=V_m(r/r_m)^b. Used only in calculation of central pressure   ***
!
	b=2.0
!
!   *** Set level from which parcels lifted   ***
!
	NK=1
!
!   *** Factor to reduce gradient wind to 10 m wind
!
	VREDUC=0.8
!
!   ***   Normalize certain quantities   ***
!
	SSTK=SST+273.15
	ES0=6.112*EXP(17.67*SST/(243.5+SST))
	DO 40 I=1,N
	 R(I)=R_IN(I)*0.001
	 T(I)=T_IN(I)+273.15
   40	CONTINUE
!
!   ***   Default values   ***
!
      VMAX=0.0
	PMIN=PSL
	IFL=1
!
	NP=0
	PM=950.0
!
!   ***   Find environmental CAPE ***
!
      TP=T(NK)
      RP=R(NK)
      PP=P(NK)
      CALL CAPE(TP,RP,PP,T,R,P,NA,N,SIG,CAPEA,TOA,IFLAG)
     !write(6,*) "ENV",TP,RP,PP,CAPEA,IFLAG
     !do I=1,NA
     ! write(6,*) "   ",I,T(I),R(I),P(I)
     !enddo
      IF(IFLAG.NE.1)IFL=2
!
!   ***   Begin iteration to find mimimum pressure   ***
!
  100 CONTINUE
!
!   ***  Find CAPE at radius of maximum winds   ***
!
      TP=T(NK)
      PP=MIN(PM,1000.0)
      RP=0.622*R(NK)*PSL/(PP*(0.622+R(NK))-R(NK)*PSL)
      CALL CAPE(TP,RP,PP,T,R,P,NA,N,SIG,CAPEM,TOM,IFLAG)
      IF(IFLAG.NE.1)IFL=2
      RAT=SSTK/TOM
      IF(IDISS.EQ.0)RAT=1.0
!
!  ***  Find saturation CAPE at radius of maximum winds   ***
!
      TP=SSTK
      PP=MIN(PM,1000.0)
      RP=0.622*ES0/(PP-ES0)
      CALL CAPE(TP,RP,PP,T,R,P,NA,N,SIG,CAPEMS,TOMS,IFLAG)
      IF(IFLAG.NE.1)IFL=2
!
!  ***  Initial estimate of minimum pressure   ***
!
      RS0=RP
      TV1=T(1)*(1.+R(1)/0.622)/(1.+R(1))
	TVAV=0.5*(TV1+SSTK*(1.+RS0/0.622)/(1.+RS0))
!	CAT=0.5*CKCD*RAT*(CAPEMS-CAPEM)
	CAT=CAPEM-CAPEA+0.5*CKCD*RAT*(CAPEMS-CAPEM)
	CAT=MAX(CAT,0.0)
	PNEW=PSL*EXP(-CAT/(287.04*TVAV))
!
!   ***  Test for convergence   ***
!
	IF(ABS(PNEW-PM).GT.0.2)THEN
	 PM=PNEW
	 NP=NP+1
	 IF(NP.GT.1000.OR.PM.LT.400.0)THEN
	  PMIN=PSL
	  IFL=0
	  GOTO 900
	 END IF
	 GOTO 100
	ELSE
	 CATFAC=0.5*(1.+1./b)
!	 CAT=CKCD*RAT*CATFAC*(CAPEMS-CAPEM)
	 CAT=CAPEM-CAPEA+CKCD*RAT*CATFAC*(CAPEMS-CAPEM)
	 CAT=MAX(CAT,0.0)
	 PMIN=PSL*EXP(-CAT/(287.04*TVAV))
	END IF
  900	CONTINUE
	FAC=MAX(0.0,(CAPEMS-CAPEM))
	VMAX=VREDUC*SQRT(CKCD*RAT*FAC)
!
!   ***  Renormalize sounding arrays   ***
!
	DO 910 I=1,N
	 R(I)=R(I)*1000.0
	 T(I)=T(I)-273.15
  910	CONTINUE
!
	RETURN
	END
!
      SUBROUTINE CAPE(TP,RP,PP,T,R,P,ND,N,SIG,CAPED,TOB,IFLAG)
      implicit none
!
!     This subroutine calculates the CAPE of a parcel with pressure PP (mb),
!       temperature TP (K) and mixing ratio RP (gm/gm) given a sounding
!       of temperature (T in K) and mixing ratio (R in gm/gm) as a function
!       of pressure (P in mb). ND is the dimension of the arrays T,R and P,
!       while N is the actual number of points in the sounding. CAPED is
!       the calculated value of CAPE and TOB is the temperature at the
!       level of neutral buoyancy.  IFLAG is a flag
!       integer. If IFLAG = 1, routine is successful; if it is 0, routine did
!       not run owing to improper sounding (e.g.no water vapor at parcel level).
!       IFLAG=2 indicates that routine did not converge.
!
      REAL(4),intent(in):: TP,RP,PP,T(ND),R(ND),P(ND),SIG
      INTEGER(4),intent(in) :: ND,N
      INTEGER(4),intent(out) :: IFLAG
      REAL(4),   intent(out) :: CAPED
      REAL(4) ::  NA, TVRDIF(100), TOB, CPD,CPV,CL,CPVMCL,RV,RD,EPS,ALV0,TPC,ESP,EVP,RH,ALV,S
      REAL(4) ::  CHI,PLCL,TG,RG,TLVR,TJC,ES,SL,EM,SG,AP,TGNEW,TC,ENEW,RMEAN,PA,PFAC,PMA,PAT,PINB
      INTEGER(4) :: NCMAX,J,JMIN,NC,INB
!
!   ***   Default values   ***
!
      CAPED=0.0
      TOB=T(1)
      IFLAG=1
!
!   ***   Check that sounding is suitable    ***
!
      IF(RP.LT.1.0E-6.OR.TP.LT.200.0)THEN
       IFLAG=0
       RETURN
      END IF
!
!   ***   Assign values of thermodynamic constants     ***
!
      CPD=1005.7
      CPV=1870.0
!      CL=4190.0
      CL=2500.0
      CPVMCL=CPV-CL
      RV=461.5
      RD=287.04
      EPS=RD/RV
      ALV0=2.501E6
!
!   ***  Define various parcel quantities, including reversible   ***
!   ***                       entropy, S.                         ***
!
      TPC=TP-273.15
      ESP=6.112*EXP(17.67*TPC/(243.5+TPC))
      EVP=RP*PP/(EPS+RP)
      RH=EVP/ESP
	RH=MIN(RH,1.0)
      ALV=ALV0+CPVMCL*TPC
      S=(CPD+RP*CL)*LOG(TP)-RD*LOG(PP-EVP)+ &
     &   ALV*RP/TP-RP*RV*LOG(RH)
!
!   ***  Find lifted condensation pressure, PLCL   ***
!
	CHI=TP/(1669.0-122.0*RH-TP)
	PLCL=PP*(RH**CHI)
!
!   ***  Begin updraft loop   ***
!
	NCMAX=0
	DO J=1,N
	 TVRDIF(J)=0.0
	END DO
!
	JMIN=1E6
	DO 200 J=1,N
!
!    ***   Don't bother lifting parcel above 60 mb and skip sections of sounding below parcel level  ***
!
      IF(P(J).LT.59.0.OR.P(J).GE.PP)GOTO 200
!
	JMIN=MIN(JMIN,J)
!
!    ***  Parcel quantities below lifted condensation level   ***
!
	 IF(P(J).GE.PLCL)THEN
	  TG=TP*(P(J)/PP)**(RD/CPD)
	  RG=RP
!
!   ***   Calculate buoyancy   ***
!
	  TLVR=TG*(1.+RG/EPS)/(1.+RG)
	  TVRDIF(J)=TLVR-T(J)*(1.+R(J)/EPS)/(1.+R(J))
	 ELSE
!
!   ***  Parcel quantities above lifted condensation level  ***
!
	  TG=T(J)
	  TJC=T(J)-273.15
	  ES=6.112*EXP(17.67*TJC/(243.5+TJC))
	  RG=EPS*ES/(P(J)-ES)
!
!   ***  Iteratively calculate lifted parcel temperature and mixing   ***
!   ***                ratio for reversible ascent                    ***
!
	  NC=0
  120	  CONTINUE
	  NC=NC+1
!
!   ***  Calculate estimates of the rates of change of the entropy    ***
!   ***           with temperature at constant pressure               ***
!
	  ALV=ALV0+CPVMCL*(TG-273.15)
	  SL=(CPD+RP*CL+ALV*ALV*RG/(RV*TG*TG))/TG
	  EM=RG*P(J)/(EPS+RG)
	  SG=(CPD+RP*CL)*LOG(TG)-RD*LOG(P(J)-EM)+  &
     &      ALV*RG/TG
	  IF(NC.LT.3)THEN
	   AP=0.3
	  ELSE
	   AP=1.0
	  END IF
	  TGNEW=TG+AP*(S-SG)/SL
!
!   ***   Test for convergence   ***
!
	  IF(ABS(TGNEW-TG).GT.0.001)THEN
	   TG=TGNEW
	   TC=TG-273.15
	   ENEW=6.112*EXP(17.67*TC/(243.5+TC))
!
!   ***   Bail out if things get out of hand   ***
!
	   IF(NC.GT.500.OR.ENEW.GT.(P(J)-1.0))THEN
            IFLAG=2
            RETURN
	   END IF
	   RG=EPS*ENEW/(P(J)-ENEW)
	   GOTO 120
	  END IF
	  NCMAX=MAX(NC,NCMAX)
!
!   *** Calculate buoyancy   ***
!
        RMEAN=SIG*RG+(1.-SIG)*RP
	  TLVR=TG*(1.+RG/EPS)/(1.+RMEAN)
	  TVRDIF(J)=TLVR-T(J)*(1.+R(J)/EPS)/(1.+R(J))
	 END IF
  200	CONTINUE
!
!  ***  Begin loop to find NA, PA, and CAPE from reversible ascent ***
!
	NA=0.0
	PA=0.0
!
!   ***  Find maximum level of positive buoyancy, INB    ***
!
	INB=1
	DO 550 J=N,JMIN,-1
	 IF(TVRDIF(J).GT.0.0)INB=MAX(INB,J)
  550	CONTINUE
	IF(INB.EQ.1)RETURN
!
!   ***  Find positive and negative areas and CAPE  ***
!
	IF(INB.GT.1)THEN
	 DO 600 J=JMIN+1,INB
	  PFAC=RD*(TVRDIF(J)+TVRDIF(J-1))*(P(J-1)-P(J))/(P(J)+P(J-1))
	  PA=PA+MAX(PFAC,0.0)
	  NA=NA-MIN(PFAC,0.0)
  600	 CONTINUE
!
!   ***   Find area between parcel pressure and first level above it ***
!
	PMA=(PP+P(JMIN))
	PFAC=RD*(PP-P(JMIN))/PMA
	PA=PA+PFAC*MAX(TVRDIF(JMIN),0.0)
	NA=NA-PFAC*MIN(TVRDIF(JMIN),0.0)
!
!   ***   Find residual positive area above INB and TO  ***
!
       PAT=0.0
       TOB=T(INB)
       IF(INB.LT.N)THEN
        PINB=(P(INB+1)*TVRDIF(INB)-P(INB)*TVRDIF(INB+1))/  &
     &   (TVRDIF(INB)-TVRDIF(INB+1))
        PAT=RD*TVRDIF(INB)*(P(INB)-PINB)/(P(INB)+PINB)
	  TOB=(T(INB)*(PINB-P(INB+1))+T(INB+1)*(P(INB)-PINB))/  &
     &    (P(INB)-P(INB+1))
       END IF
!
!   ***   Find CAPE  ***
!
	 CAPED=PA+PAT-NA
	 CAPED=MAX(CAPED,0.0)
	END IF
!
	RETURN
	END


text:
c     Subroutine to process all lines in single file
c     AJ_Kettle, 24Feb2019
c     13Mar2020: modified to read combined USAF files

c     Subroutine structure:
c     -get_params_singleline2
c     -find_minmax_string_vec
c     -find_nelem_present
c     -convert_variables2
c       -convert_single_string_to_float

      SUBROUTINE process_single_file5(
     +   s_dir_usaf_filelist,s_single_stnlist,
     +   f_ndflag,
     +   l_rgh_datalines,l_rgh_char,
     +   i_linelength_min,i_linelength_max,
     +   i_commacnt_min,i_commacnt_max,

     +   l_data,l_c1,l_c2,
     +   s_vec_platformid,s_vec_networktype,
     +   s_vec_ncdc_ob_time,s_vec_reporttypecode,
     +   s_vec_latitude,s_vec_longitude,s_vec_platformheight,
     +   s_vec_winddirection_deg,s_vec_winddirection_qc,
     +   s_vec_windspeed_ms,s_vec_windspeed_qc,
     +   s_vec_airtemperature_c,s_vec_airtemperature_qc,
     +   s_vec_dewpointtemperature_c,s_vec_dewpointtemperature_qc,
     +   s_vec_sealevelpressure_hpa,s_vec_sealevelpressure_qc,
     +   s_vec_stationpressure_hpa,s_vec_stationpressure_qc,
     +   f_vec_latitude,f_vec_longitude,f_vec_platformheight,
     +   f_vec_windspeed_ms,
     +   f_vec_sealevelpressure_hpa,f_vec_stationpressure_hpa)

      IMPLICIT NONE
c************************************************************************
c     Declare variables passed into program

      CHARACTER(LEN=*)    :: s_dir_usaf_filelist
      CHARACTER(LEN=*)    :: s_single_stnlist

      REAL                :: f_ndflag

      INTEGER             :: l_rgh_datalines
      INTEGER             :: l_rgh_char

      INTEGER             :: i_linelength_min
      INTEGER             :: i_linelength_max
      INTEGER             :: i_commacnt_min
      INTEGER             :: i_commacnt_max
c*****
c     Variables used within program

      INTEGER             :: i,j,k,ii,jj,kk
      INTEGER             :: io

      INTEGER             :: ii_top
      INTEGER             :: ii_body
      INTEGER             :: ii_data
      INTEGER             :: i_flag_body

      INTEGER             :: i_flag_sep
      INTEGER             :: i_index_sep

      CHARACTER(l_rgh_char):: s_linget
      CHARACTER(l_rgh_char):: s_header

      CHARACTER(300)::    s_pathandname

      INTEGER             :: i_linelength

      INTEGER             :: i_commacnt
      INTEGER             :: i_comma_pos(200)

      INTEGER             :: l_data
      INTEGER             :: l_c1
      INTEGER             :: l_c2

      CHARACTER(LEN=20)   :: s_platformid
      CHARACTER(LEN=20)   :: s_networktype
      CHARACTER(LEN=14)   :: s_ncdc_ob_time
      CHARACTER(LEN=6)    :: s_reporttypecode
      CHARACTER(LEN=15)   :: s_latitude
      CHARACTER(LEN=15)   :: s_longitude
      CHARACTER(LEN=15)   :: s_platformheight

      CHARACTER(LEN=l_c1) :: s_winddirection_deg      !10
      CHARACTER(LEN=l_c2) :: s_winddirection_qc       !10
      CHARACTER(LEN=l_c1) :: s_windspeed_ms           !10
      CHARACTER(LEN=l_c2) :: s_windspeed_qc           !10
      CHARACTER(LEN=l_c1) :: s_airtemperature_c       !10
      CHARACTER(LEN=l_c2) :: s_airtemperature_qc      !10
      CHARACTER(LEN=l_c1) :: s_dewpointtemperature_c  !10
      CHARACTER(LEN=l_c2) :: s_dewpointtemperature_qc !10
      CHARACTER(LEN=l_c1) :: s_sealevelpressure_hpa   !10
      CHARACTER(LEN=l_c2) :: s_sealevelpressure_qc    !10
      CHARACTER(LEN=l_c1) :: s_stationpressure_hpa    !10
      CHARACTER(LEN=l_c2) :: s_stationpressure_qc     !10

      CHARACTER(LEN=20)   :: s_vec_platformid(l_rgh_datalines)
      CHARACTER(LEN=20)   :: s_vec_networktype(l_rgh_datalines)
      CHARACTER(LEN=14)   :: s_vec_ncdc_ob_time(l_rgh_datalines)
      CHARACTER(LEN=6)    :: s_vec_reporttypecode(l_rgh_datalines)
      CHARACTER(LEN=15)   :: s_vec_latitude(l_rgh_datalines)
      CHARACTER(LEN=15)   :: s_vec_longitude(l_rgh_datalines)
      CHARACTER(LEN=15)   :: s_vec_platformheight(l_rgh_datalines)
      CHARACTER(LEN=l_c1) :: s_vec_winddirection_deg(l_rgh_datalines)
      CHARACTER(LEN=l_c2) :: s_vec_winddirection_qc(l_rgh_datalines)
      CHARACTER(LEN=l_c1) :: s_vec_windspeed_ms(l_rgh_datalines)
      CHARACTER(LEN=l_c2) :: s_vec_windspeed_qc(l_rgh_datalines)
      CHARACTER(LEN=l_c1) :: s_vec_airtemperature_c(l_rgh_datalines)
      CHARACTER(LEN=l_c2) :: s_vec_airtemperature_qc(l_rgh_datalines)
      CHARACTER(LEN=l_c1)::s_vec_dewpointtemperature_c(l_rgh_datalines)
      CHARACTER(LEN=l_c2)::s_vec_dewpointtemperature_qc(l_rgh_datalines)
      CHARACTER(LEN=l_c1) :: s_vec_sealevelpressure_hpa(l_rgh_datalines)
      CHARACTER(LEN=l_c2) :: s_vec_sealevelpressure_qc(l_rgh_datalines)
      CHARACTER(LEN=l_c1) :: s_vec_stationpressure_hpa(l_rgh_datalines)
      CHARACTER(LEN=l_c2) :: s_vec_stationpressure_qc(l_rgh_datalines)

      REAL                :: f_vec_latitude(l_rgh_datalines)
      REAL                :: f_vec_longitude(l_rgh_datalines)
      REAL                :: f_vec_platformheight(l_rgh_datalines)
      REAL                :: f_vec_winddirection_deg(l_rgh_datalines)
      REAL                :: f_vec_windspeed_ms(l_rgh_datalines)
      REAL                :: f_vec_airtemperature_c(l_rgh_datalines)
      REAL                ::f_vec_dewpointtemperature_c(l_rgh_datalines)
      REAL                :: f_vec_sealevelpressure_hpa(l_rgh_datalines)
      REAL                :: f_vec_stationpressure_hpa(l_rgh_datalines)

      INTEGER             :: i_minmax_strlen_platformid(2)
      INTEGER             :: i_minmax_strlen_networktype(2)
      INTEGER             :: i_minmax_strlen_ncdc_ob_time(2)
      INTEGER             :: i_minmax_strlen_reporttypecode(2)
      INTEGER             :: i_minmax_strlen_latitude(2)
      INTEGER             :: i_minmax_strlen_longitude(2)
      INTEGER             :: i_minmax_strlen_platformheight(2)
      INTEGER             :: i_minmax_strlen_winddirection(2)
      INTEGER             :: i_minmax_strlen_windspeed(2)
      INTEGER             :: i_minmax_strlen_sealevelpressure(2)
      INTEGER             :: i_minmax_strlen_stationpressure(2)

      INTEGER             :: i_cnt_strlen_platformid
      INTEGER             :: i_cnt_strlen_networktype
      INTEGER             :: i_cnt_strlen_ncdc_ob_time
      INTEGER             :: i_cnt_strlen_reporttypecode
      INTEGER             :: i_cnt_strlen_latitude
      INTEGER             :: i_cnt_strlen_longitude
      INTEGER             :: i_cnt_strlen_platformheight
      INTEGER             :: i_cnt_strlen_windspeed
      INTEGER             :: i_cnt_strlen_sealevelpressure
      INTEGER             :: i_cnt_strlen_stationpressure
c************************************************************************
c      print*,'just entered process_single_file4'

c      print*,'s_fileuse=',TRIM(s_fileuse)

      ii=0
      ii_top=0
      ii_body=0
      ii_data=0
      i_flag_body=0

      i_flag_sep=0

      i_linelength_min=5000
      i_linelength_max=-5000
      i_commacnt_min=5000
      i_commacnt_max=-5000

      s_pathandname=TRIM(s_dir_usaf_filelist)//TRIM(s_single_stnlist)

      OPEN(UNIT=1,FILE=TRIM(s_pathandname),FORM='formatted',
     +  STATUS='OLD',ACTION='READ')

c     Read in header
      READ(1,1002,IOSTAT=io) s_header

      DO

c      Read data line
       READ(1,1002,IOSTAT=io) s_linget
1002   FORMAT(a4000)
       IF (io .GT. 0) THEN
        WRITE(*,*) 'Check input. Something went wrong'
        GOTO 100
       ELSE IF (io .LT. 0) THEN
c        WRITE(*,*) 'end of file reached'
        GOTO 100
       ELSE

        ii=ii+1

cc       Augment running index of top information
c        IF (i_flag_sep.EQ.0) THEN
c         ii_top=ii_top+1
c        ENDIF

cc       Augment running index of body after separator found
c        IF (i_flag_sep.EQ.1) THEN
c         ii_body=ii_body+1
c        ENDIF

cc       Print out header line
c        IF (ii_body.EQ.1) THEN
c         s_header=s_linget
cc         print*,'header line=',TRIM(s_linget)
cc         print*,LEN_TRIM(s_linget)
c        ENDIF

cc       Query data lines
c        IF (ii_body.GE.2) THEN
         ii_data=ii_data+1

c        Test if data lines are within storage
         IF (ii_data.GT.l_rgh_datalines) THEN
          print*,'emergency stop; data lines exceed storage'
          STOP 'process_single_file5'
         ENDIF

         i_linelength=LEN_TRIM(s_linget)
         i_linelength_min=MIN(i_linelength,i_linelength_min)
         i_linelength_max=MAX(i_linelength,i_linelength_max)

c        Find number commas in lines
         i_commacnt=0
         DO j=1,i_linelength
          IF (s_linget(j:j).EQ.',') THEN
           i_commacnt=i_commacnt+1
           i_comma_pos(i_commacnt)=j

c          Condito for crashing program: too many commas
           IF (i_commacnt.GT.199) THEN
            print*,'too many commas'
            STOP 'process_single_file5'
           ENDIF
          ENDIF
         ENDDO

         i_commacnt_min=MIN(i_commacnt,i_commacnt_min)
         i_commacnt_max=MAX(i_commacnt,i_commacnt_max)

c         print*,'ii=',ii
c         print*,'i_linelength minmax',i_linelength_min,i_linelength_max
c         print*,'i_commacnt minmax',  i_commacnt_min,i_commacnt_max
c         print*,'s_linget=',TRIM(s_linget)
c         print*,'i_linelength=',i_linelength

c        Export bad line if number of commas is too few
         IF (i_commacnt.LT.164) THEN
          print*,'ii=',ii
          print*,'s_linget=',TRIM(s_linget)
          print*,'i_linelength=',i_linelength
          print*,'i_commacnt=',i_commacnt
          print*,'number commas below threshold 164 stnpress'

c         16Jul2019: removed to avoid problems with parallel run
c          CALL export_bad_line(ii,
c     +      s_directory_outfile_generalstats,s_subdirectory,
c     +      i_linelength,i_commacnt,
c     +      l_rgh_char,s_linget)

c          STOP 'process_single_file4'
         ENDIF

c        Get parameters from single line
c         print*,'ii=',ii
c         print*,'s_linget=',TRIM(s_linget)

         CALL get_params_singleline20210527(
     +     l_rgh_char,s_linget,
     +     i_commacnt,i_comma_pos,
     +     l_c1,l_c2,
     +     s_platformid,s_networktype,s_ncdc_ob_time,
     +     s_reporttypecode,s_latitude,s_longitude,s_platformheight,
     +     s_winddirection_deg,s_winddirection_qc,
     +     s_windspeed_ms,s_windspeed_qc,
     +     s_airtemperature_c,s_airtemperature_qc,
     +     s_dewpointtemperature_c,s_dewpointtemperature_qc,
     +     s_sealevelpressure_hpa,s_sealevelpressure_qc,
     +     s_stationpressure_hpa,s_stationpressure_qc)

c        Crash program if lat or lon fields empty
         IF (LEN_TRIM(s_latitude).EQ.0.OR.LEN_TRIM(s_longitude).EQ.0)
     +      THEN
          print*,'emergency stop, lat or lon blank'
          print*,'s_latitude=',TRIM(s_latitude)
          print*,'s_longitude=',TRIM(s_longitude)
          STOP 'process_single_file5'
         ENDIF

c        Archive strings from single line
         s_vec_platformid(ii_data)          =s_platformid
         s_vec_networktype(ii_data)         =s_networktype
         s_vec_ncdc_ob_time(ii_data)        =s_ncdc_ob_time
         s_vec_reporttypecode(ii_data)      =s_reporttypecode
         s_vec_latitude(ii_data)            =s_latitude
         s_vec_longitude(ii_data)           =s_longitude
         s_vec_platformheight(ii_data)      =s_platformheight
         s_vec_winddirection_deg(ii_data)   =s_winddirection_deg
         s_vec_winddirection_qc(ii_data)    =s_winddirection_qc
         s_vec_windspeed_ms(ii_data)        =s_windspeed_ms
         s_vec_windspeed_qc(ii_data)        =s_windspeed_qc
         s_vec_airtemperature_c(ii_data)    =s_airtemperature_c
         s_vec_airtemperature_qc(ii_data)   =s_airtemperature_qc
         s_vec_dewpointtemperature_c(ii_data)=s_dewpointtemperature_c
         s_vec_dewpointtemperature_qc(ii_data)=s_dewpointtemperature_qc
         s_vec_sealevelpressure_hpa(ii_data)=s_sealevelpressure_hpa
         s_vec_sealevelpressure_qc(ii_data) =s_sealevelpressure_qc
         s_vec_stationpressure_hpa(ii_data) =s_stationpressure_hpa
         s_vec_stationpressure_qc(ii_data)  =s_stationpressure_qc

c         print*,'just after get_params_singleline2'

c         print*,'s_linget=',TRIM(s_linget)
c         print*,'ii_data',ii_data,
c     +     'latitude='//TRIM(s_latitude),
c     +     'longitude='//TRIM(s_longitude)

c         CALL SLEEP(1)

c         print*,'l_c1,l_c2=',l_c1,l_c2
c         print*,'cleared vector assignments'
c         STOP 'process_single_file4'

c        ENDIF

cc       Locate dashed separator
c        IF (s_linget(1:3).EQ.'---') THEN
c         i_flag_sep=1
c         i_index_sep=ii
c        ENDIF

       ENDIF
      ENDDO                !end of line-counting loop

100   CONTINUE

      CLOSE(UNIT=1)

      l_data=ii_data

c      print*,'ii=',ii
c      print*,'ii_top,ii_body=',ii_top,ii_body,ii_top+ii_body
      print*,'ii_data=',ii_data
c      print*,'i_index_sep=',i_index_sep

c      print*,'i_linelength_min/max',
c     +  i_linelength_min,i_linelength_max
c      print*,'i_commacnt_min/max',
c     +  i_commacnt_min,i_commacnt_max

c      STOP 'process_single_file4'
c************************************************************************
c************************************************************************
c     Condition for crashing program
c     27Mar2019 removed condition because it tripped at WMO06759 i=1588
      IF (i_commacnt_min.NE.i_commacnt_max) THEN
        print*,'i_commacnt_min/max',
     +     i_commacnt_min,i_commacnt_max
c        STOP 'process_single_file4'
      ENDIF
c*****
c     Call subroutine to find min/max of string vectors
      CALL find_minmax_string_vec(
     +  l_rgh_datalines,l_data,
     +  s_vec_platformid,s_vec_networktype,
     +  s_vec_ncdc_ob_time,s_vec_reporttypecode,
     +  s_vec_latitude,s_vec_longitude,s_vec_platformheight,
     +  s_vec_windspeed_ms,
     +  s_vec_sealevelpressure_hpa,s_vec_stationpressure_hpa,

     +  i_minmax_strlen_platformid,i_minmax_strlen_networktype,
     +  i_minmax_strlen_ncdc_ob_time,i_minmax_strlen_reporttypecode,
     +  i_minmax_strlen_latitude,i_minmax_strlen_longitude,
     +  i_minmax_strlen_platformheight,i_minmax_strlen_windspeed,
     +  i_minmax_strlen_sealevelpressure,
     +  i_minmax_strlen_stationpressure)
c*****
c     Find number elements represented

      CALL find_nelem_present(
     +  l_rgh_datalines,l_data,
     +  s_vec_platformid,s_vec_networktype,
     +  s_vec_ncdc_ob_time,s_vec_reporttypecode,
     +  s_vec_latitude,s_vec_longitude,s_vec_platformheight,
     +  s_vec_windspeed_ms,
     +  s_vec_sealevelpressure_hpa,s_vec_stationpressure_hpa,

     +  i_cnt_strlen_platformid,i_cnt_strlen_networktype,
     +  i_cnt_strlen_ncdc_ob_time,i_cnt_strlen_reporttypecode,
     +  i_cnt_strlen_latitude,i_cnt_strlen_longitude,
     +  i_cnt_strlen_platformheight,i_cnt_strlen_windspeed,
     +  i_cnt_strlen_sealevelpressure,i_cnt_strlen_stationpressure)

c*****
c     Convert parameters to float

      CALL convert_variables2(f_ndflag,
     +  l_rgh_datalines,l_data,
     +  s_vec_latitude,s_vec_longitude,s_vec_platformheight,
     +  s_vec_winddirection_deg,s_vec_windspeed_ms,
     +  s_vec_airtemperature_c,s_vec_dewpointtemperature_c,
     +  s_vec_sealevelpressure_hpa,s_vec_stationpressure_hpa,

     +  f_vec_latitude,f_vec_longitude,f_vec_platformheight,
     +  f_vec_winddirection_deg,f_vec_windspeed_ms,
     +  f_vec_airtemperature_c,f_vec_dewpointtemperature_c,
     +  f_vec_sealevelpressure_hpa,f_vec_stationpressure_hpa)

c     +  s_vec_windspeed_ms,
c     +  s_vec_sealevelpressure_hpa,s_vec_stationpressure_hpa,

c     +  f_vec_latitude,f_vec_longitude,f_vec_platformheight,
c     +  f_vec_windspeed_ms,
c     +  f_vec_sealevelpressure_hpa,f_vec_stationpressure_hpa)
c*****
c     Count number of filled elements

c      CALL count_filled_fields(f_ndflag,
c     +  l_rgh_datalines,l_data,
c     +  s_vec_latitude,s_vec_longitude,s_vec_platformheight,
c     +  s_vec_windspeed_ms,
c     +  s_vec_sealevelpressure_hpa,s_vec_stationpressure_hpa,

c*****
c      print*,'just leaving process_single_file4'

      RETURN
      END

ext: f

max_stars_repo_name: glamod/glamod-nuim

max_issues_repo_name: glamod/glamod-nuim

max_forks_repo_name: glamod/glamod-nuim


text:
!----- GPL ---------------------------------------------------------------------
!
!  Copyright (C)  Stichting Deltares, 2011-2016.
!
!  This program is free software: you can redistribute it and/or modify
!  it under the terms of the GNU General Public License as published by
!  the Free Software Foundation version 3.
!
!  This program is distributed in the hope that it will be useful,
!  but WITHOUT ANY WARRANTY; without even the implied warranty of
!  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!  GNU General Public License for more details.
!
!  You should have received a copy of the GNU General Public License
!  along with this program.  If not, see <http://www.gnu.org/licenses/>.
!
!  contact: delft3d.support@deltares.nl
!  Stichting Deltares
!  P.O. Box 177
!  2600 MH Delft, The Netherlands
!
!  All indications and logos of, and references to, "Delft3D" and "Deltares"
!  are registered trademarks of Stichting Deltares, and remain the property of
!  Stichting Deltares. All rights reserved.
!
!-------------------------------------------------------------------------------
!  $Id: wstmod.f90 6649 2016-10-12 12:09:03Z jeuke_ml $
!  $HeadURL: https://svn.oss.deltares.nl/repos/delft3d/tags/6686/src/utils_lgpl/hyd_to_waq_part_lib/packages/hyd_to_waq_part_lib/src/wstmod.f90 $

      module wstmod

      ! module contains everything for the wastloads discription
      ! created June 2004 by Jan van Beek

      implicit none

      integer, parameter, private               :: NAME_SIZE  =  256      ! size of descriptive names
      integer, parameter, private               :: MAX_NUM    =    5      ! allocated per bunch

      ! wateload types

      integer, parameter                        :: DLWQ_WASTE_NORMAL = 1  ! normal wasteload
      integer, parameter                        :: DLWQ_WASTE_INLET  = 2  ! inlet wasteload
      integer, parameter                        :: DLWQ_WASTE_OUTLET = 3  ! outlet wasteload
      integer, parameter                        :: DLWQ_WASTE_WALK   = 4  ! walking wasteload

      type t_wasteload
         character(len=NAME_SIZE)               :: name                   ! name of wasteload
         character(len=NAME_SIZE)               :: waqtype                ! waqtype of wasteload
         integer                                :: m                      ! m coordinate
         integer                                :: n                      ! n coordinate
         integer                                :: k                      ! k coordinate
         integer                                :: type                   ! type of wasteload
      end type t_wasteload

      type t_wasteload_coll
         type(t_wasteload), pointer             :: wasteload_pnts(:)      ! pointer to the wasteloads
         integer                                :: maxsize                ! maximum size of the current array
         integer                                :: cursize                ! filled up to this size
         integer                                :: no_flow                ! total number of flows
         logical                                :: l_seconds              ! if time is in seconds or ddhhmmss
      end type t_wasteload_coll

      contains

      ! function to find a wasteload name in a collection, case sensitive at the moment

      function wasteload_coll_find( wasteload_coll, name ) result ( iret )

         type(t_wasteload_coll)                 :: wasteload_coll         ! collection of wasteloads
         character(LEN=*)                       :: name                   ! name of wasteload to be found
         integer                                :: iret                   ! result index in collection or 0 if not found

         integer                                :: i                      ! loop counter

         iret = 0
         do i = 1 , wasteload_coll%cursize
            if ( wasteload_coll%wasteload_pnts(i)%name .eq. name ) then
               iret = i
               return
            endif
         end do

      end function wasteload_coll_find

      ! function to add to a collection of wasteloads (copy)

      function wasteload_coll_add( wasteload_coll , wasteload ) result ( cursize )

         type(t_wasteload_coll)                 :: wasteload_coll         ! collection of wasteloads
         type(t_wasteload)                      :: wasteload              ! wasteload to be added
         integer                                :: cursize                ! return value the new current collection size
                                                                          ! and the index of the added wasteload

         type(t_wasteload), pointer             :: wasteload_pnts(:)      ! pointer for the resize operation
         integer                                :: i                      ! loop counter

         if ( wasteload_coll%cursize .eq. wasteload_coll%maxsize ) then

            ! resize, allocate new array

            allocate ( wasteload_pnts ( wasteload_coll%maxsize + MAX_NUM ) )

            ! copy the wasteloads into the new array

            do i = 1 , wasteload_coll%maxsize
               wasteload_pnts(i) = wasteload_coll%wasteload_pnts(i)   ! copies the wasteloads
            enddo

            ! deallocate the old array and attach the new array to the collection

            if ( wasteload_coll%maxsize .ne. 0 ) deallocate ( wasteload_coll%wasteload_pnts )
            wasteload_coll%wasteload_pnts => wasteload_pnts
            wasteload_coll%maxsize = wasteload_coll%maxsize + MAX_NUM

         endif

         wasteload_coll%cursize = wasteload_coll%cursize + 1
         wasteload_coll%wasteload_pnts(wasteload_coll%cursize ) = wasteload
         cursize = wasteload_coll%cursize
         return

      end function wasteload_coll_add

      end module wstmod

max_stars_repo_name: liujiamingustc/phd

max_issues_repo_name: liujiamingustc/phd

max_forks_repo_name: liujiamingustc/phd


text:
: keydown   over #5 rshift #1 + #2 lshift + @ swap #31 and #1 swap lshift and 0<> ;
: klast     kblast keydown ;
: kstate    kbstate keydown ;
: kdelta    >r  r@ kstate #1 and  r> klast #1 and  - ;
: pressed   kdelta #1 = ;
: released  kdelta #-1 = ;

max_stars_repo_name: RogerLevy/RamenEngine

max_issues_repo_name: Alim-Oezdemir/RamenEngine

max_forks_repo_name: Alim-Oezdemir/RamenEngine


text:
	subroutine m70wti (fcb, func, time, button, x, y)
c
c     This routine waits for the appropriate interupt from the
c     Unibus M70 interface card, then returns button and cursor
c     information.
c
c
c     parameter descriptions:
c
c        fcb   is a system info. array.
c
c        func  is interpreted:
c                 0 ==> wait for button push.
c                 1 ==> wait for cursor move.
c                 2 ==> wait for button push or cursor move.
c
	integer fcb(*), func, time, button, x, y
c
	include 'fcbu.inc'
	external io$_rewindoff
	integer*4 sys$qiow, mask, chan, iosb(2)
	integer*2 chan2(2)
	equivalence (chan, chan2)
c
	chan2(1) = fcb(fcb_u_m70_chan)	! get M70 channel
	chan2(2) = fcb(fcb_u_m70_chan+1)
c
	if (func.eq.0) then
	  mask = '0400'x		! wait for button
	elseif (func.eq.1) then
	  mask = '0800'x		! wait for trackball
	elseif (func.eq.2) then
	  mask = '0C00'x		! wait for button or trackball
	else
	  mask = '0C00'x
	endif
c
	status = sys$qiow (, %val(chan), io$_rewindoff,
	1		iosb,,,%val(mask),,,,,)
c
c  Get button word and X-Y position of cursor
c
	call rbutn (fcb, button, x, y)
c
	return
	end

max_stars_repo_name: ysBach/irafdocgen

max_issues_repo_name: ysBach/irafdocgen

max_forks_repo_name: ysBach/irafdocgen


text:
!! Copyright (c) 1994 Unicomp, Inc.  All rights reserved.
!
!! Developed at Unicomp, Inc.
!
!! Permission to USE, copy, modify, and distribute this
!! software is freely granted, provided that this notice
!! is preserved.

!! NAG F90

module USR_SMlib_IO

   implicit none

!!  Default input and output units:

   integer, parameter :: DEFAULT_INPUT_UNIT = 5
   integer, parameter :: DEFAULT_OUTPUT_UNIT = 6

!!  Number and value of preconnected units

   integer, parameter :: NUMBER_OF_PRECONNECTED_UNITS = 3
   integer, parameter :: PRECONNECTED_UNITS (NUMBER_OF_PRECONNECTED_UNITS) = &
                         (/ 0, 5, 6 /)

!!  Values returned to IOSTAT for end of record and end of file

   integer, parameter :: END_OF_RECORD = -2
   integer, parameter :: END_OF_FILE = -1

!!  Largest allowed unit number (or a large number, if none)

   integer, parameter :: MAX_UNIT_NUMBER = 1000

contains

   function new_unit ()  result (result)

!!  Returns a unit number of a unit that exists and is not connected

      integer :: result
      logical :: exists, opened
      integer :: ios

      do result = 1, max_unit_number
         if (result == DEFAULT_INPUT_UNIT .or. &
             result == DEFAULT_OUTPUT_UNIT) cycle
         if (any (result == PRECONNECTED_UNITS)) cycle
         inquire (unit = result,  &
                  exist = exists,  &
                  opened = opened,  &
                  iostat = ios)
         if (exists .and. .not. opened .and. ios == 0) return
      end do

      result = -1

   end function new_unit

end module USR_SMlib_IO

max_stars_repo_name: wawiesel/tapack

max_issues_repo_name: wawiesel/tapack

max_forks_repo_name: wawiesel/tapack


text:
A shortstory about an alien in Davis, CA:

Originally from a small peaceful town in the warm Pacific coast state of http://www.history.com/topics/mexico/nayarit/ Nayarit, Mexico, Im told I was brought to the southern California coast to be united with family kicking and screaming when I was 3 years old. But being given no other choice, I eventually learned to adjust to city living in the States and flourished in the English language and adopted an American lifestyle while attending K12. As far as I knew, http://www.citydata.com/city/OxnardCalifornia.html/ Oxnard, CA had been and would always be my hometown. By senior year, I had participated in Honors courses, played Varsity Basketball, been an officer of a school club, and was more than ready to join fellow classmates in the process of choosing majors and applying to Universities. But I learned that was not quite what life had in store for me. By this time, the http://travel.state.gov/content/visas/english/visit/visitor.html/ B2 visitor visa with which my grandparents brought me here had expired and the http://travel.state.gov/content/visas/english/immigrate/types/family/familypreference.html#1/ F4 (I130) familybased immigrant visa petition entered by a US citizen brother on my behalf had been processing for 9 years with https://egov.uscis.gov/cris/processingTimesDisplay.do/ no given date of completion to look forward to. I was an illegal http://legaldictionary.thefreedictionary.com/alien/ alien. I was stifled by the knowledge that I could not at this point follow my high school dreams of attending any of the http://www.universityofcalifornia.edu/ UCs Id been admitted to or of accepting an invitation to pursue runway modeling abroad. So I did the only thing I could do I ran. I packed what I could along with a couple of scholarships, jumped in a car with a cousin, and thus arrived in http://www.citydata.com/city/ElkGroveCalifornia.html/ Elk Grove, CA two days after graduating from http://www.huenemehigh.us/index.php/ Hueneme High School in the summer of 2003. Needing some semblance of control in my own life, I soon moved on my own to http://www.citydata.com/city/SacramentoCalifornia.html/ Sacramento, CA where I thankfully found fulltime work as a server, went to school little by little at the local community college without any kind of financial aid, stayed out of trouble, and established my roots as an adult as best I could while living in the shadows for 4 years.

Early in 2007, I finally received word from http://www.uscis.gov/ USCIS and was granted an appointment to appear in http://photos.state.gov/libraries/ciudadjuarez/231771/PDFs/Packet4dec13.pdf/ Ciudad Juarez, Mexico for final processing of my immigrant visa petition. Hoping for the best, I was unprepared for the worst when I arrived I was informed that my visa petition had been denied on grounds of my residency in the US since my visitors visa expired when I was 9 years old. I was not permitted to cross the border back into the US and thus was exiled from the only home I knew. I was taken in by my godparents in http://www.visitmexico.com/meetings/toluca/ Toluca, Mexico and was surprised at how hard it was to transition into life in my country of origin. Bureaucratic red tape made it impossible to translate my high school transcripts with which to find employment on my own, so for the first 6 months my godparents invited me to tag along on family vacations throughout Mexico. Such places of archaeological interest as the https://www.google.com/search?qtenochtitlan&sourcelnms&tbmisch&saX&eiv9lrU9noLNGpyATiuYLICQ&sqi2&ved0CAYQ_AUoAQ&biw1366&bih753#facrc_&imgdii_&imgrcAZzLdqlfX1M95M%253A%3B6YYzGy0chTF3M%3Bhttp%253A%252F%252Fwww.solucionpolitica.net%252Fwpcontent%252Fuploads%252F2012%252F07%252Ftenochtitlan.jpg%3Bhttp%253A%252F%252Fwww.solucionpolitica.net%252Ffundacionmexicotenochtitlancelebrasu687aniversario%252F%3B1024%3B768/ Aztec pyramids of Tenochtitlan and those of natural beauty like the deep blue waterfalls of https://www.google.com/search?qtenochtitlan&sourcelnms&tbmisch&saX&eiv9lrU9noLNGpyATiuYLICQ&sqi2&ved0CAYQ_AUoAQ&biw1366&bih753#qPuente+de+Dios&tbmisch&facrc_&imgdii_&imgrczLLwBSdyOXwjuM%253A%3BWsz2fAR6HGK1qM%3Bhttp%253A%252F%252Fwww.tamasopo.com.mx%252FGaleria%252520nueva%252520de%252520Tamasopo.JPG%3Bhttp%253A%252F%252Fwww.tamasopo.com.mx%252Fgaleria2.html%3B1600%3B1200/ Puente de Dios in San Luis Potosi succeeded in distracting me for a time. When my godparents were finally able to secure me a job as a receptionist at one of the busiest hotels in Toluca, I also found that though I looked Mexican enough and spoke Spanish well enough, the slang and social references and popular customs eluded me to the point of ridicule amongst my coworkers. I did the best I could as I busied myself working 12 hours a day, 6 days a week, for the equivalent of $400 USD a month, but I knew Mexico would never be home and I longed for California with all my heart. 1.5 years later, just as I was giving up all hope, I received notice of the attorney generals approval of a waiver on the http://www.uscis.gov/sites/default/files/files/form/i601instr.pdf/ grounds of inadmissibility on my behalf due to my extenuating circumstances and otherwise good citizenship, not unlike many later http://dreamact.info/ DREAM Act petitioners. I was granted http://www.uscis.gov/tools/glossary/permanentresidentalien/ LPR status and was finally allowed to go home.

I returned to Sacramento to continue where Id left off and finally move forward with my life. I was welcomed back to my original serving job and branched out with a second job at http://blogs.wsj.com/bankruptcy/2012/10/25/lessonsfromthedeathofcircuitcity/ Circuit City for full time employment, but the http://www.csmonitor.com/Business/2013/0908/TimelineontheGreatRecession/ economic recession I found the country in didnt allow that to last for long. I found that I needed more help as an LPR than I ever did while living in the shadows and became indebted financially as well as to many friends and family members for their help. As I returned to school with a change to a more lucrative major, I accepted an invitation to reunite with another cousin in http://www.citydata.com/city/SuisunCityCalifornia.html/ Suisun City, CA rentfree for 9 months in order to afford payments on my first car while I commuted to work and school in Sacramento. I then took a break and followed a friend to work in http://www.citydata.com/city/LasVegasNevada.html/ Las Vegas, NV for 3 months. I returned and settled down again in Sacramento for another 3 years during which I helped one and eventually a second younger sister to follow me from So Cal and get settled much in the same way I had done almost 10 years earlier. I also took opportunities to volunteer at the http://www.ucdmc.ucdavis.edu/medicalcenter/index.html/ UC Davis Medical Centers https://www.ucdmc.ucdavis.edu/children/clinical_services/NICU// Neonatal I.C.U. and http://www.ucdmc.ucdavis.edu/emergency/ E.R. departments for career exploration, to pursue a natural love for nature in roadtripping/tidepooling/hiking around the Northern California and Oregon coast, and to volunteer with the http://www.deltaconservancy.ca.gov/aboutdeltaconservancy/ Delta Conservancy in Waterway Clean Up and conservation efforts. In the Spring of 2013, I finally graduated from http://www.crc.losrios.edu/ Cosumnes River College with Biology, PreNursing option A.S. & Social Science A.A. and prepared to apply to Bachelors of Science in Nursing programs, as well as for http://www.uscis.gov/uscitizenship/citizenshipthroughnaturalization/ Naturalization to US Citizen status next year.

Throughout my 28 years, Ive followed family, my own sense of adventure and need for independence, government orders, and friends in my wanderings within Mexico and California. Ive learned to disregard the http://idioms.thefreedictionary.com/keep+up+with+the+Joneses keeping up with the Joneses blues as I watched fellows from my https://www.facebook.com/HuenemeHighSchoolClassOf2003Reunions/ Class of 2003 accomplish goals, dreams, and lifes milestones far ahead of me. Ive worked hard, played little, and waited patiently for my time to come. When it did, unfulfilled high school dreams faded or transformed and gave way to new dreams, forged with a new sense of hope. Ive taken some time to take care of me and do things I had not the opportunity to do before. Amongst those things, Ive been lucky in finding my other half and Im happy to say that this time Ive followed my heart to http://www.citydata.com/city/DavisCalifornia.html/ Davis, CA. I havent been here long, but Im excited to explore what this town has to offer and what I can do for it as I share my life with Daniel Pfister my love and we pursue our dreams hand in hand.

(written/edited as practice to attain Wiki skills for a potential employer feel free to help :)

max_stars_repo_name: voflo/Search

max_issues_repo_name: voflo/Search

max_forks_repo_name: voflo/Search


text:
c
c     version 1.2  Mar 1992
c
c     ----------------------------------------------------------------
c     |                                                              |
c     |             ********  T X M A T C H   ********               |
c     |                                                              |
c     |         select picks with the same positions for two         |
c     |           different phases from two "tx.in" files            |
c     |                                                              |
c     |                   Written by C. A. Zelt                      |
c     |                                                              |
c     |                Geological Survey of Canada                   |
c     |                  Ottawa, Canada K1A 0Y3                      |
c     |                                                              |
c     ----------------------------------------------------------------
c
      npickm=0
      npicku=0
c
      write(6,55)
55    format(/'Enter file 1 and file 2 phases to match')
      read(5,*) ip1,ip2
c
      open(11, file='tx1.in', status='old')
      open(12, file='tx2.in', status='old')
      open(13, file='tx1.out')
      open(14, file='tx2.out')
c
100   read(11,*,end=999) x1,t1,u1,i1
c
      if(i1.eq.-1) go to 999
      if(i1.eq.0) then
        x1shot=x1
        t1shot=t1
        isw=0
        go to 100
      end if
      if(i1.ne.ip1) go to 100
c
101   read(12,*,end=998) x2,t2,u2,i2
      if(i2.eq.-1) go to 998
      if(i2.eq.0) then
        x2shot=x2
        t2shot=t2
        u2shot=u2
        i2shot=i2
        if(abs(x2shot-x1shot).lt..001.and.abs(t2shot-t1shot).lt..001)
     +  then
          iflag=1
        else
          iflag=0
        end if
        go to 101
      end if
      if(i2.ne.ip2) go to 101
      if(iflag.eq.1) then
        if(abs(x2-x1).lt..001) then
          if(isw.eq.0) then
            write(13,5) x2shot,t2shot,u2shot,i2shot
            write(14,5) x2shot,t2shot,u2shot,i2shot
            isw=1
            nshot=nshot+1
          end if
          write(13,5) x1,t1,u1,i1
          write(14,5) x2,t2,u2,i2
          npickm=npickm+1
          go to 100
        end if
      end if
      go to 101
c
998   npicku=npicku+1
      rewind(12)
      go to 100
999   write(13,5) x1,t1,u1,i1
      write(14,5) x1,t1,u1,i1
5     format(3f10.3,i10)
c
      write(*,*, fmt="('number of picks matched:   ',i6)") npickm
      write(*,*, fmt="('number of picks unmatched: ',i6)") npicku
      write(*,*, fmt="('number of shots:           ',i6/)") nshot
c
      stop
      end


text:
*----------------------------------------------------------------------*
      subroutine svd_drv(ndim,mat,singval)
*----------------------------------------------------------------------*
*     split matrix into non-coupling blocks and perform the svd
*     on these individual blocks
*     matrix should be symmetric
*
*     matthias, dec 2010
*----------------------------------------------------------------------*
      implicit none

      include 'stdunit.h'
      include 'routes.h'
      include 'def_filinf.h'
      include 'ifc_memman.h'

      integer, parameter ::
     &     ntest = 00
      real(8), parameter ::
     &     coupl_thresh = 1d-14 ! threshold for coupling elements
      integer, intent(in) ::
     &     ndim
      real(8), intent(out) ::
     &     singval(ndim)
      real(8), intent(inout) ::
     &     mat(ndim,ndim)
      real(8) ::
     &     dum1,dum2,umat(ndim,ndim)
      real(8), pointer ::
     &     mat_tmp(:,:),mat_tmpv(:), !mat_tmpv is only a handler for the mat_tmp memory as memman can only allocate one dimensional vectors
     &     sing_tmp(:), wrk(:)
      integer ::
     &     cur_els(ndim), nels, iels, iidim
      logical ::
     &     decomposed(ndim)

      integer ::
     &     nrot, idx, lwrk, info, idx2, luinp, iostatus,
     &     ifree
      integer, external ::
     &     idxlist

      umat = 0d0
      singval = 0d0
      decomposed = .false.

      do iidim = 1, ndim
        ifree = mem_setmark('svd_drv')
        if (decomposed(iidim)) cycle

        ! find all coupling elements
        nels = 1
        cur_els(1) = iidim
        iels = 1
        ! find all elements which couple with elements on list
        do while(iels.le.nels)
          do idx2 = 1, ndim
            ! skip if already on list
            if (idxlist(idx2,cur_els,nels,1).gt.0) cycle
            if (abs(mat(cur_els(iels),idx2)).gt.coupl_thresh) then
              nels = nels + 1
              cur_els(nels) = idx2
            end if
          end do
          decomposed(cur_els(iels)) = .true.
          iels = iels + 1
        end do
        ! zero element? do nothing
        if (nels.eq.1.and.mat(cur_els(1),cur_els(1)).lt.coupl_thresh)
     &       cycle
c dbg
c        print *,'cur_els: ',cur_els(1:nels)
c dbgend

! assemble block

        ifree = mem_alloc_real(mat_tmpv, nels*nels,'mat_tmp')
        ifree = mem_alloc_real(sing_tmp, nels,'sing_tmp')
        mat_tmp(1:nels,1:nels) => mat_tmpv
        do idx = 1, nels
          do idx2 = 1, nels
            mat_tmp(idx,idx2) = mat(cur_els(idx),cur_els(idx2))
          end do
        end do

        ! call singular value decomposer
        lwrk=max(1024,nels**2)
        ifree = mem_alloc_real(wrk, lwrk,'wrk')
        info = 0
        call dgesvd('O','N',nels,nels,
     &       mat_tmp,nels,sing_tmp,
     &       dum1,1,dum2,1,
     &       wrk,lwrk,info)
        if (info.ne.0) then
          write(lulog,*) 'WARNING in svd_drv: SVD in trouble'
        end if

        ! copy singular values and eigenvector
        do idx = 1, nels
          singval(cur_els(idx)) = sing_tmp(idx)
          do idx2 = 1, nels
            umat(cur_els(idx),cur_els(idx2)) = mat_tmp(idx,idx2)
          end do
        end do

        ifree = mem_flushmark()
      end do

      mat = umat

      return
      end

max_stars_repo_name: ak-ustutt/GeCCo-public

max_issues_repo_name: ak-ustutt/GeCCo-public

max_forks_repo_name: ak-ustutt/GeCCo-public


text:
C
C       $Id: vvinit.f,v 1.17 2008-07-27 00:17:35 haley Exp $
C
C                Copyright (C)  2000
C        University Corporation for Atmospheric Research
C                All Rights Reserved
C
C The use of this Software is governed by a License Agreement.
C
C-----------------------------------------------------------------------
C
      SUBROUTINE VVINIT (U,LU,V,LV,P,LP,M,N,WRK,LW)
C
C Argument dimensions.
C
      DIMENSION       U(LU,N)    ,V(LV,N)    ,P(LP,N)
      DIMENSION       WRK(LW)
C
C Input parameters
C
C U,V   - 2-d arrays holding the component values of a vector field
C LU,LV - The first dimensions of the U and V arrays, respectively
C P     - A 2-d array containing a scalar data field. The contents
C         of this array may be used to color the vectors.
C LP    - The first dimension of the P array
C M     - The first data dimension (must be less than or equal to
C         MIN(LU,LV) or, if the P array is used, MIN(LU,LV,LP))
C N     - The second data dimension
C WRK   - an internally used work array
C LW    - dimension of the work array (must be at least 2*M*N)
C
C Output parameters:
C
C None
C
C ---------------------------------------------------------------------
C
C NOTE:
C Since implicit typing is used for all real and integer variables
C a consistent length convention has been adopted to help clarify the
C significance of the variables encountered in the code for this
C utility. All local variable and subroutine parameter identifiers
C are limited to 1,2,or 3 characters. Four character names identify
C members of common blocks. Five and 6 character variable names
C denote PARAMETER constants or subroutine or function names.
C
C Declare the VV common blocks.
C
C IPLVLS - Maximum number of color threshold level values
C IPAGMX - Maximum number of area groups allowed in the area map
C
      PARAMETER (IPLVLS = 256, IPAGMX = 64)
C
C Integer and real common block variables
C
C
      COMMON /VVCOM/
     +                IUD1       ,IVD1       ,IPD1       ,IXDM       ,
     +                IYDN       ,VLOM       ,VHIM       ,ISET       ,
     +                VRMG       ,VRLN       ,VFRC       ,IXIN       ,
     +                IYIN       ,ISVF       ,UUSV       ,UVSV       ,
     +                UPSV       ,IMSK       ,ICPM       ,UVPS       ,
     +                UVPL       ,UVPR       ,UVPB       ,UVPT       ,
     +                UWDL       ,UWDR       ,UWDB       ,UWDT       ,
     +                UXC1       ,UXCM       ,UYC1       ,UYCN       ,
     +                NLVL       ,IPAI       ,ICTV       ,WDLV       ,
     +                UVMN       ,UVMX       ,PMIN       ,PMAX       ,
     +                RVMN       ,RVMX       ,RDMN       ,RDMX       ,
     +                ISPC       ,RVMD       ,IPLR       ,IVST       ,
     +                IVPO       ,ILBL       ,IDPF       ,IMSG       ,
     +                ICLR(IPLVLS)           ,TVLU(IPLVLS)
C
C Arrow size/shape parameters
C
        COMMON / VVARO /
     +                HDSZ       ,HINF       ,HANG       ,IAST       ,
     +                HSIN       ,HCOS       ,FAMN       ,FAMX       ,
     +                UVMG       ,FAIR       ,FAWR       ,FAWF       ,
     +                FAXR       ,FAXF       ,FAYR       ,FAYF       ,
     +                AROX(8)    ,AROY(8)    ,FXSZ       ,FYSZ       ,
     +                FXRF       ,FXMN       ,FYRF       ,FYMN       ,
     +                FWRF       ,FWMN       ,FIRF       ,FIMN       ,
     +                AXMN       ,AXMX       ,AYMN       ,AYMX       ,
     +                IACM       ,IAFO       ,WBAD       ,WBTF       ,
     +                WBCF       ,WBDF       ,WBSC
C
C
C Text related parameters
C
        COMMON /VVTXP /
     +                FCWM    ,ICSZ    ,
     +                FMNS    ,FMNX    ,FMNY    ,IMNP    ,IMNC  ,
     +                FMXS    ,FMXX    ,FMXY    ,IMXP    ,IMXC  ,
     +                FZFS    ,FZFX    ,FZFY    ,IZFP    ,IZFC  ,
     +                FILS    ,FILX    ,FILY    ,IILP    ,IILC  ,
     +                FLBS    ,ILBC

C
C Character variable declartions
C
      CHARACTER*160 CSTR
      PARAMETER (IPCHSZ=36)
      CHARACTER*(IPCHSZ)  CMNT,CMXT,CZFT,CLBT,CILT
C
C Text string parameters
C
      COMMON /VVCHAR/ CSTR,CMNT,CMXT,CZFT,CLBT,CILT
C
      SAVE /VVCOM/, /VVARO/, /VVTXP/, /VVCHAR/
C
C The mapping common block: made available to user mapping routines
C
      COMMON /VVMAP/
     +                IMAP       ,
     +                XVPL       ,XVPR       ,YVPB       ,YVPT       ,
     +                WXMN       ,WXMX       ,WYMN       ,WYMX       ,
     +                XLOV       ,XHIV       ,YLOV       ,YHIV       ,
     +                SXDC       ,SYDC       ,NXCT       ,NYCT       ,
     +                RLEN       ,LNLG       ,INVX       ,INVY       ,
     +                ITRT       ,IWCT       ,FW2W       ,FH2H       ,
     +                DVMN       ,DVMX       ,RBIG       ,IBIG
C
      SAVE /VVMAP/
C
C Math constants
C
      PARAMETER (PDTOR  = 0.017453292519943,
     +           PRTOD  = 57.2957795130823,
     +           P1XPI  = 3.14159265358979,
     +           P2XPI  = 6.28318530717959,
     +           P1D2PI = 1.57079632679489,
     +           P5D2PI = 7.85398163397448)
C
C --------------------------------------------------------------------
C
C Do a call forcing a BLOCKDATA to be loaded from a binary library.
C
      CALL VVDATA
C
C ---------------------------------------------------------------------
C
C Write the array sizes into the common block
C
      IUD1=LU
      IVD1=LV
      IPD1=LP
C
C Error if M > LU or M > LV
C
      IXDM=M
      IYDN=N
      IF (IXDM.GT.IUD1 .OR. IXDM.GT.IVD1) THEN
         CSTR(1:45)='VVINIT - U AND/OR V ARRAY DIMENSIONS EXCEEDED'
         CALL SETER (CSTR(1:45),1,1)
         RETURN
      END IF
      IF (RVMD .GT. 0.0 .AND. LW .LT. 2*IXDM*IYDN) THEN
         CSTR(1:50)='VVINIT - WORK ARRAY TOO SMALL FOR THINNING VECTORS'
         CALL SETER (CSTR(1:45),1,1)
         RETURN
      END IF
C
C Initialize and transfer some arguments to local variables.
C
      IBIG = I1MACH(9)
      RBIG = R1MACH(2)
C
C Decide what the range of values in X and Y should be.
C
      IF (UXC1.EQ.UXCM) THEN
        XLOV=1.
        XHIV=REAL(MAX(IXDM, 2))
      ELSE
        XLOV=UXC1
        XHIV=UXCM
      END IF
C
      IF (UYC1.EQ.UYCN) THEN
        YLOV=1.
        YHIV=REAL(MAX(IYDN,2))
      ELSE
        YLOV=UYC1
        YHIV=UYCN
      END IF
C
      IXIN = MAX(IXIN,1)
      IYIN = MAX(IYIN,1)
C
      NXCT = IXDM/IXIN
      NYCT = IYDN/IYIN
C
C If the user has done a SET call, retrieve the arguments; if he hasn't
C done a SET call, do it for him.
C
      IF (ISET .EQ .0) THEN
C
        CALL GETSET (XVPL,XVPR,YVPB,YVPT,WXMN,WXMX,WYMN,WYMX,LNLG)
C
      ELSE
C
        LNLG=1
C
        IF (UWDL.EQ.UWDR) THEN
          WXMN=XLOV
          WXMX=XHIV
        ELSE
          WXMN=UWDL
          WXMX=UWDR
        END IF
C
        IF (UWDB.EQ.UWDT) THEN
          WYMN=YLOV
          WYMX=YHIV
        ELSE
          WYMN=UWDB
          WYMX=UWDT
        END IF
C
C Determine the viewport based on the setting of the viewport
C shape and viewport extent parameters
C
        IF (UVPS.LT.0.) THEN
          AR=ABS(UVPS)
        ELSE IF (UVPS.EQ.0.) THEN
          AR=(UVPR-UVPL)/(UVPT-UVPB)
        ELSE IF (UVPS.LE.1.) THEN
          AR=ABS((WXMX-WXMN)/(WYMX-WYMN))
          IF (MIN(AR,1./AR).LT.UVPS) AR=(UVPR-UVPL)/(UVPT-UVPB)
        ELSE
          AR=ABS((WXMX-WXMN)/(WYMX-WYMN))
          IF (MAX(AR,1./AR).GT.UVPS) AR=1.
        END IF
C
        IF (AR.LT.(UVPR-UVPL)/(UVPT-UVPB)) THEN
          XVPL=.5*(UVPL+UVPR)-.5*(UVPT-UVPB)*AR
          XVPR=.5*(UVPL+UVPR)+.5*(UVPT-UVPB)*AR
          YVPB=UVPB
          YVPT=UVPT
        ELSE
          XVPL=UVPL
          XVPR=UVPR
          YVPB=.5*(UVPB+UVPT)-.5*(UVPR-UVPL)/AR
          YVPT=.5*(UVPB+UVPT)+.5*(UVPR-UVPL)/AR
        END IF
C
        CALL SET (XVPL,XVPR,YVPB,YVPT,WXMN,WXMX,WYMN,WYMX,LNLG)
C
      END IF
C
C Calculate fraction of VP width to fractional size factor.
C Calculate fraction of VP height to fractional size factor.
C These are for convenience.
C
      FW2W = XVPR - XVPL
      FH2H = YVPT - YVPB
C
C Swap window rectangle if it is inverted, but keep track
C This makes it easier to exclude out-of-bounds points in the
C projection mapping routines
C
      INVX=0
      INVY=0
      IF (WXMN .GT. WXMX) THEN
         T=WXMN
         WXMN=WXMX
         WXMX=T
         INVX=1
      END IF
      IF (WYMN .GT. WYMX) THEN
         T=WYMN
         WYMN=WYMX
         WYMX=T
         INVY=1
      END IF
C
C Calculate a default value for the maximum vector length in
C fractional coordinates. Set minimum length to 0.0.
C
      IXIN = MAX(IXIN, 1)
      IYIN = MAX(IYIN, 1)
      I=MAX(IXDM/IXIN, 1)
      J=MAX(IYDN/IYIN, 1)
C
      TX=(XVPR-XVPL)/REAL(I)
      TY=(YVPT-YVPB)/REAL(J)
      DVMX=SQRT((TX*TX+TY*TY)/2.0)
      DVMN = 0.0
C
C Determine the equivalent maximum length in the user coordinate
C system. Note that if the X axis is non-uniform or if the Y axis
C has a different scaling from the X-Axis, this value may not be
C meaningful.
C
      RLEN=DVMX*(WXMX-WXMN)/(XVPR-XVPL)
C
C Find the maximum and minimum magnitude vectors.
C If there are no special values use simpler code for efficiency.
C
      UVMN=RBIG
      UVMX=0.0
      IF (ISVF .LE. 0) THEN
C
         DO 111 J=1,IYDN,IYIN
            DO 110 I=1,IXDM,IXIN
               TU=U(I,J)
               TV=V(I,J)
               IF (IPLR.NE.0) THEN
                  UVMN = MIN(UVMN, TU)
                  UVMX = MAX(UVMX, TU)
               ELSE
                  UVMN = MIN(UVMN, TU*TU+TV*TV)
                  UVMX = MAX(UVMX, TU*TU+TV*TV)
               END IF
 110        CONTINUE
 111     CONTINUE
C
      ELSE
C
         DO 116 J=1,IYDN,IYIN
            DO 115 I=1,IXDM,IXIN
               TU=U(I,J)
               TV=V(I,J)
C
               IF (ISVF .EQ. 1) THEN
                  IF (TU .EQ. UUSV) GO TO 115
               ELSE IF (ISVF .EQ. 2) THEN
                  IF (TV .EQ. UVSV) GO TO 115
               ELSE IF (ISVF .EQ. 3) THEN
                  IF (TU .EQ. UUSV .OR. TV .EQ. UVSV) THEN
                     GO TO 115
                  END IF
               ELSE IF (ISVF .EQ. 4) THEN
                  IF (TU .EQ. UUSV .AND. TV .EQ. UVSV) THEN
                     GO TO 115
                  END IF
               END IF
C
               IF (IPLR.NE.0) THEN
                  UVMN = MIN(UVMN, TU)
                  UVMX = MAX(UVMX, TU)
               ELSE
                  UVMN = MIN(UVMN, TU*TU+TV*TV)
                  UVMX = MAX(UVMX, TU*TU+TV*TV)
               END IF
C
 115        CONTINUE
 116     CONTINUE
C
      END IF
C
      IF (IPLR .LE. 0) THEN
         UVMX = SQRT(UVMX)
         UVMN = SQRT(UVMN)
      END IF
C
C Decide how to use the scalar array.
C
      IF (ABS(ICTV).EQ.1) THEN
         PMIN=UVMN
         PMAX=UVMX
         IP=1
      ELSE IF (ICTV.NE.0) THEN
C
C Check for error condition
C
         IF (IPD1 .LT. IXDM) THEN
            CSTR(1:31)='VVINIT - SCALAR ARRAY TOO SMALL'
            CALL SETER (CSTR(1:31),2,1)
            RETURN
         END IF
C
C Find the max and min scalar values when necessary
C
         IF (ICTV.GE.2 .OR. PMIN.EQ.PMAX) THEN
            PMIN=RBIG
            PMAX=-RBIG
            IP=0
            DO 121 J=1,IYDN,IYIN
               DO 120 I=1,IXDM,IXIN
                  IF (ISPC .GE. 0 .AND. P(I,J) .EQ. UPSV) GO TO 120
                  IP=1
                  PMIN = MIN(PMIN, P(I,J))
                  PMAX = MAX(PMAX, P(I,J))
 120           CONTINUE
 121        CONTINUE
         END IF
C
      END IF
C
C Determine the threshold values. If there were no non-special
C values in the P array, then set all threshold values to the
C maximum possible.
C
      IF (NLVL.GT.0 .AND. ICTV.GT.0) THEN
C
         IF (IP .NE. 0) THEN
            TV=(PMAX-PMIN)/NLVL
            DO 125 I=1,NLVL,1
               TVLU(I)=PMIN+REAL(I)*TV
 125        CONTINUE
         ELSE
            DO 130 I=1,NLVL,1
               TVLU(I)=RBIG
 130        CONTINUE
         END IF
C
      END IF
C
C Copy min/max vector magnitude and NDC size to read-only params.
C
      RVMN=UVMN
      RVMX=UVMX
      RDMN=DVMN
      RDMX=DVMX
C
C Done.
C
      RETURN
C
      END

max_stars_repo_name: tenomoto/ncl

max_issues_repo_name: tenomoto/ncl

max_forks_repo_name: tenomoto/ncl


! { dg-do run }
!
! PR fortran/50981
! Test the handling of optional, polymorphic and non-polymorphic arguments
! to elemental procedures.
!
! Original testcase by Tobias Burnus <burnus@net-b.de>

implicit none
type t
  integer :: a
end type t

type t2
  integer, allocatable :: a
  integer, allocatable :: a2(:)
  integer, pointer :: p => null()
  integer, pointer :: p2(:) => null()
end type t2

type(t), allocatable :: ta, taa(:)
type(t), pointer :: tp, tpa(:)
class(t), allocatable :: ca, caa(:)
class(t), pointer :: cp, cpa(:)

type(t2) :: x

integer :: s, v(2)

tp => null()
tpa => null()
cp => null()
cpa => null()

! =============== sub1 ==================
! SCALAR COMPONENTS: Non alloc/assoc

s = 3
v = [9, 33]

call sub1 (s, x%a, .false.)
call sub1 (v, x%a, .false.)
!print *, s, v
if (s /= 3) STOP 1
if (any (v /= [9, 33])) STOP 2

call sub1 (s, x%p, .false.)
call sub1 (v, x%p, .false.)
!print *, s, v
if (s /= 3) STOP 3
if (any (v /= [9, 33])) STOP 4

! SCALAR COMPONENTS: alloc/assoc

allocate (x%a, x%p)
x%a = 4
x%p = 5
call sub1 (s, x%a, .true.)
call sub1 (v, x%a, .true.)
!print *, s, v
if (s /= 4*2) STOP 5
if (any (v /= [4*2, 4*2])) STOP 6

call sub1 (s, x%p, .true.)
call sub1 (v, x%p, .true.)
!print *, s, v
if (s /= 5*2) STOP 7
if (any (v /= [5*2, 5*2])) STOP 8

! ARRAY COMPONENTS: Non alloc/assoc

v = [9, 33]

call sub1 (v, x%a2, .false.)
!print *, v
if (any (v /= [9, 33])) STOP 9

call sub1 (v, x%p2, .false.)
!print *, v
if (any (v /= [9, 33])) STOP 10

! ARRAY COMPONENTS: alloc/assoc

allocate (x%a2(2), x%p2(2))
x%a2(:) = [84, 82]
x%p2    = [35, 58]

call sub1 (v, x%a2, .true.)
!print *, v
if (any (v /= [84*2, 82*2])) STOP 11

call sub1 (v, x%p2, .true.)
!print *, v
if (any (v /= [35*2, 58*2])) STOP 12

! =============== sub_t ==================
! SCALAR DT: Non alloc/assoc

s = 3
v = [9, 33]

call sub_t (s, ta, .false.)
call sub_t (v, ta, .false.)
!print *, s, v
if (s /= 3) STOP 13
if (any (v /= [9, 33])) STOP 14

call sub_t (s, tp, .false.)
call sub_t (v, tp, .false.)
!print *, s, v
if (s /= 3) STOP 15
if (any (v /= [9, 33])) STOP 16

call sub_t (s, ca, .false.)
call sub_t (v, ca, .false.)
!print *, s, v
if (s /= 3) STOP 17
if (any (v /= [9, 33])) STOP 18

call sub_t (s, cp, .false.)
call sub_t (v, cp, .false.)
!print *, s, v
if (s /= 3) STOP 19
if (any (v /= [9, 33])) STOP 20

! SCALAR COMPONENTS: alloc/assoc

allocate (ta, tp, ca, cp)
ta%a = 4
tp%a = 5
ca%a = 6
cp%a = 7

call sub_t (s, ta, .true.)
call sub_t (v, ta, .true.)
!print *, s, v
if (s /= 4*2) STOP 21
if (any (v /= [4*2, 4*2])) STOP 22

call sub_t (s, tp, .true.)
call sub_t (v, tp, .true.)
!print *, s, v
if (s /= 5*2) STOP 23
if (any (v /= [5*2, 5*2])) STOP 24

call sub_t (s, ca, .true.)
call sub_t (v, ca, .true.)
!print *, s, v
if (s /= 6*2) STOP 25
if (any (v /= [6*2, 6*2])) STOP 26

call sub_t (s, cp, .true.)
call sub_t (v, cp, .true.)
!print *, s, v
if (s /= 7*2) STOP 27
if (any (v /= [7*2, 7*2])) STOP 28

! ARRAY COMPONENTS: Non alloc/assoc

v = [9, 33]

call sub_t (v, taa, .false.)
!print *, v
if (any (v /= [9, 33])) STOP 29

call sub_t (v, tpa, .false.)
!print *, v
if (any (v /= [9, 33])) STOP 30

call sub_t (v, caa, .false.)
!print *, v
if (any (v /= [9, 33])) STOP 31

call sub_t (v, cpa, .false.)
!print *, v
if (any (v /= [9, 33])) STOP 32

deallocate(ta, tp, ca, cp)

! ARRAY COMPONENTS: alloc/assoc

allocate (taa(2), tpa(2))
taa(1:2)%a = [44, 444]
tpa(1:2)%a = [55, 555]
allocate (caa(2), source=[t(66), t(666)])
allocate (cpa(2), source=[t(77), t(777)])

select type (caa)
type is (t)
  if (any (caa(:)%a /= [66, 666])) STOP 33
end select

select type (cpa)
type is (t)
  if (any (cpa(:)%a /= [77, 777])) STOP 34
end select

call sub_t (v, taa, .true.)
!print *, v
if (any (v /= [44*2, 444*2])) STOP 35

call sub_t (v, tpa, .true.)
!print *, v
if (any (v /= [55*2, 555*2])) STOP 36

call sub_t (v, caa, .true.)
!print *, v
if (any (v /= [66*2, 666*2])) STOP 37

call sub_t (v, cpa, .true.)
!print *, v
if (any (v /= [77*2, 777*2])) STOP 38

deallocate (taa, tpa, caa, cpa)

contains

  elemental subroutine sub1 (x, y, alloc)
    integer, intent(inout) :: x
    integer, intent(in), optional :: y
    logical, intent(in) :: alloc
    if (alloc .neqv. present (y)) &
      x = -99
    if (present(y)) &
      x = y*2
  end subroutine sub1

  elemental subroutine sub_t(x, y, alloc)
    integer, intent(inout) :: x
    type(t), intent(in), optional :: y
    logical, intent(in) :: alloc
    if (alloc .neqv. present (y)) &
      x = -99
    if (present(y)) &
      x = y%a*2
  end subroutine sub_t

end


text:
C +-======-+
C  Copyright (c) 2003-2007 United States Government as represented by
C  the Admistrator of the National Aeronautics and Space Administration.
C  All Rights Reserved.
C
C  THIS OPEN  SOURCE  AGREEMENT  ("AGREEMENT") DEFINES  THE  RIGHTS  OF USE,
C  REPRODUCTION,  DISTRIBUTION,  MODIFICATION AND REDISTRIBUTION OF CERTAIN
C  COMPUTER SOFTWARE ORIGINALLY RELEASED BY THE UNITED STATES GOVERNMENT AS
C  REPRESENTED BY THE GOVERNMENT AGENCY LISTED BELOW ("GOVERNMENT AGENCY").
C  THE UNITED STATES GOVERNMENT, AS REPRESENTED BY GOVERNMENT AGENCY, IS AN
C  INTENDED  THIRD-PARTY  BENEFICIARY  OF  ALL  SUBSEQUENT DISTRIBUTIONS OR
C  REDISTRIBUTIONS  OF THE  SUBJECT  SOFTWARE.  ANYONE WHO USES, REPRODUCES,
C  DISTRIBUTES, MODIFIES  OR REDISTRIBUTES THE SUBJECT SOFTWARE, AS DEFINED
C  HEREIN, OR ANY PART THEREOF,  IS,  BY THAT ACTION, ACCEPTING IN FULL THE
C  RESPONSIBILITIES AND OBLIGATIONS CONTAINED IN THIS AGREEMENT.
C
C  Government Agency: National Aeronautics and Space Administration
C  Government Agency Original Software Designation: GSC-15354-1
C  Government Agency Original Software Title:  GEOS-5 GCM Modeling Software
C  User Registration Requested.  Please Visit http://opensource.gsfc.nasa.gov
C  Government Agency Point of Contact for Original Software:
C  			Dale Hithon, SRA Assistant, (301) 286-2691
C
C +-======-+
       module ESMF_CFIOBaseMod
       include "netcdf.inc"
       end module ESMF_CFIOBaseMod

max_stars_repo_name: joeylamcy/gchp

max_stars_repo_licenses: NCSA Apache-2.0 MIT

max_issues_repo_name: joeylamcy/gchp

max_issues_repo_licenses: NCSA Apache-2.0 MIT

max_forks_repo_name: joeylamcy/gchp

max_forks_repo_licenses: NCSA Apache-2.0 MIT


! ****************************************************************************
! ... Module for variable types
! ... COSMO Project
! ... Quim Ballabrera, March 2017
! ****************************************************************************

MODULE types

IMPLICIT NONE

PRIVATE
PUBLIC sp,dp

INTEGER, PARAMETER                 :: sp = KIND(0.0)    ! Single precision
INTEGER, PARAMETER                 :: dp = KIND(0.0D0)  ! Double precision

END MODULE types


text:
FUNCTION AIRTOVAC(lam)

  !convert wavelengths from air to vac
  !see Morton (1991 Ap.J. Suppl. 77, 119)
  !this code was adapted from the IDL routine airtovac.pro

  USE sps_vars; USE sps_utils, ONLY : locate
  IMPLICIT NONE

  INTEGER :: vv,nn
  REAL(SP), DIMENSION(:), INTENT(in) :: lam
  REAL(SP), DIMENSION(SIZE(lam)) :: airtovac,fact,sigma2

  !------------------------------------------------------!

  nn = SIZE(lam)

  !Convert to wavenumber squared
  sigma2(1:nn) = (1E4/lam(1:nn) )**2.

  !Compute conversion factor
  fact(1:nn) = 1.D0 + 6.4328D-5 + 2.94981D-2/(146.d0 - sigma2(1:nn)) + &
       2.5540D-4/( 41.d0 - sigma2(1:nn))

  !no conversion for wavelengths <2000A
  IF (lam(1).LT.2000.D0) THEN
     vv = MIN(MAX(locate(lam(1:nn),2000.D0),1),nn)
     fact(1:vv) = 1.0
  ENDIF

  !Convert wavelength
  airtovac(1:nn) = lam(1:nn)*fact(1:nn)

END FUNCTION AIRTOVAC

!---------------------------------------------------------------!
!---------------------------------------------------------------!
!---------------------------------------------------------------!

FUNCTION VACTOAIR(lam)

  !convert wavelengths from air to vac
  !see Morton (1991 Ap.J. Suppl. 77, 119)
  !this code was adapted from the IDL routine vactoair.pro

  USE sps_vars; USE sps_utils, ONLY : locate
  IMPLICIT NONE

  INTEGER :: vv,nn
  REAL(SP), DIMENSION(:), INTENT(in) :: lam
  REAL(SP), DIMENSION(SIZE(lam)) :: vactoair,fact

  !------------------------------------------------------!

  nn = SIZE(lam)

  fact(1:nn) = 1.D0 + 2.735182D-4 + 131.4182/lam(1:nn)**2 + &
       2.76249D8/lam(1:nn)**4

  !no conversion for wavelengths <2000A
   IF (lam(1).LT.2000.D0) THEN
     vv = MIN(MAX(locate(lam(1:nn),2000.D0),1),nn)
     fact(1:vv) = 1.0
  ENDIF

  !Convert wavelengths
  vactoair(1:nn) = lam(1:nn)/fact(1:nn)

END FUNCTION VACTOAIR

max_stars_repo_name: christopherlovell/fsps

max_issues_repo_name: christopherlovell/fsps

max_forks_repo_name: christopherlovell/fsps


text:
ccc_____________________________________________________________________
ccc
      subroutine align_SEED(npts,nch,ix,ierr)

ccc   starting from a set of data blocks in array next_block,
ccc   find first possible common start time, initialize
ccc   array ix and starting record number irec  ...

      include 'seed.inc'
      integer ix(0:nch,*),nch,npts,ich,i,ierr,read_ms,k,
     &   irec_1(2,nchmx),irec_l(2,nchmx),start_rec,i0
      logical check_again
      double precision tdiff
      record /INT_TIME/exptime

ccc   find start, end sample numbers of (already input) next data block
ccc   (stored in next_blk)
      write(0,*) 'doy0 ==',doy0
      start_rec = 0
      do ich = 1,nch
         call rec_num(hdr(ich).begtime,doy0,samp_rate,irec_1(1,ich))
         start_rec = max(start_rec,irec_1(1,ich))
         call rec_num(hdr(ich).endtime,doy0,samp_rate,irec_l(1,ich))
      end do

      write(0,*) 'In align_seed : irec_1',(irec_1(1,k),k=1,5)
      write(0,*) 'In align_seed : irec_l',(irec_l(1,k),k=1,5)

ccc   make sure data blocks all contain a common starting record
ccc   read additional blocks for "early" channels until this is
ccc   achieved
      check_again = .true.
      do while(check_again)
         check_again = .false.
         do ich = 1,nch
            do while(start_rec .gt. irec_l(1,ich))
               check_again = .true.
               ngot(ich) = read_ms(hdr(ich),next_blk(1,ich),nblkmx,
     &         inunits(ich))
               if(ngot(ich).le.0) then
ccc               end of file hit while searching for alignement
                  ierr = -1
                  return
               end if
               call rec_num(hdr(ich).begtime,doy0,samp_rate,
     &               irec_1(1,ich))
               start_rec = max(start_rec,irec_1(1,ich))
               call rec_num(hdr(ich).endtime,doy0,samp_rate,
     &                   irec_l(1,ich))
            end do
         end do
      end do

ccc   JUne 1997  add code to allow for offsets in sampling ...
      do ich = 1,nch
        dt(ich,1) = float(irec_1(2,ich))/TICKS_PER_SEC
        if(samp_rate.lt.0) then
           dt(ich,1) = dt(ich,1)*abs(samp_rate)
        else
           dt(ich,1) = dt(ich,1)/samp_rate
        endif
        dt_chng_rec(ich,1) = maxint
        dt_chng_n(ich) = 1
        dt_chng_i(ich) = 1
      enddo
      write(0,*) 'dt', (dt(ich,1),ich=1,nch)

c      irec = start_rec - 1
       write(0,*) 'start_rec = ',start_rec
      irec = start_rec
ccc   initialize ix, nix
      do ich = 1,nch
         i0 = start_rec - irec_1(1,ich)
         do i = i0+1,ngot(ich)
            ix(ich,i) = next_blk(i-i0,ich)
         end do
         nix(ich) = ngot(ich) - i0
      end do
      write(0,*) 'In align_seed : nix',(nix(k),k=1,5)

ccc   load next data block for each channel
      do ich = 1,nch
ccc      first save last header block
         call time_interval(ngot(ich),hdr(ich).sample_rate,second,
     &         ticks)
         call add_time(hdr(ich).begtime,second,ticks,exptime)
         ngot(ich) = read_ms(hdr(ich),next_blk(1,ich),nblkmx,
     &                      inunits(ich))
c         write(0,*) 'ngot(ich),ich',ngot(ich),ich
         if(ngot(ich) .le. 0 ) then
            iend(ich) = -1
            ngot(ich) = 0
         else
ccc         compute gap between blocks (if any) for channel ich
            iend(ich) = 0
            ngap(ich) = nint(tdiff(hdr(ich).begtime,exptime)/
     &              (TICKS_PER_SEC*ddr))
            if(ngap(ich) .gt. 0 ) then
               call rec_num(hdr(ich).begtime,doy0,samp_rate,
     &                     irec_1(1,ich))
               dt_chng_rec(ich,dt_chng_n(ich)) = irec_1(1,ich)
               dt_chng_n(ich) = dt_chng_n(ich) + 1
               dt(ich,dt_chng_n(ich)) =
     &                  float(irec_1(2,ich))/TICKS_PER_SEC
               dt_chng_rec(ich,dt_chng_n(ich)) = maxint
            endif
         endif
      end do
c      write(0,*) 'in align seeed : ngot',(ngot(k),k=1,nch)
c      write(0,*) 'in align seeed : ngap',(ngap(k),k=1,nch)
      npts = 0
c      write(0,*) 'IEND',(iend(k),k=1,nch)
      return
      end

max_stars_repo_name: simpeg-research/iris-mt-scratch

max_issues_repo_name: simpeg-research/iris-mt-scratch

max_forks_repo_name: simpeg-research/iris-mt-scratch