Skip to content

Commit

Permalink
address comments
Browse files Browse the repository at this point in the history
  • Loading branch information
apchoiCMD committed Nov 14, 2024
1 parent 9ea280f commit b3d78ed
Show file tree
Hide file tree
Showing 3 changed files with 255 additions and 1 deletion.
1 change: 0 additions & 1 deletion scripts/old/exgdas_global_marine_analysis_vrfy.py
Original file line number Diff line number Diff line change
@@ -1,4 +1,3 @@
#!/usr/bin/env python3
################################################################################
# UNIX Script Documentation Block
# . .
Expand Down
210 changes: 210 additions & 0 deletions utils/soca/fig_gallery/exgdas_global_marine_analysis_vrfy_manual.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,210 @@
import os
import numpy as np
import gen_eva_obs_yaml
import marine_eva_post
import diag_statistics
from multiprocessing import Process
from soca_vrfy import statePlotter, plotConfig
import subprocess

comout = os.getenv('COM_OCEAN_ANALYSIS')
com_ice_history = os.getenv('COM_ICE_HISTORY_PREV')
com_ocean_history = os.getenv('COM_OCEAN_HISTORY_PREV')
cyc = os.getenv('cyc')
RUN = os.getenv('RUN')

bcyc = str((int(cyc) - 3) % 24).zfill(2)
gcyc = str((int(cyc) - 6) % 24).zfill(2)
grid_file = os.path.join(comout, f'{RUN}.t'+bcyc+'z.ocngrid.nc')
layer_file = os.path.join(comout, f'{RUN}.t'+cyc+'z.ocninc.nc')

# for eva
diagdir = os.path.join(comout, 'diags')
HOMEgfs = os.getenv('HOMEgfs')

# Get flags from environment variables (set in the bash driver)
run_ensemble_analysis = os.getenv('RUN_ENSENBLE_ANALYSIS', 'OFF').upper() == 'ON'
run_bkgerr_analysis = os.getenv('RUN_BACKGROUND_ERROR_ANALYSIS', 'OFF').upper() == 'ON'
run_bkg_analysis = os.getenv('RUN_BACKGROUND_ANALYSIS', 'OFF').upper() == 'ON'
run_increment_analysis = os.getenv('RUN_INCREMENT_ANLYSIS', 'OFF').upper() == 'ON'

# Initialize an empty list for the main config
configs = [plotConfig(grid_file=grid_file,
data_file=os.path.join(comout, f'{RUN}.t'+cyc+'z.ocnana.nc'),
variables_horiz={'ave_ssh': [-1.8, 1.3],
'Temp': [-1.8, 34.0],
'Salt': [32, 40]},
colormap='nipy_spectral',
comout=os.path.join(comout, 'vrfy', 'ana')), # ocean surface analysis
plotConfig(grid_file=grid_file,
data_file=os.path.join(comout, f'{RUN}.t'+cyc+'z.iceana.nc'),
variables_horiz={'aice_h': [0.0, 1.0],
'hi_h': [0.0, 4.0],
'hs_h': [0.0, 0.5]},
colormap='jet',
projs=['North', 'South', 'Global'],
comout=os.path.join(comout, 'vrfy', 'ana'))] # sea ice analysis

# Define each config and add to main_config if its flag is True
if run_ensemble_analysis:
config_ens = [plotConfig(grid_file=grid_file,
data_file=os.path.join(comout, f'{RUN}.t{cyc}z.ocn.recentering_error.nc'),
variables_horiz={'ave_ssh': [-1, 1]},
colormap='seismic',
comout=os.path.join(comout, 'vrfy', 'recentering_error')), # recentering error
plotConfig(grid_file=grid_file,
data_file=os.path.join(comout, f'{RUN}.t{cyc}z.ocn.ssh_steric_stddev.nc'),
variables_horiz={'ave_ssh': [0, 0.8]},
colormap='gist_ncar',
comout=os.path.join(comout, 'vrfy', 'bkgerr', 'ssh_steric_stddev')), # ssh steric stddev
plotConfig(grid_file=grid_file,
data_file=os.path.join(comout, f'{RUN}.t{cyc}z.ocn.ssh_unbal_stddev.nc'),
variables_horiz={'ave_ssh': [0, 0.8]},
colormap='gist_ncar',
comout=os.path.join(comout, 'vrfy', 'bkgerr', 'ssh_unbal_stddev')), # ssh unbal stddev
plotConfig(grid_file=grid_file,
data_file=os.path.join(comout, f'{RUN}.t{cyc}z.ocn.ssh_total_stddev.nc'),
variables_horiz={'ave_ssh': [0, 0.8]},
colormap='gist_ncar',
comout=os.path.join(comout, 'vrfy', 'bkgerr', 'ssh_total_stddev')), # ssh total stddev
plotConfig(grid_file=grid_file,
data_file=os.path.join(comout, f'{RUN}.t{cyc}z.ocn.steric_explained_variance.nc'),
variables_horiz={'ave_ssh': [0, 1]},
colormap='seismic',
comout=os.path.join(comout, 'vrfy', 'bkgerr', 'steric_explained_variance'))] # steric explained variance
configs.extend(config_ens)

