Skip to content
Snippets Groups Projects
Commit 5a50a51f authored by Evgeny Mortikov's avatar Evgeny Mortikov
Browse files

obl output update

parent 0fff4552
No related branches found
No related tags found
No related merge requests found
......@@ -77,7 +77,7 @@ program obl_main
!1 - pacanowski-philander, 2 - pacanowski-philander+,
!3 - k-epsilon explicit, 4 - k-epsilon semiimplicit, 5 - inmom
integer, parameter :: output_mode = 1 ! 1 -- netcdf, 2 -- ascii, 3 -- tecplot
integer, parameter :: output_mode = 3 ! 1 -- netcdf, 2 -- ascii, 3 -- tecplot
integer, parameter :: obl_setup = 1 ! 1 - Kato-Phillips, 2 - Papa station
......@@ -422,6 +422,9 @@ program obl_main
call push_value_to_tseries(output_mld, mld)
call push_value_to_tseries(output_mld_ref, lab_mld)
call push_value_to_tseries(output_tau_x, flux_u_surf_res)
call push_value_to_tseries(output_tau_y, flux_v_surf_res)
call push_value_to_tseries(output_time, time_current / 3600.0)
endif
enddo
......@@ -431,46 +434,31 @@ program obl_main
output_Salin%data(:,1:output_Salin%num) = output_Salin%data(:,1:output_Salin%num) + Salin_ref
output_TinC%data(:,1:output_TinC%num) = output_TinC%data(:,1:output_TinC%num) + Theta_ref - 273.15
output_tau_x%data(1:output_tau_x%num) = output_tau_x%data(1:output_tau_x%num) * Rho_ref
output_tau_y%data(1:output_tau_y%num) = output_tau_y%data(1:output_tau_y%num) * Rho_ref
endif
if (output_mode.eq.1) then
write(*, *) ' >> writing netcdf output ...'
! Writing 2D arrays (variables with depth and time)
call write_netcdf
! Writing 1D arrays (scalar variables over time)
!call write_1d_real_nc(Tau_x_surf, 'output.nc', meta_tau_u)
!call write_1d_real_nc(Tau_y_surf, 'output.nc', meta_tau_v)
endif
if (output_mode.eq.2) then
write(*, *) ' >> writing ascii output ...'
! Writing 2D arrays (variables with depth and time) to ASCII files
call write_ascii
! Writing 1D arrays (scalar variables over time) to ASCII files
!call write_1d_real_ascii(Tau_x_surf, 'output_tau_u.txt', meta_tau_u)
!call write_1d_real_ascii(Tau_y_surf, 'output_tau_v.txt', meta_tau_v)
endif
if (output_mode.eq.3) then
write(*, *) ' >> writing tecplot output ...'
! time slice output
call write_tecplot(grid%z, output_time%data)
! Writing 1D arrays (scalar variables over time) to Tecplot files
!call write_1d_real_ascii(Tau_x_surf, 'output_tau_u.txt', meta_tau_u)
!call write_1d_real_ascii(Tau_y_surf, 'output_tau_v.txt', meta_tau_v)
endif
!close(199)
!close(20)
! deallocation
!> deallocation
call deallocate_state_vec
deallocate(dz)
......
......@@ -37,6 +37,8 @@ module obl_output
type(oblTimeSeries), public :: output_mld
type(oblTimeSeries), public :: output_mld_ref
type(oblTImeSeries), public :: output_tau_x, output_tau_y
contains
......@@ -98,6 +100,9 @@ module obl_output
call write_1d_real_nc(output_mld%data(1:output_mld%num), 'output.nc', meta_mld)
call write_1d_real_nc(output_mld_ref%data(1:output_mld_ref%num), 'output.nc', meta_lab_mld)
call write_1d_real_nc(output_tau_x%data(1:output_tau_x%num), 'output.nc', meta_tau_u)
call write_1d_real_nc(output_tau_y%data(1:output_tau_y%num), 'output.nc', meta_tau_v)
end subroutine write_netcdf
! --------------------------------------------------------------------------------
......@@ -147,6 +152,11 @@ module obl_output
call write_1d_real_ascii(output_mld_ref%data(1:output_mld_ref%num), &
trim(path) // trim(meta_lab_mld%name) // '.txt', meta_lab_mld)
call write_1d_real_ascii(output_tau_x%data(1:output_tau_x%num), &
trim(path) // trim(meta_tau_u%name) // '.txt', meta_tau_u)
call write_1d_real_ascii(output_tau_y%data(1:output_tau_y%num), &
trim(path) // trim(meta_tau_v%name) // '.txt', meta_tau_v)
end subroutine write_ascii
! --------------------------------------------------------------------------------
......@@ -200,6 +210,11 @@ module obl_output
call write_1d_real_plt(output_mld_ref%data(1:output_mld_ref%num), tcoord(1:output_mld_ref%num), &
trim(path) // trim(meta_lab_mld%name) // '.plt', meta_lab_mld)
call write_1d_real_plt(output_tau_x%data(1:output_tau_x%num), tcoord(1:output_tau_x%num), &
trim(path) // trim(meta_tau_u%name) // '.plt', meta_tau_u)
call write_1d_real_plt(output_tau_y%data(1:output_tau_y%num), tcoord(1:output_tau_y%num), &
trim(path) // trim(meta_tau_v%name) // '.plt', meta_tau_v)
end subroutine write_tecplot
! --------------------------------------------------------------------------------
......@@ -232,6 +247,9 @@ module obl_output
call deallocate_tseries(output_mld)
call deallocate_tseries(output_mld_ref)
call deallocate_tseries(output_tau_x)
call deallocate_tseries(output_tau_y)
end subroutine output_cleanup
end module
\ No newline at end of file
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment