This subroutine sets the open face lengths at selected points to restrict passages to their observed widths from a list read from a file.
808 type(dyn_horgrid_type),
intent(inout) :: G
809 type(param_file_type),
intent(in) :: param_file
810 type(unit_scale_type),
optional,
intent(in) :: US
813 character(len=120),
pointer,
dimension(:) :: lines => null()
814 character(len=120) :: line
815 character(len=200) :: filename, chan_file, inputdir, mesg
816 character(len=40) :: mdl =
"reset_face_lengths_list"
817 real,
pointer,
dimension(:,:) :: &
818 u_lat => null(), u_lon => null(), v_lat => null(), v_lon => null()
819 real,
pointer,
dimension(:) :: &
820 u_width => null(), v_width => null()
827 logical :: found_u, found_v
828 logical :: unit_in_use
829 integer :: ios, iounit, isu, isv
830 integer :: last, num_lines, nl_read, ln, npt, u_pt, v_pt
831 integer :: i, j, isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
832 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
833 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
835 call calltree_enter(trim(mdl)//
"(), MOM_shared_initialization.F90")
837 call get_param(param_file, mdl,
"CHANNEL_LIST_FILE", chan_file, &
838 "The file from which the list of narrowed channels is read.", &
839 default=
"MOM_channel_list")
840 call get_param(param_file, mdl,
"INPUTDIR", inputdir, default=
".")
841 inputdir = slasher(inputdir)
842 filename = trim(inputdir)//trim(chan_file)
843 call log_param(param_file, mdl,
"INPUTDIR/CHANNEL_LIST_FILE", filename)
844 call get_param(param_file, mdl,
"CHANNEL_LIST_360_LON_CHECK", check_360, &
845 "If true, the channel configuration list works for any "//&
846 "longitudes in the range of -360 to 360.", default=.true.)
848 if (is_root_pe())
then
850 if (.not.file_exists(filename))
call mom_error(fatal, &
851 " reset_face_lengths_list: Unable to open "//trim(filename))
855 INQUIRE(iounit,opened=unit_in_use) ;
if (.not.unit_in_use)
exit
857 if (iounit >= 512)
call mom_error(fatal, &
858 "reset_face_lengths_list: No unused file unit could be found.")
861 open(iounit, file=trim(filename), access=
'SEQUENTIAL', &
862 form=
'FORMATTED', action=
'READ', position=
'REWIND', iostat=ios)
863 if (ios /= 0)
call mom_error(fatal, &
864 "reset_face_lengths_list: Error opening "//trim(filename))
867 call read_face_length_list(iounit, filename, num_lines, lines)
870 len_lon = 360.0 ;
if (g%len_lon > 0.0) len_lon = g%len_lon
871 len_lat = 180.0 ;
if (g%len_lat > 0.0) len_lat = g%len_lat
873 call broadcast(num_lines, root_pe())
875 if (num_lines > 0)
then
876 allocate (lines(num_lines))
877 if (num_lines > 0)
then
878 allocate(u_lat(2,num_lines)) ; u_lat(:,:) = -1e34
879 allocate(u_lon(2,num_lines)) ; u_lon(:,:) = -1e34
880 allocate(u_width(num_lines)) ; u_width(:) = -1e34
882 allocate(v_lat(2,num_lines)) ; v_lat(:,:) = -1e34
883 allocate(v_lon(2,num_lines)) ; v_lon(:,:) = -1e34
884 allocate(v_width(num_lines)) ; v_width(:) = -1e34
888 if (is_root_pe())
then
889 call read_face_length_list(iounit, filename, nl_read, lines)
890 if (nl_read /= num_lines) &
891 call mom_error(fatal,
'reset_face_lengths_list : Found different '// &
892 'number of valid lines on second reading of '//trim(filename))
893 close(iounit) ; iounit = -1
897 call broadcast(lines, 120, root_pe())
903 found_u = .false.; found_v = .false.
904 isu = index(uppercase(line),
"U_WIDTH" );
if (isu > 0) found_u = .true.
905 isv = index(uppercase(line),
"V_WIDTH" );
if (isv > 0) found_v = .true.
910 read(line(isu+8:),*) u_lon(1:2,u_pt), u_lat(1:2,u_pt), u_width(u_pt)
911 if (is_root_pe())
then
913 if ((abs(u_lon(1,u_pt)) > len_lon) .or. (abs(u_lon(2,u_pt)) > len_lon)) &
914 call mom_error(warning,
"reset_face_lengths_list : Out-of-bounds "//&
915 "u-longitude found when reading line "//trim(line)//
" from file "//&
917 if ((abs(u_lat(1,u_pt)) > len_lat) .or. (abs(u_lat(2,u_pt)) > len_lat)) &
918 call mom_error(warning,
"reset_face_lengths_list : Out-of-bounds "//&
919 "u-latitude found when reading line "//trim(line)//
" from file "//&
922 if (u_lat(1,u_pt) > u_lat(2,u_pt)) &
923 call mom_error(warning,
"reset_face_lengths_list : Out-of-order "//&
924 "u-face latitudes found when reading line "//trim(line)//
" from file "//&
926 if (u_lon(1,u_pt) > u_lon(2,u_pt)) &
927 call mom_error(warning,
"reset_face_lengths_list : Out-of-order "//&
928 "u-face longitudes found when reading line "//trim(line)//
" from file "//&
930 if (u_width(u_pt) < 0.0) &
931 call mom_error(warning,
"reset_face_lengths_list : Negative "//&
932 "u-width found when reading line "//trim(line)//
" from file "//&
935 elseif (found_v)
then
937 read(line(isv+8:),*) v_lon(1:2,v_pt), v_lat(1:2,v_pt), v_width(v_pt)
938 if (is_root_pe())
then
940 if ((abs(v_lon(1,v_pt)) > len_lon) .or. (abs(v_lon(2,v_pt)) > len_lon)) &
941 call mom_error(warning,
"reset_face_lengths_list : Out-of-bounds "//&
942 "v-longitude found when reading line "//trim(line)//
" from file "//&
944 if ((abs(v_lat(1,v_pt)) > len_lat) .or. (abs(v_lat(2,v_pt)) > len_lat)) &
945 call mom_error(warning,
"reset_face_lengths_list : Out-of-bounds "//&
946 "v-latitude found when reading line "//trim(line)//
" from file "//&
949 if (v_lat(1,v_pt) > v_lat(2,v_pt)) &
950 call mom_error(warning,
"reset_face_lengths_list : Out-of-order "//&
951 "v-face latitudes found when reading line "//trim(line)//
" from file "//&
953 if (v_lon(1,v_pt) > v_lon(2,v_pt)) &
954 call mom_error(warning,
"reset_face_lengths_list : Out-of-order "//&
955 "v-face longitudes found when reading line "//trim(line)//
" from file "//&
957 if (v_width(v_pt) < 0.0) &
958 call mom_error(warning,
"reset_face_lengths_list : Negative "//&
959 "v-width found when reading line "//trim(line)//
" from file "//&
968 do j=jsd,jed ;
do i=isdb,iedb
969 lat = g%geoLatCu(i,j) ; lon = g%geoLonCu(i,j)
970 if (check_360)
then ; lon_p = lon+len_lon ; lon_m = lon-len_lon
971 else ; lon_p = lon ; lon_m = lon ;
endif
974 if (((lat >= u_lat(1,npt)) .and. (lat <= u_lat(2,npt))) .and. &
975 (((lon >= u_lon(1,npt)) .and. (lon <= u_lon(2,npt))) .or. &
976 ((lon_p >= u_lon(1,npt)) .and. (lon_p <= u_lon(2,npt))) .or. &
977 ((lon_m >= u_lon(1,npt)) .and. (lon_m <= u_lon(2,npt)))) )
then
979 g%dy_Cu(i,j) = g%mask2dCu(i,j) * min(g%dyCu(i,j), max(u_width(npt), 0.0))
980 if (j>=g%jsc .and. j<=g%jec .and. i>=g%isc .and. i<=g%iec)
then
981 if ( g%mask2dCu(i,j) == 0.0 )
then
982 write(*,
'(A,2F8.2,A,4F8.2,A)')
"read_face_lengths_list : G%mask2dCu=0 at ",lat,lon,
" (",&
983 u_lat(1,npt), u_lat(2,npt), u_lon(1,npt), u_lon(2,npt),
") so grid metric is unmodified."
985 write(*,
'(A,2F8.2,A,4F8.2,A5,F9.2,A1)') &
986 "read_face_lengths_list : Modifying dy_Cu gridpoint at ",lat,lon,
" (",&
987 u_lat(1,npt), u_lat(2,npt), u_lon(1,npt), u_lon(2,npt),
") to ",g%dy_Cu(i,j),
"m"
993 g%areaCu(i,j) = g%dxCu(i,j)*g%dy_Cu(i,j)
995 if (g%areaCu(i,j) > 0.0) g%IareaCu(i,j) = g%mask2dCu(i,j) / g%areaCu(i,j)
998 do j=jsdb,jedb ;
do i=isd,ied
999 lat = g%geoLatCv(i,j) ; lon = g%geoLonCv(i,j)
1000 if (check_360)
then ; lon_p = lon+len_lon ; lon_m = lon-len_lon
1001 else ; lon_p = lon ; lon_m = lon ;
endif
1004 if (((lat >= v_lat(1,npt)) .and. (lat <= v_lat(2,npt))) .and. &
1005 (((lon >= v_lon(1,npt)) .and. (lon <= v_lon(2,npt))) .or. &
1006 ((lon_p >= v_lon(1,npt)) .and. (lon_p <= v_lon(2,npt))) .or. &
1007 ((lon_m >= v_lon(1,npt)) .and. (lon_m <= v_lon(2,npt)))) )
then
1008 g%dx_Cv(i,j) = g%mask2dCv(i,j) * min(g%dxCv(i,j), max(v_width(npt), 0.0))
1009 if (i>=g%isc .and. i<=g%iec .and. j>=g%jsc .and. j<=g%jec)
then
1010 if ( g%mask2dCv(i,j) == 0.0 )
then
1011 write(*,
'(A,2F8.2,A,4F8.2,A)')
"read_face_lengths_list : G%mask2dCv=0 at ",lat,lon,
" (",&
1012 v_lat(1,npt), v_lat(2,npt), v_lon(1,npt), v_lon(2,npt),
") so grid metric is unmodified."
1014 write(*,
'(A,2F8.2,A,4F8.2,A5,F9.2,A1)') &
1015 "read_face_lengths_list : Modifying dx_Cv gridpoint at ",lat,lon,
" (",&
1016 v_lat(1,npt), v_lat(2,npt), v_lon(1,npt), v_lon(2,npt),
") to ",g%dx_Cv(i,j),
"m"
1022 g%areaCv(i,j) = g%dyCv(i,j)*g%dx_Cv(i,j)
1023 g%IareaCv(i,j) = 0.0
1024 if (g%areaCv(i,j) > 0.0) g%IareaCv(i,j) = g%mask2dCv(i,j) / g%areaCv(i,j)
1027 if (num_lines > 0)
then
1028 deallocate(u_lat) ;
deallocate(u_lon) ;
deallocate(u_width)
1029 deallocate(v_lat) ;
deallocate(v_lon) ;
deallocate(v_width)
1032 call calltree_leave(trim(mdl)//
'()')