if run_bkgerr_analysis:
config_bkgerr = [plotConfig(grid_file=grid_file,
layer_file=layer_file,
data_file=os.path.join(comout, os.path.pardir, os.path.pardir,
'bmatrix', 'ocean', f'{RUN}.t'+cyc+'z.ocean.bkgerr_stddev.nc'),
lats=np.arange(-60, 60, 10),
lons=np.arange(-280, 80, 30),
variables_zonal={'Temp': [0, 2],
'Salt': [0, 0.2],
'u': [0, 0.2],
'v': [0, 0.2]},
variables_meridional={'Temp': [0, 2],
'Salt': [0, 0.2],
'u': [0, 0.2],
'v': [0, 0.2]},
variables_horiz={'Temp': [0, 2],
'Salt': [0, 0.2],
'u': [0, 0.2],
'v': [0, 0.2],
'ave_ssh': [0, 0.1]},
colormap='jet',
comout=os.path.join(comout, 'vrfy', 'bkgerr'))] # ocn bkgerr stddev
configs.extend(config_bkgerr)

if run_bkg_analysis:
config_bkg = [plotConfig(grid_file=grid_file,
data_file=os.path.join(com_ice_history, f'{RUN}.ice.t{gcyc}z.inst.f006.nc'),
variables_horiz={'aice_h': [0.0, 1.0],
'hi_h': [0.0, 4.0],
'hs_h': [0.0, 0.5]},
colormap='jet',
projs=['North', 'South', 'Global'],
comout=os.path.join(comout, 'vrfy', 'bkg')), # sea ice background
plotConfig(grid_file=grid_file,
layer_file=layer_file,
data_file=os.path.join(com_ocean_history, f'{RUN}.ocean.t{gcyc}z.inst.f006.nc'),
lats=np.arange(-60, 60, 10),
lons=np.arange(-280, 80, 30),
variables_zonal={'Temp': [-1.8, 34.0],
'Salt': [32, 40]},
variables_meridional={'Temp': [-1.8, 34.0],
'Salt': [32, 40]},
variables_horiz={'ave_ssh': [-1.8, 1.3],
'Temp': [-1.8, 34.0],
'Salt': [32, 40]},
colormap='nipy_spectral',
comout=os.path.join(comout, 'vrfy', 'bkg'))]
configs.extend(config_bkg)

if run_increment_analysis:
config_incr = [plotConfig(grid_file=grid_file,
layer_file=layer_file,
data_file=os.path.join(comout, f'{RUN}.t'+cyc+'z.ocninc.nc'),
lats=np.arange(-60, 60, 10),
lons=np.arange(-280, 80, 30),
variables_zonal={'Temp': [-0.5, 0.5],
'Salt': [-0.1, 0.1]},
variables_horiz={'Temp': [-0.5, 0.5],
'Salt': [-0.1, 0.1],
'ave_ssh': [-0.1, 0.1]},
variables_meridional={'Temp': [-0.5, 0.5],
'Salt': [-0.1, 0.1]},
colormap='seismic',
comout=os.path.join(comout, 'vrfy', 'incr')), # ocean increment
plotConfig(grid_file=grid_file,
data_file=os.path.join(comout, f'{RUN}.t'+cyc+'z.ice.incr.nc'),
lats=np.arange(-60, 60, 10),
variables_horiz={'aice_h': [-0.2, 0.2],
'hi_h': [-0.5, 0.5],
'hs_h': [-0.1, 0.1]},
colormap='seismic',
projs=['North', 'South'],
comout=os.path.join(comout, 'vrfy', 'incr'))] # sea ice increment
configs.extend(config_incr)


# plot marine analysis vrfy

