Scripts for WW3 data pre/post processing #6
Comments
|
Please continue the discussion here. Thanks, |
|
@aliabdolali, I can run your example case and reproduce the results, however, when I use the GFS wave spec which I am using in my configuration, it produces the following error for the ascii read ([IN] = swden_ww3_read_ascii('../obc/gfswave.22108.spec');): could you please help me fix the error for those files? |
|
please pull wave-tools and try again. |
|
Hi @aliabdolali , Thank you, that error disappeared, however, the results seem there is a mismatch in the reading, please take a look at the figure, |
|
You need to be careful about spectral resolution, WW3 freq increment is logarithmic (the larger the frequeny, the wider the spectral band), and in GFS, it starts with 0.035 with 1.07. |
|
when you read netcdf files, it automatically calculates hs, tp, ..., so you can compare them to make sure the mismatch is acceptable or not. |
|
will continue updating here as needed as @yunfangsun builds scripts/notebooks for WW3 model pre- and post-processing. |
|
Update from Yunfang: updating scripts to prep namelist generation |
|
Hi I am using schism analytical wwm test case (Test_WWM_Analytical) to test ufs-coastal (schism2ww3). Then I fake both mesh into geographic (lon/lat), above CMEPS errors are gone, but I got error from WW3 like the following: @aliabdolali told me it is due to the precision issue in WW3. So, I don’t have solution for this case and turn to test with real case (suggested by @josephzhang8 ). To @yunfangsun, do you know how to introduce open boundary into WW3 for real case? |
|
@danishyo I could not see the error from screenshot. Is it possible to attach actual log? CMEPs does not need to know anything about the grid but if there is an issue with the grid or the definition in the SCHSIM side, that could trigger issue. You could setup the same case by removing mediator to see what happens. |
|
@uturuncoglu Sorry for the figure, I didn't notice it's blurry after uploading. |
|
@danishyo Okay. It seems you are trying to couple SCHSIM with WW3. Right? I think your are using cartesian grid in WW3 side. Right? It seems that you are using |
|
@uturuncoglu This case is using Cartesian setting for both schism & WW3. It is set (ics=1) Cartesian in schism with its grid setup (hgrid.gr3), as well as WW3 is set GRID%COORD = 'CART' in ww3_grid.nml when generating mod_def.ww3. For spherical case setup, you can see /scratch1/06923/hyu05/FV3/ufs-coastal/ww3_couple_test_cpp |
|
@danishyo Okay. I am not sure but there might be some restriction in CMEPS side and may not support Cartesian case and it would be hard to bring that capability. But, You could try to remove mediator from your run sequence and config files and just using NUOPC connectors and we could test the case like it. I could try to modify your cartesian case next week and let you know. Is this the cartesian one, /scratch1/06923/hyu05/FV3/ufs-coastal/ww3_couple_test? |
|
@uturuncoglu Thanks! Yes, that's the right one. |
|
@uturuncoglu Could you give a hint how to setup ufs.configure by only using NUOPC connectors? |
|
@danishyo Try following ufs.configure: BTW, i am not sure why you set coupling interval to 1 (in the run sequence). Do you want to couple each component every second? Also, since we remove the mediator, this will use 40 processor in total. So, you need to adjust your job submission script. |
|
@danishyo BTW, if both component are using same grid you could use following run sequence. |
|
@uturuncoglu Both setting show errors in PET log but different. The 2nd one shows You can find them at the following path: |
|
@danishyo Okay. If you are getting error in the first one then it means that CMEPS is not the issue in here. Maybe there is some kind of configuration issue. You could also try to remove |
|
@uturuncoglu I checked OCN and WW3 grid and they are the same from their original ascii files. By the way, could you show me how to generate data nc files like "wind_atm_fin_ch_time_vec_ESMFmesh.nc & wind_atm_fin_ch_time_vec_STR_fixed.nc" from scratch? |
@aliabdolali shared a private repo (https://github.com/erdc/wave-tools) containing pre/post scripts for spectral wave models.
@saeed-moghimi-noaa
Attached is one example of directional spec interpolation from ascii to netcdf.
clear all
addpath /Users/rdchlaa9/Desktop/Denise_Grid/wave-tools/bin/matlab/
% read asii spec
[IN] = swden_ww3_read_ascii('bnd_ndbc.spec.spc');
%% user defined interpolation parameters
% input data %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
nfreq=32; % number of frequencies
nDir=36; % number of Directions
inc=1.1; % frequency increment
f0=0.050;
filenameIN='bnd_ndbc.spec_ww3.nc';% name of netcdf file (boundary)
filenameOUT='bnd_ndbc.spec_interp_ww3.nc';% name of netcdf file (boundary)
testcase='JXFL57'; % test vase
pointID='JXFL57'; % id (length should be 16)
coordinate='spherical'; % coordinate : Spherical, cartesian
visualize='true';
%----------------------------------------------------------%%----------------------------------------------------------%
deltatheta=360/nDir; % DeltaDir
%----------------------------------------------------------%
% convert ascii spec to nc spec
[filename] = write_directional_spectra_nc(filenameIN,testcase,...
IN.pointID,IN.lat,IN.lon,IN.dpt,IN.wnd,IN.wnddir,IN.cur,IN.curdir,...
IN.time,IN.frequency,IN.direction,IN.efth,coordinate);
% interpolate input spec to target spectral resolution
[IN,OUT]= spec_interp_nc(filenameIN,filenameOUT,nfreq,f0,inc,nDir,...
coordinate,visualize);
This script is using swden_ww3_read_ascii:
write_directional_spectra_nc is as follow:
the interpolation
function[IN,OUT]=spec_interp_nc(filenameIN,filenameOUT,nfreq,f0,inc,nDir,coordinate,visualize)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This program reads all netcdf file contents %
% [email protected] March 06, 2024 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%% INPUT %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% ncf: name of input netcdf file
% nfreq=35; % number of frequencies
% nDir=36; % number of Directions
% inc=1.1; % frequency increment
% f0=0.038; % first freq
% visualize='true';
% filenameIN='JXFL57_swan.nc';% name of netcdf file (boundary)
% coordinate='spherical'; % coordinate : Spherical, cartesian
%%%%%%%%%%%%%%%%%%% OUTPUT %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% filenameOUT='JXFL57_ww3.nc';% name of netcdf file (boundary)
%%%%%%%%%%%%%%%%%%%%%%%% Dependency %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% convert_time
% read_nc
% dictionary_wave
%%%%%%%%%%%%%%%%%%% example %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%[IN,OUT]= spec_interp_nc('in.nc','out.nc',35,0.038,1.1,36,...
% 'spherical','true')
%-------------------------------------------------------------------------
%-------------------------------------------------------------------------%
testcase=['interpolated from ',filenameIN]; % test vase
deltatheta=360/nDir; % DeltaDir
% Dir0=5; % first direction (deg)
% theta0=0; % first Dir
%-------------------------------------------------------------------------%
% read input
IN=read_nc(filenameIN);
% sort direction and rearrange efth
[IN.direction,id]=sort(IN.direction);
IN.efth=IN.efth(id,:,:);
%-------------------------------------------------------------------------%
% prepare freq and direction
f=f0;
for i=1:nfreq-1
f(end+1)=f(end)inc;
end
dir0(1,:)=90:-deltatheta:0;
dir0(end+1:nDir)=360-deltatheta+dir0(end):-deltatheta:90+deltatheta-dir0(end);
dir(1,:)=pidir0/180;
[dirtmpO,ftmpO]=meshgrid(dir0,f);
%-------------------------------------------------------------------------%
% create output structure
OUT=IN;
OUT=rmfield(OUT,'efth');
OUT.pointID=[[OUT.station_name{:}],repmat(' ', [1, 16-strlength([OUT.station_name{:}])])];
%-------------------------------------------------------------------------%
% make a tmp IN and add two columns on both sides of directions
INtmp=IN;
INtmp=rmfield(INtmp,'direction');
INtmp=rmfield(INtmp,'efth');
imax=find(IN.direction==max(IN.direction(:)));
imin=find(IN.direction==min(IN.direction(:)));
INtmp.direction=[IN.direction(imax)-360; IN.direction;IN.direction(imin)+360];
INtmp.efth=[IN.efth(imax,:,:); IN.efth; IN.efth(imin,:,:)];
%-------------------------------------------------------------------------%
% interpolation
[dirtmp,ftmp]=meshgrid(INtmp.direction,INtmp.frequency);
for i=1:length(IN.time)
DENStmp(:,:)=INtmp.efth(:,:,i);
DENStmp2=DENStmp';
tmp=reshape(interp2(dirtmp,ftmp,DENStmp',dirtmpO(:),ftmpO(:)),size(dirtmpO));
OUT.efth(:,:,1,i)=tmp';
end
%-------------------------------------------------------------------------%
display (['Generating ', filenameOUT,' ...'])
%dump into netcdf
[filename] = write_directional_spectra_nc(filenameOUT,testcase,...
OUT.pointID,OUT.lat,OUT.lon,OUT.dpt,OUT.wnd,OUT.wnddir,...
OUT.cur,OUT.curdir,OUT.time,f,dir0,OUT.efth,coordinate);
%%
%%
tf=strcmp(visualize,'true');
if tf==1
display (['Plot Generation ', [filename,'.gif'],' ...'])
width=500; % Width of figure for movie [pixels]
height=800; % Height of figure of movie [pixels]
left=200; % Left margin between figure and screen edge [pixels]
bottom=200; % Bottom margin between figure and screen edge [pixels]
h=figure('Color','w');
set(gcf,'Position', [left bottom width height])
for k=1:length(OUT.time)
k
clf
sss1=subplot(2,1,1);
clear SPEC1
clear SPECC
SPEC1(:,:)=IN.efth(:,:,k);
SPECC=SPEC1;
[~,cc] = polarPcolor(IN.frequency',[IN.direction; IN.direction(1)]',...
[SPECC; SPECC(1,:)],'Nspokes',36,'ncolor',10,'labelR','f (Hz)','Rscale','log');
ylabel(cc,['INPUT - ',datestr(IN.time(k))],'FontSize',14);
sss2=subplot(2,1,2);
clear SPECC
SPECC=OUT.efth(:,:,1,k);
[~,cc] = polarPcolor(f,[dir0 dir0(1)],...
[SPECC; SPECC(1,:)],'Nspokes',36,'ncolor',10,'labelR','f (Hz)','Rscale','log');
ylabel(cc,['OUTPUT - ',datestr(OUT.time(k))],'FontSize',14);
set(sss1, 'Position', [.08 .55 .88 .4]);
set(sss2, 'Position', [.08 .08 .88 .4]);
frame = getframe(h);
im = frame2im(frame);
[imind,cm] = rgb2ind(im,256);
end
end
%----------------------------------------------------------%
display (['Finished'])
The text was updated successfully, but these errors were encountered: