import os, pickle, io
import numpy as np
import xarray as xr
import matplotlib.pyplot as plt
import datetime
import re
native_outputs_path='/proj/bolinc/users/x_abaro/MIMICA_radiation/COMBLE_NEW_SIMULATIONS/'
converted_outputs_path='/proj/bolinc/users/x_abaro/MIMICA_radiation/COMBLE_CONVERTED/'
converted_outputs_path_time_mod='/proj/bolinc/users/x_abaro/MIMICA_radiation/COMBLE_CONVERTED_time_mod/'
#COMBLE_part2='/proj/mimica-chalmers-comble/users/x_abaro/COMBLE/'
Output_statistic = {
'ps': { 'std_name': 'surface_pressure', 'var_id': 'ps', 'units': 'Pa', 'dim': 'time'},
'ts': { 'std_name': 'surface_temperature', 'var_id': 'ts', 'units': 'K' , 'dim': 'time'},
'ustar': { 'std_name': 'surface_friction_velocity', 'var_id': 'ustar', 'units': 'm s-1', 'dim': 'time'},
'z0' : { 'std_name': 'surface_roughness_length_for_momentum_in_air', 'var_id': 'z0', 'units': 'm', 'dim': 'time'},
'z0h': { 'std_name': 'surface_roughness_length_for_heat_in_air', 'var_id': 'z0h', 'units': 'm', 'dim': 'time'},
'z0q': { 'std_name': 'surface_roughness_length_for_humidity_in_air', 'var_id': 'z0q', 'units': 'm', 'dim': 'time'},
'hfss': { 'std_name': 'surface_upward_sensible_heat_flux', 'var_id': 'hfss', 'units': 'W m-2', 'dim': 'time'},
'hfls': { 'std_name': 'surface_upward_latent_heat_flux', 'var_id': 'hfls', 'units': 'W m-2', 'dim': 'time'},
'ol' : { 'std_name': 'obukhov_length', 'var_id': 'ol', 'units': 'm', 'dim': 'time'},
'lwpc': { 'std_name': 'atmosphere_mass_content_of_liquid_cloud_water', 'var_id': 'lwpc', 'units': 'kg m-2', 'dim': 'time'},
'lwpr': { 'std_name': 'atmosphere_mass_content_of_rain_water', 'var_id': 'lwpr', 'units': 'kg m-2', 'dim': 'time'},
'iwp' : { 'std_name': 'atmosphere_mass_content_of_ice_water', 'var_id': 'iwp', 'units': 'kg m-2', 'dim': 'time'},
'clt' : { 'std_name': 'cloud_area_fraction', 'var_id': 'clt', 'units': '1', 'dim': 'time'},
'od' : { 'std_name': 'optical_depth', 'var_id': 'od', 'units': '1', 'dim': 'time'},
'odlc' : { 'std_name': 'optical_depth_of_liquid_cloud', 'var_id': 'odlc', 'units': '1', 'dim': 'time'},
'pr' : { 'std_name': 'precipitation_flux_at_surface', 'var_id': 'pr', 'units': 'kg m-2 s-1', 'dim': 'time'},
'pri' : { 'std_name': 'precipitation_flux_at_surface_in_ice_phase', 'var_id': 'pri', 'units': 'kg m-2 s-1', 'dim': 'time'},
'rlut' : { 'std_name': 'toa_outgoing_longwave_flux', 'var_id': 'rlut', 'units': 'W m-2', 'dim': 'time'},
'rlds' : { 'std_name': 'surface_downwelling_longwave_flux', 'var_id': 'rlds', 'units': 'W m-2', 'dim': 'time'},
'rlus' : { 'std_name': 'surface_upwelling_longwave_flux', 'var_id': 'rlus', 'units': 'W m-2', 'dim': 'time'},
'ssaf' : { 'std_name': 'surface_sea_spray_number_flux', 'var_id': 'ssaf', 'units': 'm-2 s-1', 'dim': 'time'},
'pa' : { 'std_name': 'air_pressure', 'var_id': 'pa', 'units': 'Pa', 'dim': 'time, height'},
'ua' : { 'std_name': 'eastward_wind', 'var_id': 'ua', 'units': 'm s-1', 'dim': 'time, height'},
'va' : { 'std_name': 'northward_wind', 'var_id': 'va', 'units': 'm s-1', 'dim': 'time, height'},
'rhoa' : { 'std_name': 'air_dry_density', 'var_id': 'rhoa', 'units': 'kg m-3', 'dim': 'time, height'},
'ta' : { 'std_name': 'air_temperature', 'var_id': 'ta', 'units': 'K', 'dim': 'time, height'},
'qv' : { 'std_name': 'water_vapor_mixing_ratio', 'var_id': 'qv', 'units': 'kg kg-1', 'dim': 'time, height'},
'hur' : { 'std_name': 'relative_humidity', 'var_id': 'hur', 'units': '1', 'dim': 'time, height'},
'huri' : { 'std_name': 'relative_humidity_over_ice', 'var_id': 'huri', 'units': '1', 'dim': 'time, height'},
'theta' : { 'std_name': 'air_potential_temperature', 'var_id': 'theta', 'units': 'K', 'dim': 'time, height'},
'tke_res': { 'std_name': 'specific_turbulent_kinetic_energy_resolved', 'var_id': 'tke_res', 'units': 'm2 s-2', 'dim': 'time, height'},
'tke_sgs': { 'std_name': 'specific_turbulent_kinetic_energy_sgs', 'var_id': 'tke_sgs', 'units': 'm2 s-2', 'dim': 'time, height'},
'qlc': { 'std_name': 'mass_mixing_ratio_of_cloud_liquid_water_in_air','var_id': 'qlc', 'units': 'kg kg-1', 'dim': 'time, height'},
'qlr': { 'std_name': 'mass_mixing_ratio_of_rain_water_in_air', 'var_id': 'qlr', 'units': 'kg kg-1', 'dim': 'time, height'},
'qic': { 'std_name': 'mass_mixing_ratio_of_cloud_ice_in_air', 'var_id': 'qic', 'units': 'kg kg-1', 'dim': 'time, height'},
'qis': { 'std_name': 'mass_mixing_ratio_of_snow_in_air', 'var_id': 'qis', 'units': 'kg kg-1', 'dim': 'time, height'},
'qig': { 'std_name': 'mass_mixing_ratio_of_graupel_in_air', 'var_id': 'qig', 'units': 'kg kg-1', 'dim': 'time, height'},
'nlc': { 'std_name': 'number_of_liquid_cloud_droplets_in_air', 'var_id': 'nlc', 'units': 'kg-1', 'dim': 'time, height'},
'nlr': { 'std_name': 'number_of_rain_drops_in_air', 'var_id': 'nlr', 'units': 'kg-1', 'dim': 'time, height'},
'nic': { 'std_name': 'number_of_cloud_ice_crystals_in_air', 'var_id': 'nic', 'units': 'kg-1', 'dim': 'time, height'},
'nis': { 'std_name': 'number_of_snow_crystals_in_air', 'var_id': 'nis', 'units': 'kg-1', 'dim': 'time, height'},
'nig': { 'std_name': 'number_of_graupel_crystals_in_air', 'var_id': 'nig', 'units': 'kg-1', 'dim': 'time, height'},
'na1': { 'std_name': 'number_of_total_aerosol_mode1', 'var_id': 'na1', 'units': 'kg-1', 'dim': 'time, height'},
'na2': { 'std_name': 'number_of_total_aerosol_mode2', 'var_id': 'na2', 'units': 'kg-1', 'dim': 'time, height'},
'na3': { 'std_name': 'number_of_total_aerosol_mode3', 'var_id': 'na3', 'units': 'kg-1', 'dim': 'time, height'},
'nlcic': { 'std_name': 'number_of_liquid_cloud_droplets_in_cloud', 'var_id': 'nlcic', 'units': 'kg-1', 'dim': 'time, height'},
'niic' : { 'std_name': 'number_of_ice_crystals_in_cloud', 'var_id': 'niic', 'units': 'kg-1', 'dim': 'time, height'}
}
Output_2D = {
'time' : {'std_name': 'time' , 'var_id': 'time' , 'units': 's'},
'x' : {'std_name': 'x' , 'var_id': 'x' , 'units': 'm'},
'y' : {'std_name': 'y' , 'var_id': 'y' , 'units': 'm'},
'hfss' : {'std_name': 'surface_upward_sensible_heat_flux' , 'var_id': 'hfss' , 'units': 'W m-2' , 'dimensions': ['time', 'x', 'y']},
'hfls' : {'std_name': 'surface_upward_latent_heat_flux' , 'var_id': 'hfls' , 'units': 'W m-2' , 'dimensions': ['time', 'x', 'y']},
'ustar': {'std_name': 'surface_friction_velocity' , 'var_id': 'ustar' , 'units': 'm s-1' , 'dimensions': ['time', 'x', 'y']},
'us' : {'std_name': 'surface_eastward_wind' , 'var_id': 'us' , 'units': 'm s-1' , 'dimensions': ['time', 'x', 'y']},
'vs' : {'std_name': 'surface_northward_wind' , 'var_id': 'vs' , 'units': 'm s-1' , 'dimensions': ['time', 'x', 'y']},
'pr' : {'std_name': 'precipitation_flux_at_surface' , 'var_id': 'pr' , 'units': 'kg m-2 s-1' , 'dimensions': ['time', 'x', 'y']},
'lwp' :{'std_name': 'atmosphere_mass_content_of_liquid_water' , 'var_id': 'lwp' , 'units': 'kg m-2' , 'dimensions': ['time', 'x', 'y']},
'iwp' :{'std_name': 'atmosphere_mass_content_of_ice' , 'var_id': 'iwp' , 'units': 'kg m-2' , 'dimensions': ['time', 'x', 'y']},
'opt' :{'std_name': 'atmosphere_optical_thickness' , 'var_id': 'opt' , 'units': ' ' , 'dimensions': ['time', 'x', 'y']},
'alb' :{'std_name': 'pseudo-albedo' , 'var_id': 'alb' , 'units': ' ' , 'dimensions': ['time', 'x', 'y']},
'olr11':{'std_name': 'toa_outgoing_longwave_flux_atmospheric_window', 'var_id': 'olr11' , 'units': 'K' , 'dimensions': ['time', 'x', 'y']}
}
Output_3D = {
'time': {'std_name': 'time' , 'var_id': 'time', 'units': 's'},
'x' : {'std_name': 'x' , 'var_id': 'x' , 'units': 'm'},
'y' : {'std_name': 'y' , 'var_id': 'y' , 'units': 'm'},
'zf' : {'std_name': 'height' , 'var_id': 'zf' , 'units': 'm', 'dimensions': ['time']},
'pa' : {'std_name': 'air_pressure' , 'var_id': 'pa' , 'units': 'Pa', 'dimensions': ['time', 'height']},
'ua' : {'std_name': 'eastward_wind' , 'var_id': 'ua' , 'units': 'm s^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'va' : {'std_name': 'northward_wind' , 'var_id': 'va' , 'units': 'm s^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'wa' : {'std_name': 'upward_air_velocity' , 'var_id': 'wa' , 'units': 'm s^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'ta' : {'std_name': 'air_temperature' , 'var_id': 'ta' , 'units': 'K', 'dimensions': ['time', 'height', 'x', 'y']},
'rhoa': {'std_name': 'air_dry_density' , 'var_id': 'rhoa', 'units': 'kg m^-3', 'dimensions': ['time', 'height', 'x', 'y']},
'qv' : {'std_name': 'water_vapor_mixing_ratio' , 'var_id': 'qv', 'units': 'kg kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'qlc' : {'std_name': 'mass_mixing_ratio_of_cloud_liquid_water_in_air' , 'var_id': 'qlc', 'units': 'kg kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'qlr' : {'std_name': 'mass_mixing_ratio_of_rain_water_in_air' , 'var_id': 'qlr', 'units': 'kg kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'qic' : {'std_name': 'mass_mixing_ratio_of_cloud_ice_in_air' , 'var_id': 'qic', 'units': 'kg kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'qis' : {'std_name': 'mass_mixing_ratio_of_snow_in_air' , 'var_id': 'qis', 'units': 'kg kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'qig' : {'std_name': 'mass_mixing_ratio_of_graupel_in_air' , 'var_id': 'qig', 'units': 'kg kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'na' : {'std_name': 'number_of_dry_aerosol_in_air' , 'var_id': 'na', 'units': 'kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'nlc' : {'std_name': 'number_of_cloud_droplets_in_air' , 'var_id': 'nlc', 'units': 'kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'nlr' : {'std_name': 'number_of_rain_droplets_in_air' , 'var_id': 'nlr', 'units': 'kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'nic' : {'std_name': 'number_of_cloud_ice_crystals_in_air' , 'var_id': 'nic', 'units': 'kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'nis' : {'std_name': 'number_of_snow_crystals_in_air' , 'var_id': 'nis', 'units': 'kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'nig' : {'std_name': 'number_of_graupel_crystals_in_air' , 'var_id': 'nig', 'units': 'kg^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'eps' : {'std_name': 'dissipation_rate_of_turbulent_kinetic_energy' , 'var_id': 'eps', 'units': 'm^2 s^-3', 'dimensions': ['time', 'height', 'x', 'y']},
'vmlc': {'std_name': 'mass_weighted_fall_speed_of_liquid_cloud_water' , 'var_id': 'vmlc', 'units': 'm s^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'vmlr': {'std_name': 'mass_weighted_fall_speed_of_rain' , 'var_id': 'vmlr', 'units': 'm s^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'vmic': {'std_name': 'mass_weighted_fall_speed_of_ice_cloud' , 'var_id': 'vmic', 'units': 'm s^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'vmis': {'std_name': 'mass_weighted_fall_speed_of_snow' , 'var_id': 'vmis', 'units': 'm s^-1', 'dimensions': ['time', 'height', 'x', 'y']},
'vmig': {'std_name': 'mass_weighted_fall_speed_of_graupel' , 'var_id': 'vmig', 'units': 'm s^-1', 'dimensions': ['time', 'height', 'x', 'y']}
}
def convert_TS_to_xarray(TS_output):
data=np.genfromtxt(TS_output, skip_header=1)
# Read the variable names from the first line
with open(TS_output, 'r') as file:
header = file.readline().strip()
variable_names = header.split()[2:] # Exclude the # symbol and the 'time'
# Extract time and data for each variable
time = data[:, 0]
variables_data = data[:,1:]
# Create an xarray dataset
ds = xr.Dataset()
# Add time as a dimension
ds['time'] = time
# Add variables as data variables
for ii, var_name in enumerate(variable_names):
ds[var_name] = (['time'], variables_data[:,ii])
file.close()
return(ds)
def Output_add_attributes(output_data,output_name,mixed_phase):
if mixed_phase:
Output_statistic_attr={'title': 'Mean output for '+output_name,
'reference': 'COMBLE',
'author': 'Alejandro Baro Perez (alejandro.baro.perez@misu.su.se)',
'version': 'Created on 2026-05-27',
'format_version': 'DEPHY SCM format version 1',}
else:
Output_statistic_attr={'case': 'COMBLE_CAO','title': 'Mean output for '+output_name,
'reference': 'COMBLE',
'author': 'Alejandro Baro Perez (alejandro.baro.perez@misu.su.se)',
'version': 'Created on 2026-05-27',
'format_version': 'DEPHY SCM format version 1',}
output_data.attrs=Output_statistic_attr
return(output_data)
def convert_time_new(output_xarray,yyear,mmonth,dday,hhour,time_format):
time_arr=output_xarray.time.values
if(time_format=='days'):
t_origin = datetime.datetime(year=yyear, month=mmonth, day=dday, hour=hhour)
date_list = [datetime.timedelta(days=x) + t_origin for x in time_arr[:].astype(np.float64)]
if(time_format=='seconds'):
t_origin = datetime.datetime(year=yyear, month=mmonth, day=dday, hour=hhour)
date_list = [datetime.timedelta(seconds=x) + t_origin for x in time_arr[:].astype(np.float64)]
# Set seconds and microseconds to zero
date_list = [dt.replace(second=0, microsecond=0) for dt in date_list]
date_list=np.asarray(date_list)
output_xarray=output_xarray.assign(time=np.asarray(date_list))
return(output_xarray)
#=============================================================================================================
# OUTPUT STATISTIC
# domain-mean quantities and statistics every 10 minutes (time-averages reported at the end of each period)
def create_stat_output(MIMICA_prof,MIMICA_prof_cl,MIMICA_surf,MIMICA_TS,output_name,mixed_phase):
#Ensuring that the time values are the same for prof,prof_cl,surf and TS outputs
MIMICA_TS=MIMICA_TS.assign_coords(time=MIMICA_surf.time.values)
MIMICA_prof=MIMICA_prof.assign_coords(time=MIMICA_surf.time.values)
MIMICA_prof_cl=MIMICA_prof_cl.assign_coords(time=MIMICA_surf.time.values)
Statistical_output = xr.Dataset()
#Adding the surface pressure (it is constant in MIMICA)
psurf=99544.5 #[Pa]
ps_array=np.full(MIMICA_surf.time.shape,psurf)
ps_dataset = xr.Dataset({"ps": ("time",ps_array,)},coords={"time": MIMICA_surf.time.values},)
Statistical_output['ps']=ps_dataset.ps
Statistical_output.ps.attrs['standard_name']='surface_pressure'
Statistical_output.ps.attrs['units']='Pa'
# Friction velocity, surface fluxes
DEPHY_name=['ts','ustar', 'z0','z0h','hfss','hfls', 'ol','rlut', 'rlus' , 'rlds']
MIMICA_name=['SST','U_star', 'Z0_m','z0_h','SHF','LHF','L_MO','OLR', 'Surface_LW_up','Surface_LW_down']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
Statistical_output[dephy]=MIMICA_TS[mimica]
Statistical_output[dephy].attrs['standard_name']= Output_statistic[dephy]['std_name']
Statistical_output[dephy].attrs['units']=Output_statistic[dephy]['units']
# Water paths
DEPHY_name=['lwpc','lwpr','iwp']
MIMICA_name=['LWPc','LWPr', 'IWP']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
Statistical_output[dephy]=MIMICA_TS[mimica]
Statistical_output[dephy].attrs['standard_name']= Output_statistic[dephy]['std_name']
Statistical_output[dephy].attrs['units']=Output_statistic[dephy]['units']
#Cloud area fraction
Statistical_output['clt']=MIMICA_TS['C_COL']
Statistical_output['clt'].attrs['standard_name']= Output_statistic['clt']['std_name']
Statistical_output['clt'].attrs['units']=Output_statistic['clt']['units']
#Optical depths
Statistical_output['od'] =(MIMICA_surf['OPTC']+MIMICA_surf['OPTI']+MIMICA_surf['OPTR']).mean(dim=('Y','X'))
Statistical_output['od'].attrs['standard_name']= Output_statistic['od']['std_name']
Statistical_output['od'].attrs['units']=Output_statistic['od']['units']
Statistical_output['odlc'] =MIMICA_surf['OPTC'].mean(dim=('Y','X'))
Statistical_output['odlc'].attrs['standard_name']= Output_statistic['odlc']['std_name']
Statistical_output['odlc'].attrs['units']=Output_statistic['odlc']['units']
# Surface flux of aerosols MIMICA
if 'Surf_naer_flux' in MIMICA_surf:
Statistical_output['ssaf']=MIMICA_surf['Surf_naer_flux'].mean(dim=('Y','X'))
Statistical_output['ssaf'].attrs['standard_name']= Output_statistic['ssaf']['std_name']
Statistical_output['ssaf'].attrs['units']= Output_statistic['ssaf']['units']
# U,V,temperature,potential temperature,water vapor mixing fraction(ratio)
DEPHY_name=['pa', 'ua', 'va', 'rhoa','ta','theta','tke_res','tke_sgs']
MIMICA_name=['P_tot','U', 'V','RHO','T','PT' , 'TKE' , 'TKE_sgs']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
Statistical_output[dephy]=MIMICA_prof[mimica]
Statistical_output[dephy].attrs['standard_name']= Output_statistic[dephy]['std_name']
Statistical_output[dephy].attrs['units']= Output_statistic[dephy]['units']
## Precipitation flux. Change from precipitation rate (MIMICA) to precipitation flux
#Considering that
# density of water = 997 kg/m3
# density of ice = 917 kg/m3
# precipitation rates in MIMICA are given in = mm/d = 1e-3/86400 [m/s]
# Mulitplying the density of water or ice by the precipitation rates in MIMICA we can get the precipitation fluxes
# in units [kg m-2 s-1]
# Then precipitation rates in MIMICA should be multiplied by
# a factor = 0.997/86400 if water
# a factor = 0.917/86400 if ice
water_factor= 0.997/86400
ice_factor= 0.917/86400
DEPHY_name=['pr','pri']
MIMICA_name=['R_rate','R_ice_rate']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
if mixed_phase:
#Assuming that the liquid phase precipitation in negligible in the mixed phase case
Statistical_output[dephy]=ice_factor*MIMICA_TS[mimica]
else:
Statistical_output[dephy]=water_factor*MIMICA_TS[mimica]
Statistical_output[dephy].attrs['standard_name']= Output_statistic[dephy]['std_name']
Statistical_output[dephy].attrs['units']= Output_statistic[dephy]['units']
# Relative humidity, and relative humidity over ice:
DEPHY_name=['hur','huri']
MIMICA_name=['RH','RHi']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
Statistical_output[dephy]=1e-2*MIMICA_prof[mimica]
Statistical_output[dephy].attrs['standard_name']= Output_statistic[dephy]['std_name']
Statistical_output[dephy].attrs['units']= Output_statistic[dephy]['units']
# Mixing fraction from MIMICA are directly mass mixing ratio in DEPHY
DEPHY_name=[ 'qv','qlc','qlr','qic','qis','qig']
MIMICA_name=['Qv','Qc', 'Qr', 'Qi', 'Qs', 'Qg']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
if mimica in MIMICA_prof:
Statistical_output[dephy]=MIMICA_prof[mimica]/MIMICA_prof['RHO']
Statistical_output[dephy].attrs['standard_name']= Output_statistic[dephy]['std_name']
Statistical_output[dephy].attrs['units']= Output_statistic[dephy]['units']
# Hydrometeors number density:
DEPHY_name=['nlc','nlr','nic','nis','nig']
MIMICA_name=['Nc','Nr','Ni','Ns','Ng']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
if mimica in MIMICA_prof:
Statistical_output[dephy]=MIMICA_prof[mimica]/MIMICA_prof['RHO']
Statistical_output[dephy].attrs['standard_name']= Output_statistic[dephy]['std_name']
Statistical_output[dephy].attrs['units']= Output_statistic[dephy]['units']
# Hydrometeors number density in cloud:
DEPHY_name=['nlcic','niic']
MIMICA_name=['Nc','Ni']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
if mimica in MIMICA_prof_cl:
#Statistical_output[dephy]=MIMICA_prof_cl[mimica]/MIMICA_prof_cl['RHO']
Statistical_output[dephy]=save_division(MIMICA_prof_cl[mimica],MIMICA_prof_cl['RHO'])
Statistical_output[dephy].attrs['standard_name']= Output_statistic[dephy]['std_name']
Statistical_output[dephy].attrs['units']= Output_statistic[dephy]['units']
# Number of total aerosol in air #For later !activated + unactivated???
if 'N_AERO1' in MIMICA_prof:
DEPHY_name=['na1','na2','na3']
MIMICA_name=['N_AERO1','N_AERO2','N_AERO3']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
Statistical_output[dephy]=MIMICA_prof[mimica]
Statistical_output[dephy].attrs['units']= Output_statistic[dephy]['units']
# Renaming vertical coordinate
Statistical_output=Statistical_output.rename({'Z': 'zf'})
Statistical_output.zf.attrs['standard_name']='height'
Statistical_output.zf.attrs['comment']='altitude of mid level points above sea surface'
# Adding the global attributes to the file
Statistical_output=Output_add_attributes(Statistical_output,output_name,mixed_phase)
# Changing time format
Statistical_output=convert_time_new(Statistical_output,2020,3,12,22,'seconds')
return Statistical_output
def create_2D_output(MIMICA_surf,output_name,mixed_phase):
D2_output = xr.Dataset()
# Friction velocity, surface fluxes
DEPHY_name=['ustar','hfss','hfls']
MIMICA_name=['U_star','SHF','LHF']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
D2_output[dephy]=MIMICA_surf[mimica]
D2_output[dephy].attrs['standard_name']= Output_2D[dephy]['std_name']
D2_output[dephy].attrs['units']= Output_2D[dephy]['units']
# Liquid water path, ice water path
DEPHY_name=['lwp','iwp']
MIMICA_name=['LWP','IWP']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
D2_output[dephy]=MIMICA_surf[mimica]
D2_output[dephy].attrs['standard_name']= Output_2D[dephy]['std_name']
D2_output[dephy].attrs['units']= Output_2D[dephy]['units']
# Precipitation flux. Change from precipitation rate (MIMICA) to precipitation flux
water_factor= 0.997/86400
ice_factor= 0.917/86400
if mixed_phase:
D2_output['pr']=ice_factor*MIMICA_surf['R_rate']
else:
D2_output['pr']=water_factor*MIMICA_surf['R_rate']
D2_output['pr'].attrs['standard_name']= Output_2D['pr']['std_name']
D2_output['pr'].attrs['units']= Output_2D['pr']['units']
# Atmosphere optical thickness
OPT= MIMICA_surf['OPTC']+MIMICA_surf['OPTI']+MIMICA_surf['OPTR']
D2_output['opt'] = OPT
D2_output['opt'].attrs['standard_name']= Output_2D['opt']['std_name']
D2_output['opt'].attrs['units']= Output_2D['opt']['units']
#Pseudo-albedo
D2_output['alb'] = OPT/(OPT+13.)
D2_output['alb'].attrs['standard_name']= Output_2D['alb']['std_name']
D2_output['alb'].attrs['units']= Output_2D['alb']['units']
# Outgoing longwave radiation in atmospheric window
D2_output['olr11'] = MIMICA_surf['OLR_11mum']
D2_output['olr11'].attrs['starndard_name']= Output_2D['olr11']['std_name']
D2_output['olr11'].attrs['units']= Output_2D['olr11']['units']
D2_output['olr11'].attrs['description']= 'Outgoing longwave flux in window containing 10.8 micrometers (wavelength range 10-12.5 micrometers)'
# Surface winds u,v
D2_output['us'] = MIMICA_surf['U_zlevel_1']
D2_output['us'].attrs['starndard_name']= Output_2D['us']['std_name']
D2_output['us'].attrs['units']= Output_2D['us']['units']
D2_output['vs'] = MIMICA_surf['V_zlevel_1']
D2_output['vs'].attrs['starndard_name']= Output_2D['vs']['std_name']
D2_output['vs'].attrs['units']= Output_2D['vs']['units']
# Changing time format
D2_output=convert_time_new(D2_output,2020,3,12,22,'seconds')
# Adding the global attributes to the file
D2_output=Output_add_attributes(D2_output,output_name,mixed_phase)
return D2_output
def create_3D_output(MIMICA_3D,output_name,mixed_phase):
D3_output = xr.Dataset()
# U,V,temperature,potential temperature,water vapor mixing fraction(ratio)
DEPHY_name=[ 'pa','ua','va','wa','ta','rhoa']
MIMICA_name=['P_tot','U', 'V','W', 'T','RHO']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
D3_output[dephy]=MIMICA_3D[mimica]
D3_output[dephy].attrs['standard_name']= Output_3D[dephy]['std_name']
D3_output[dephy].attrs['units']= Output_3D[dephy]['units']
# Mixing fraction from MIMICA are directly mass mixing ratio in DEPHY
DEPHY_name=[ 'qv','qlc','qlr']
MIMICA_name=['Qv','Qc', 'Qr' ]
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
D3_output[dephy]=MIMICA_3D[mimica]
D3_output[dephy].attrs['standard_name']= Output_3D[dephy]['std_name']
D3_output[dephy].attrs['units']= Output_3D[dephy]['units']
# Mixing fraction from MIMICA are directly mass mixing ratio in DEPHY [ice phase]
if 'Qi' in MIMICA_3D:
DEPHY_name=[ 'qic','qis','qig']
MIMICA_name=['Qi', 'Qs', 'Qg' ]
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
D3_output[dephy]=MIMICA_3D[mimica]
D3_output[dephy].attrs['standard_name']= Output_3D[dephy]['std_name']
D3_output[dephy].attrs['units']= Output_3D[dephy]['units']
# Hydrometeors number density: #Question here
DEPHY_name=['nlc','nlr']
MIMICA_name=['Nc','Nr']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
D3_output[dephy]=MIMICA_3D[mimica]
D3_output[dephy].attrs['standard_name']= Output_3D[dephy]['std_name']
D3_output[dephy].attrs['units']= Output_3D[dephy]['units']
# Hydrometeors number density (ice phase)
if 'Ni' in MIMICA_3D:
DEPHY_name=['nic','nis','nig']
MIMICA_name=['Ni','Ns','Ng']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
D3_output[dephy]=MIMICA_3D[mimica]
D3_output[dephy].attrs['standard_name']= Output_3D[dephy]['std_name']
D3_output[dephy].attrs['units']= Output_3D[dephy]['units']
#Check these variables. Are they mass weighted?
DEPHY_name=[ 'vmlc','vmlr','vmic','vmis','vmig']
MIMICA_name=[ 'Vpc','Vpr','Vpi','Vps','Vpg']
for dephy,mimica in zip(DEPHY_name,MIMICA_name):
if mimica in MIMICA_3D:
D3_output[dephy]=MIMICA_3D[mimica]
D3_output[dephy].attrs['standard_name']= Output_3D[dephy]['std_name']
D3_output[dephy].attrs['units']= Output_3D[dephy]['units']
# Changing time format
D3_output=convert_time_new(D3_output,2020,3,12,22,'seconds')
# Adding the global attributes to the file
D3_output=Output_add_attributes(D3_output,output_name,mixed_phase)
return D3_output
def save_division(numerator,denominator):
# Replace zeros in the denominator with NaN to avoid division by zero
safe_denominator = denominator.where(denominator != 0, np.nan)
# Perform the division
resulting_variable = numerator / safe_denominator
# Replace NaN values in the result with zero
resulting_variable = resulting_variable.fillna(0)
# Assign the resulting variable back to the dataset if needed
#COMBLE_FixN_256_new_prof_cl['resulting_variable'] = resulting_variable
return resulting_variable
def convert_time_from_minutes_to_seconds_new(file, file_time_modified):
# Open dataset
ds = xr.open_dataset(file)
# Force correct encoding
ds.time.encoding["units"] = "seconds since 2020-03-12 22:00:00"
ds.time.encoding["calendar"] = "proleptic_gregorian"
# Save modified dataset
ds.to_netcdf(file_time_modified)
return ds
Simulation with ice part 1: FixN#
#Reading MIMICA native outputs
COMBLE_FixN_256_prof= xr.open_dataset(native_outputs_path+'COMBLE_FixN_256x256/OUTPUT/profiles_tot.nc')
COMBLE_FixN_256_prof_cl=xr.open_dataset(native_outputs_path+'COMBLE_FixN_256x256/OUTPUT/profiles_cl.nc')
COMBLE_FixN_256_surf= xr.open_dataset(native_outputs_path+'COMBLE_FixN_256x256/OUTPUT/slice_surf.nc')
COMBLE_FixN_256_3D= xr.open_dataset(native_outputs_path+'COMBLE_FixN_256x256/OUTPUT/COMBLE_MPI.nc')
COMBLE_FixN_256_TS= xr.open_dataset(native_outputs_path+'COMBLE_FixN_256x256/OUTPUT/T_S.nc')
Stat_output_FixN_256=create_stat_output(COMBLE_FixN_256_prof,COMBLE_FixN_256_prof_cl,COMBLE_FixN_256_surf,
COMBLE_FixN_256_TS,'MIMICA_Lx25_dx100_FixN_alt',True)
OUTPUT_2D_FixN_256=create_2D_output(COMBLE_FixN_256_surf,'MIMICA_Lx25_dx100_FixN_2D_alt',True)
OUTPUT_3D_FixN_256=create_3D_output(COMBLE_FixN_256_3D,'MIMICA_Lx25_dx100_FixN_3D_alt',True)
Create netcdf files for FixN#
Stat_output_FixN_256.to_netcdf(path=converted_outputs_path+'MIMICA_Lx25_dx100_FixN'+'.nc',mode='w',engine="netcdf4")
OUTPUT_2D_FixN_256.to_netcdf(path=converted_outputs_path+'MIMICA_Lx25_dx100_FixN_2D'+'.nc',mode='w',engine="netcdf4")
OUTPUT_3D_FixN_256.to_netcdf(path=converted_outputs_path+'MIMICA_Lx25_dx100_FixN_3D'+'.nc',mode='w',engine="netcdf4")
Simulation without ice part 1 : FixN noice#
COMBLE_FixN_noice_256_prof= xr.open_dataset(native_outputs_path+'COMBLE_FixN_noice_256x256/OUTPUT/profiles_tot.nc')
COMBLE_FixN_noice_256_prof_cl=xr.open_dataset(native_outputs_path+'COMBLE_FixN_noice_256x256/OUTPUT/profiles_cl.nc')
COMBLE_FixN_noice_256_surf= xr.open_dataset(native_outputs_path+'COMBLE_FixN_noice_256x256/OUTPUT/slice_surf.nc')
COMBLE_FixN_noice_256_3D= xr.open_dataset(native_outputs_path+'COMBLE_FixN_noice_256x256/OUTPUT/COMBLE_MPI.nc')
COMBLE_FixN_noice_256_TS= xr.open_dataset(native_outputs_path+'COMBLE_FixN_noice_256x256/OUTPUT/T_S.nc')
Stat_output_FixN_noice_256=create_stat_output(COMBLE_FixN_noice_256_prof,COMBLE_FixN_noice_256_prof_cl,COMBLE_FixN_noice_256_surf,
COMBLE_FixN_noice_256_TS,'MIMICA_Lx25_dx100_FixN_noice',False)
OUTPUT_2D_FixN_noice_256=create_2D_output(COMBLE_FixN_noice_256_surf,'MIMICA_Lx25_dx100_FixN_noice_2D',False)
OUTPUT_3D_FixN_noice_256=create_3D_output(COMBLE_FixN_noice_256_3D,'MIMICA_Lx25_dx100_FixN_noice_3D',False)
Create netcdf files for FixN noice#
Stat_output_FixN_noice_256.to_netcdf(path=converted_outputs_path+'MIMICA_Lx25_dx100_FixN_noice'+'.nc',mode='w',engine="netcdf4")
OUTPUT_2D_FixN_noice_256.to_netcdf(path=converted_outputs_path+'MIMICA_Lx25_dx100_FixN_noice_2D'+'.nc',mode='w',engine="netcdf4")
OUTPUT_3D_FixN_noice_256.to_netcdf(path=converted_outputs_path+'MIMICA_Lx25_dx100_FixN_noice_3D'+'.nc',mode='w',engine="netcdf4")
Taking the last time step of the 3D outputs for the first iteration:#
OUTPUT_3D_FixN_256.isel(time=20).to_netcdf(path=converted_outputs_path+'MIMICA_Lx25_dx100_FixN_3D_at_18UTC'+'.nc',mode='w',engine="netcdf4")
OUTPUT_3D_FixN_noice_256.isel(time=20).to_netcdf(path=converted_outputs_path+'MIMICA_Lx25_dx100_FixN_noice_3D_at_18UTC'+'.nc',mode='w',engine="netcdf4")
Converting the time from minutes to seconds#
#Doing this only because Tim Juliano could not read properly the netcdf files after the conversion
file_names=[ 'MIMICA_Lx25_dx100_FixN.nc','MIMICA_Lx25_dx100_FixN_2D.nc', 'MIMICA_Lx25_dx100_FixN_3D.nc',
'MIMICA_Lx25_dx100_FixN_noice.nc','MIMICA_Lx25_dx100_FixN_noice_2D.nc', 'MIMICA_Lx25_dx100_FixN_noice_3D.nc' ]
for ii_file in file_names:
print("Converting time in file ",ii_file, "saving in",converted_outputs_path_time_mod+ii_file)
convert_time_from_minutes_to_seconds_new(converted_outputs_path+ii_file,converted_outputs_path_time_mod+ii_file)
Converting time in file MIMICA_Lx25_dx100_FixN_noice.nc saving in /proj/bolinc/users/x_abaro/MIMICA_radiation/COMBLE_CONVERTED_time_mod/MIMICA_Lx25_dx100_FixN_noice.nc
Converting time in file MIMICA_Lx25_dx100_FixN_noice_2D.nc saving in /proj/bolinc/users/x_abaro/MIMICA_radiation/COMBLE_CONVERTED_time_mod/MIMICA_Lx25_dx100_FixN_noice_2D.nc
Converting time in file MIMICA_Lx25_dx100_FixN_noice_3D.nc saving in /proj/bolinc/users/x_abaro/MIMICA_radiation/COMBLE_CONVERTED_time_mod/MIMICA_Lx25_dx100_FixN_noice_3D.nc
Part 2#
#COMBLE_Part2_Prog_Na_prof=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNa_new/OUTPUT/profiles_tot.nc')
#COMBLE_Part2_Prog_Na_prof_cl=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNa_new/OUTPUT/profiles_cl.nc')
#COMBLE_Part2_Prog_Na_surf=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNa_new/OUTPUT/slice_surf.nc')
#COMBLE_Part2_Prog_Na_3D=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNa_new/OUTPUT/COMBLE_MPI.nc')
#COMBLE_Part2_Prog_Na_TS=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNa_new/OUTPUT/T_S.nc')
#Stat_output_COMBLE_Part2_Prog_Na=create_stat_output(COMBLE_Part2_Prog_Na_prof,COMBLE_Part2_Prog_Na_prof_cl,COMBLE_Part2_Prog_Na_surf,
# COMBLE_Part2_Prog_Na_TS,'MIMICA_Lx25_dx100_Prog_Na_alt',True)
#OUTPUT_2D_COMBLE_Part2_Prog_Na=create_2D_output(COMBLE_Part2_Prog_Na_surf,'MIMICA_Lx25_dx100_Prog_Na_2D_alt',True)
#OUTPUT_3D_COMBLE_Part2_Prog_Na=create_3D_output(COMBLE_Part2_Prog_Na_3D,'MIMICA_Lx25_dx100_Prog_Na_3D_alt',True)
#COMBLE_Part2_Prog_NaNi_prof=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNaNi/OUTPUT/profiles_tot.nc')
#COMBLE_Part2_Prog_NaNi_prof_cl=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNaNi/OUTPUT/profiles_cl.nc')
#COMBLE_Part2_Prog_NaNi_surf=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNaNi/OUTPUT/slice_surf.nc')
#COMBLE_Part2_Prog_NaNi_3D=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNaNi/OUTPUT/COMBLE_MPI.nc')
#COMBLE_Part2_Prog_NaNi_TS=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNaNi/OUTPUT/T_S.nc')
#Stat_output_COMBLE_Part2_Prog_NaNi=create_stat_output(COMBLE_Part2_Prog_NaNi_prof,COMBLE_Part2_Prog_NaNi_prof_cl,COMBLE_Part2_Prog_NaNi_surf,
# COMBLE_Part2_Prog_NaNi_TS,'MIMICA_Lx25_dx100_Prog_NaNi_alt',True)
#OUTPUT_2D_COMBLE_Part2_Prog_NaNi=create_2D_output(COMBLE_Part2_Prog_NaNi_surf,'MIMICA_Lx25_dx100_Prog_NaNi_2D_alt',True)
#OUTPUT_3D_COMBLE_Part2_Prog_NaNi=create_3D_output(COMBLE_Part2_Prog_NaNi_3D,'MIMICA_Lx25_dx100_Prog_NaNi_3D_alt',True)
#COMBLE_Part2_Prog_Na_noice_prof=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNa_noice/OUTPUT/profiles_tot.nc')
#COMBLE_Part2_Prog_Na_noice_prof_cl=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNa_noice/OUTPUT/profiles_cl.nc')
#COMBLE_Part2_Prog_Na_noice_surf=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNa_noice/OUTPUT/slice_surf.nc')
#COMBLE_Part2_Prog_Na_noice_3D=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNa_noice/OUTPUT/COMBLE_MPI.nc')
#COMBLE_Part2_Prog_Na_noice_TS=xr.open_dataset(COMBLE_part2+'COMBLE_Part2_ProgNa_noice/OUTPUT/T_S.nc')
#Stat_output_COMBLE_Part2_Prog_Na_noice=create_stat_output(COMBLE_Part2_Prog_Na_noice_prof,COMBLE_Part2_Prog_Na_noice_prof_cl,COMBLE_Part2_Prog_Na_noice_surf,
# COMBLE_Part2_Prog_Na_noice_TS,'MIMICA_Lx25_dx100_Prog_Na_noice',False)
#OUTPUT_2D_COMBLE_Part2_Prog_Na_noice=create_2D_output(COMBLE_Part2_Prog_Na_noice_surf,'MIMICA_Lx25_dx100_Prog_Na_noice_2D',False)
#OUTPUT_3D_COMBLE_Part2_Prog_Na_noice=create_3D_output(COMBLE_Part2_Prog_Na_noice_3D,'MIMICA_Lx25_dx100_Prog_Na_noice_3D',False)
Create netcdf files#
Prog_Na#
#file_path = simulation_path+'MIMICA_Lx25_dx100_Prog_Na_alt'+'.nc'
#
#if os.path.exists(file_path):
# os.remove(file_path)
# print(f"The file {file_path} has been deleted.")
#
#Stat_output_COMBLE_Part2_Prog_Na.to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_Na_alt'+'.nc',mode='w',engine="netcdf4")
The file /proj/mimica-chalmers-comble/users/x_abaro/COMBLE/Simulations_Folder_COMBLE/MIMICA_Lx25_dx100_Prog_Na_alt.nc has been deleted.
#file_path = simulation_path+'MIMICA_Lx25_dx100_Prog_Na_2D_alt'+'.nc'
#
#if os.path.exists(file_path):
# os.remove(file_path)
# print(f"The file {file_path} has been deleted.")
#
#OUTPUT_2D_COMBLE_Part2_Prog_Na.to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_Na_2D_alt'+'.nc',mode='w',engine="netcdf4")
The file /proj/mimica-chalmers-comble/users/x_abaro/COMBLE/Simulations_Folder_COMBLE/MIMICA_Lx25_dx100_Prog_Na_2D_alt.nc has been deleted.
#file_path = simulation_path+'MIMICA_Lx25_dx100_Prog_Na_3D_alt'+'.nc'
#
#if os.path.exists(file_path):
# os.remove(file_path)
# print(f"The file {file_path} has been deleted.")
#
#OUTPUT_3D_COMBLE_Part2_Prog_Na.to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_Na_3D_alt'+'.nc',mode='w',engine="netcdf4")
#OUTPUT_3D_COMBLE_Part2_Prog_Na.isel(time=20).to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_Na_3D_alt_at_18UTC'+'.nc',mode='w',engine="netcdf4")
Prog_NaNi#
#file_path = simulation_path+'MIMICA_Lx25_dx100_Prog_NaNi_alt'+'.nc'
#
#if os.path.exists(file_path):
# os.remove(file_path)
# print(f"The file {file_path} has been deleted.")
#
#Stat_output_COMBLE_Part2_Prog_NaNi.to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_NaNi_alt'+'.nc',mode='w',engine="netcdf4")
#file_path = simulation_path+'MIMICA_Lx25_dx100_Prog_NaNi_2D_alt'+'.nc'
#
#if os.path.exists(file_path):
# os.remove(file_path)
# print(f"The file {file_path} has been deleted.")
#
#OUTPUT_2D_COMBLE_Part2_Prog_NaNi.to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_NaNi_2D_alt'+'.nc',mode='w',engine="netcdf4")
#file_path = simulation_path+'MIMICA_Lx25_dx100_Prog_NaNi_3D_alt'+'.nc'
#
#if os.path.exists(file_path):
# os.remove(file_path)
# print(f"The file {file_path} has been deleted.")
#
#OUTPUT_3D_COMBLE_Part2_Prog_NaNi.to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_NaNi_3D_alt'+'.nc',mode='w',engine="netcdf4")
#OUTPUT_3D_COMBLE_Part2_Prog_NaNi.isel(time=20).to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_NaNi_3D_alt_at_18UTC'+'.nc',mode='w',engine="netcdf4")
Prog_Na_noice#
#file_path = simulation_path+'MIMICA_Lx25_dx100_Prog_Na_noice'+'.nc'
#
#if os.path.exists(file_path):
# os.remove(file_path)
# print(f"The file {file_path} has been deleted.")
#
#Stat_output_COMBLE_Part2_Prog_Na_noice.to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_Na_noice'+'.nc',mode='w',engine="netcdf4")
#file_path = simulation_path+'MIMICA_Lx25_dx100_Prog_Na_noice_2D'+'.nc'
#
#if os.path.exists(file_path):
# os.remove(file_path)
# print(f"The file {file_path} has been deleted.")
#
#OUTPUT_2D_COMBLE_Part2_Prog_Na_noice.to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_Na_noice_2D'+'.nc',mode='w',engine="netcdf4")
#file_path = simulation_path+'MIMICA_Lx25_dx100_Prog_Na_noice_3D'+'.nc'
#
#if os.path.exists(file_path):
# os.remove(file_path)
# print(f"The file {file_path} has been deleted.")
#
#OUTPUT_3D_COMBLE_Part2_Prog_Na_noice.to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_Na_noice_3D'+'.nc',mode='w',engine="netcdf4")
#OUTPUT_3D_COMBLE_Part2_Prog_Na_noice.isel(time=20).to_netcdf(path=simulation_path+'MIMICA_Lx25_dx100_Prog_Na_noice_3D_at_18UTC'+'.nc',mode='w',engine="netcdf4")
def convert_time_from_minutes_to_seconds(file,file_time_modified):
# Open dataset
ds = xr.open_dataset(file) # Replace with actual file
# Check if 'units' attribute exists
if 'units' in ds['time'].attrs:
# Extract the reference time from the attribute
import re
match = re.search(r"since\s(.+)", ds['time'].attrs['units'])
if match:
reference_time_str = match.group(1)
reference_time = np.datetime64(reference_time_str)
else:
raise ValueError("Could not extract reference time from units attribute.")
else:
# If units attribute is missing, define a reference time manually
print("Warning: No 'units' attribute found. Setting reference time manually.")
reference_time = ds['time'].values[0] # Use the first time entry as reference
reference_time_str = str(reference_time) # Convert to string
# Ensure time is in datetime64 format
ds['time'] = reference_time + (ds['time'] - reference_time) * 60 # Convert from minutes to seconds
# Update metadata to reflect seconds instead of minutes
ds['time'].attrs['units'] = f"seconds since {reference_time_str}"
# Move 'units' from attributes to encoding
ds['time'].encoding['units'] = ds['time'].attrs.pop('units')
# Save the updated dataset
ds.to_netcdf(file_time_modified)
# Print the updated time values and attributes
#print(ds['time'].values) # Ensure values are still datetime64[ns]
#print(ds['time'].attrs) # Check that units now say "seconds since ..."
#return ds
def convert_time_from_minutes_to_seconds_new(file, file_time_modified):
# Open dataset
ds = xr.open_dataset(file)
# Force correct encoding
ds.time.encoding["units"] = "seconds since 2020-03-12 22:00:00"
ds.time.encoding["calendar"] = "proleptic_gregorian"
# Save modified dataset
ds.to_netcdf(file_time_modified)
return ds
#file_names=['MIMICA_Lx25_dx100_FixN_alt.nc','MIMICA_Lx25_dx100_FixN_noice.nc','MIMICA_Lx25_dx100_Prog_Na_alt.nc','MIMICA_Lx25_dx100_Prog_NaNi_alt.nc',
# 'MIMICA_Lx25_dx100_Prog_Na_noice.nc']
file_names=[ 'MIMICA_Lx25_dx100_FixN.nc','MIMICA_Lx25_dx100_FixN_2D.nc', 'MIMICA_Lx25_dx100_FixN_3D.nc' ]
for ii_file in file_names:
print("Converting time in file ",ii_file, "saving in",converted_outputs_path_time_mod+ii_file)
convert_time_from_minutes_to_seconds_new(converted_outputs_path+ii_file,converted_outputs_path_time_mod+ii_file)
# print(" ")
Converting time in file MIMICA_Lx25_dx100_FixN.nc saving in /proj/bolinc/users/x_abaro/MIMICA_radiation/COMBLE_CONVERTED_time_mod/MIMICA_Lx25_dx100_FixN.nc
Converting time in file MIMICA_Lx25_dx100_FixN_2D.nc saving in /proj/bolinc/users/x_abaro/MIMICA_radiation/COMBLE_CONVERTED_time_mod/MIMICA_Lx25_dx100_FixN_2D.nc
Converting time in file MIMICA_Lx25_dx100_FixN_3D.nc saving in /proj/bolinc/users/x_abaro/MIMICA_radiation/COMBLE_CONVERTED_time_mod/MIMICA_Lx25_dx100_FixN_3D.nc