def plot_marine_vrfy(config):
ocnvrfyPlotter = statePlotter(config)
ocnvrfyPlotter.plot()


# Number of processes
num_processes = len(configs)

# Create a list to store the processes
processes = []

# Iterate over configs
for config in configs[:num_processes]:
process = Process(target=plot_marine_vrfy, args=(config,))
process.start()
processes.append(process)

# Wait for all processes to finish
for process in processes:
process.join()

#######################################
# eva plots
#######################################

evadir = os.path.join(HOMEgfs, 'sorc', f'{RUN}.cd', 'ush', 'eva')
marinetemplate = os.path.join(evadir, 'marine_gdas_plots.yaml')
varyaml = os.path.join(comout, 'yaml', 'var_original.yaml')

# it would be better to refrence the dirs explicitly with the comout path
# but eva doesn't allow for specifying output directories
os.chdir(os.path.join(comout, 'vrfy'))
if not os.path.exists('preevayamls'):
os.makedirs('preevayamls')
if not os.path.exists('evayamls'):
os.makedirs('evayamls')

gen_eva_obs_yaml.gen_eva_obs_yaml(varyaml, marinetemplate, 'preevayamls')

files = os.listdir('preevayamls')
for file in files:
infile = os.path.join('preevayamls', file)
marine_eva_post.marine_eva_post(infile, 'evayamls', diagdir)

files = os.listdir('evayamls')
for file in files:
infile = os.path.join('evayamls', file)
print('running eva on', infile)
subprocess.run(['eva', infile], check=True)

#######################################
# calculate diag statistics
#######################################

# As of 11/12/2024 not working
# diag_statistics.get_diag_stats()
45 changes: 45 additions & 0 deletions utils/soca/fig_gallery/run_marine_analysis_vrfy_manual.job
Original file line number Diff line number Diff line change
@@ -0,0 +1,45 @@
#!/bin/bash
#SBATCH --job-name=marine_vrfy # Assign a name to the job (customize as needed)
#SBATCH --account=da-cpu
#SBATCH --qos=debug
#SBATCH -A da-cpu
#SBATCH --output=run_marine_vrfy_analysis.out
#SBATCH --nodes=1 # Request 1 node
#SBATCH --ntasks=40 # Request 40 total tasks (processors across nodes)
#SBATCH --partition=hera # Specify the partition (cluster) named "hera"
#SBATCH --cpus-per-task=1 # Set 1 CPU per task (equivalent to ppn=40 and tpp=1)
#SBATCH --mem=24GB # Request 24GB of memory
#SBATCH --time=00:30:00 # Set the walltime limit to 30 minutes

# Define HOMEgfs
export HOMEgfs="/scratch1/NCEPDEV/da/Mindo.Choi/workflow_11122024/global-workflow/"

# Load EVA module
module use ${HOMEgfs}sorc/gdas.cd/modulefiles
module load EVA/hera

# Set PYTHONPATH using HOMEgfs
export PYTHONPATH="${HOMEgfs}sorc/gdas.cd/ush/:\
${HOMEgfs}sorc/gdas.cd/ush/eva/:\
${HOMEgfs}sorc/gdas.cd/ush/soca/:\
$PYTHONPATH"

# Set flags to control plotConfig in the Python script
export RUN_ENSENBLE_ANALYSIS=OFF # Check if ensemble run is ON
export RUN_BACKGROUND_ERROR_ANALYSIS=ON
export RUN_BACKGROUND_ANALYSIS=ON
export RUN_INCREMENT_ANLYSIS=ON

# Define and export the environment variables
export cyc="00"
export RUN="gdas"
export PSLOT="gdas_test"
export PDY="20210827"

# Define and export environment variables with paths
export COM_OCEAN_ANALYSIS="/scratch1/NCEPDEV/da/Mindo.Choi/sandbox/marine_vrfy/gdas.20210827/00/analysis/ocean"
export COM_ICE_HISTORY_PREV="/scratch1/NCEPDEV/da/Mindo.Choi/sandbox/marine_vrfy/gdas.20210826/18/model/ice/history"
export COM_OCEAN_HISTORY_PREV="/scratch1/NCEPDEV/da/Mindo.Choi/sandbox/marine_vrfy/gdas.20210826/18/model/ocean/history"

# Excute Marine Verify Analysis
python3 ${HOMEgfs}sorc/gdas.cd/scripts/exgdas_global_marine_analysis_vrfy_manual.py

0 comments on commit b3d78ed

Please sign in to comment.