diff --git a/.github/workflows/marscalendar_test.yml b/.github/workflows/marscalendar_test.yml index 216c8cc8..4dc135b0 100644 --- a/.github/workflows/marscalendar_test.yml +++ b/.github/workflows/marscalendar_test.yml @@ -20,6 +20,10 @@ on: - 'bin/MarsCalendar.py' - 'tests/test_marscalendar.py' - '.github/workflows/marscalendar_test.yml' + # - 'AmesCAP/amescap/Script_utils.py' + # - 'AmesCAP/amescap/Ncdf_wrapper.py' + # - 'AmesCAP/amescap/FV3_utils.py' + # - 'AmesCAP/amescap/Spectral_utils.py' jobs: test: # Run on multiple OS and Python versions for comprehensive testing diff --git a/amescap/FV3_utils.py b/amescap/FV3_utils.py index 8acf91fd..60fd98fb 100644 --- a/amescap/FV3_utils.py +++ b/amescap/FV3_utils.py @@ -956,84 +956,12 @@ def cart_to_azimut_TR(u, v, mode="from"): np.sqrt(u**2 + v**2)) -def sfc_area_deg(lon1, lon2, lat1, lat2, R=3390000.): - """ - Return the surface between two sets of latitudes/longitudes:: - - S = int[R^2 dlon cos(lat) dlat] _____lat2 - \ \ - \____\lat1 - lon1 lon2 - :param lon1: longitude from set 1 [°] - :type lon1: float - :param lon2: longitude from set 2 [°] - :type lon2: float - :param lat1: latitude from set 1 [°] - :type lat1: float - :param lat2: longitude from set 2 [°] - :type lat2: float - :param R: planetary radius [m] - :type R: int - - .. note:: - qLon and Lat define the corners of the area not the grid cell center. - """ - - lat1 *= np.pi/180 - lat2 *= np.pi/180 - lon1 *= np.pi/180 - lon2 *= np.pi/180 - return ((R**2) - * np.abs(lon1 - lon2) - * np.abs(np.sin(lat1) - np.sin(lat2))) - - -def area_meridional_cells_deg(lat_c, dlon, dlat, normalize=False, R=3390000.): - """ - Return area of invidual cells for a meridional band of thickness - ``dlon`` where ``S = int[R^2 dlon cos(lat) dlat]`` with - ``sin(a)-sin(b) = 2 cos((a+b)/2)sin((a+b)/2)`` - so ``S = 2 R^2 dlon 2cos(lat)sin(dlat/2)``:: - - _________lat + dlat/2 - \ lat \ ^ - \lon + \ | dlat - \________\lat - dlat/2 v - lon - dlon/2 lon + dlon/2 - <------> - dlon - - :param lat_c: latitude of cell center [°] - :type lat_c: float - :param dlon: cell angular width [°] - :type dlon: float - :param dlat: cell angular height [°] - :type dlat: float - :param R: planetary radius [m] - :type R: float - :param normalize: if True, the sum of the output elements = 1 - :type normalize: bool - :return: ``S`` areas of the cells, same size as ``lat_c`` in [m2] - or normalized by the total area - """ - - # Initialize - area_tot = 1. - # Compute total area in a longitude band extending from - # ``lat[0] - dlat/2`` to ``lat_c[-1] + dlat/2`` - if normalize: - area_tot = sfc_area_deg(-dlon/2, dlon/2, (lat_c[0] - dlat/2), - (lat_c[-1] + dlat/2), R) - # Now convert to radians - lat_c = lat_c*np.pi/180 - dlon *= np.pi/180 - dlat *= np.pi/180 - return (2. * R**2 * dlon * np.cos(lat_c) * np.sin(dlat/2.) / area_tot) - - def area_weights_deg(var_shape, lat_c, axis = -2): """ - Return weights for averaging the variable. + Returns weights scaled so that np.mean(var*W) gives an area-weighted + average. This works because grid cells near the poles have smaller + areas than those at the equator, so they should contribute less to + a global average. Expected dimensions are: @@ -1048,12 +976,11 @@ def area_weights_deg(var_shape, lat_c, axis = -2): may be truncated on either end (e.g., ``lat = [-20 ..., 0... 50]``) but must be continuous. - :param var_shape: variable shape + :param var_shape: the shape/dimensions of your data array :type var_shape: tuple - :param lat_c: latitude of cell centers [°] + :param lat_c: latitude cell centers in degrees [°] :type lat_c: float - :param axis: position of the latitude axis for 2D and higher - dimensional arrays. The default is the SECOND TO LAST dimension + :param axis: which dimension contains latitude, default: 2nd-to-last :type axis: int :return: ``W`` weights for the variable ready for standard averaging as ``np.mean(var*W)`` [condensed form] or @@ -1083,33 +1010,75 @@ def area_weights_deg(var_shape, lat_c, axis = -2): ``np.average(var, weights=W, axis = X)`` to average over a specific axis. """ - + # Check if either the lat array or var shape is essentially + # scalar (single value) + # np.atleast_1d() ensures the input is treated as at least a 1D + # array for checking its length if (len(np.atleast_1d(lat_c)) == 1 or len(np.atleast_1d(var_shape)) == 1): - # If variable or lat is a scalar, do nothing + # If either is scalar, returns an array of 1s matching the var + # shape (no weighting needed since there's nothing to weight across) return np.ones(var_shape) else: - # Then lat has at least 2 elements + # Calculates lat spacing by taking the difference between the + # first two latitude points. Assumes uniform grid spacing dlat = lat_c[1]-lat_c[0] - # Calculate cell areas. Since it is normalized, we can use - # ``dlon = 1`` and ``R = 1`` without changing the result. Note - # that ``sum(A) = (A1 + A2 + ... An) = 1`` - A = area_meridional_cells_deg(lat_c, 1, dlat, normalize = True, R = 1) + + # Calculate cell areas + # Use dlon = 1 (arbitrary lon spacing and planet + # radius because normalization makes the absolute values + # irrelevant), then normalize + R = 1. + dlon = 1. + + # Compute total surface area for normalization + lon1_deg = -dlon/2 + lon2_deg = dlon/2 + lat1_deg = lat_c[0] - dlat/2 + lat2_deg = lat_c[-1] + dlat/2 + + # Convert to radians for area calculation + lat1_rad = lat1_deg * np.pi/180 + lat2_rad = lat2_deg * np.pi/180 + lon1_rad = lon1_deg * np.pi/180 + lon2_rad = lon2_deg * np.pi/180 + + area_tot = ((R**2) + * np.abs(lon1_rad - lon2_rad) + * np.abs(np.sin(lat1_rad) - np.sin(lat2_rad))) + + # Convert to radians + lat_c_rad = lat_c * np.pi/180 + dlon_rad = dlon * np.pi/180 + dlat_rad = dlat * np.pi/180 + + # Calculate normalized areas. Areas sum to 1 + A = (2. * R**2 * dlon_rad * np.cos(lat_c_rad) * + np.sin(dlat_rad/2.) / area_tot) + # Check if var is 1D (just lat values) if len(var_shape) == 1: - # Variable is a 1D array of size = [lat]. Easiest case - # since ``(w1 + w2 + ...wn) = sum(A) = 1`` and - # ``N = len(lat)`` + # For 1D case: multiply normalized areas by the number of + # lat points. Creates weights where sum(W) = len(lat_c), + # allowing np.mean(var*W) to give the area-weighted average W = A * len(lat_c) else: - # Generate the appropriate shape for the area ``A``, - # e.g., [time, lev, lat, lon] > [1, 1, lat, 1] - # In this case,`` N = time * lev * lat * lon`` and - # ``(w1 + w2 + ...wn) = time * lev * lon * sum(A)``, - # therefore ``N / (w1 + w2 + ...wn) = lat`` + # For multidimensional data: create a list of 1s with + # length matching the number of dims in the var + # (e.g., [1, 1, 1, 1] for 4D data). reshape_shape = [1 for i in range(0, len(var_shape))] + # Set the lat dim to the actual number of lat points + # (e.g., [1, 1, lat, 1] for 4D data where lat = third dim) reshape_shape[axis] = len(lat_c) + # Reshape the 1D area array to match the var's + # dimensionality (broadcasting-ready shape), then multiply + # by the number of lat points. Allows the weights to + # broadcast correctly across all other dims. W = A.reshape(reshape_shape)*len(lat_c) + + # Expand the weights to full var shape by multiplying with an + # array of 1s. Creates the final weight array that can be + # directly multiplied with other data for area-weighted avg. return W*np.ones(var_shape) @@ -1738,7 +1707,8 @@ def get_trend_2D(VAR, LON, LAT, type_trend="wmean"): # dimensions # Flatten array (``[10, 36, lat, lon]`` -> ``[360, lat, lon]``) - nflatten = int(np.prod(var_shape[:-2])) + prod_result = np.prod(var_shape[:-2]) + nflatten = int(prod_result.item()) if hasattr(prod_result, 'item') else int(prod_result) reshape_flat = np.append(nflatten, var_shape[-2:]) VAR = VAR.reshape(reshape_flat) diff --git a/amescap/Ncdf_wrapper.py b/amescap/Ncdf_wrapper.py index e78a52df..9443d4c5 100644 --- a/amescap/Ncdf_wrapper.py +++ b/amescap/Ncdf_wrapper.py @@ -119,20 +119,21 @@ def add_constant(self, variable_name, value, longname_txt="", # Private Definitions def _def_variable(self, variable_name, dim_array, longname_txt="", - units_txt="",datatype="f4"): + units_txt="",datatype="f4", fill_value=None): self.var_dict[variable_name] = self.f_Ncdf.createVariable(variable_name, datatype, - dim_array) + dim_array, + fill_value=fill_value) self.var_dict[variable_name].units = units_txt self.var_dict[variable_name].long_name = longname_txt - self.var_dict[variable_name].dim_name = str(dim_array) def _def_axis1D(self, variable_name, dim_array, longname_txt="", - units_txt="", cart_txt="",datatype="f8"): + units_txt="", cart_txt="",datatype="f8", fill_value=None): self.var_dict[variable_name] = self.f_Ncdf.createVariable(variable_name, datatype, - dim_array) + dim_array, + fill_value=fill_value) self.var_dict[variable_name].units = units_txt self.var_dict[variable_name].long_name = longname_txt self.var_dict[variable_name].cartesian_axis = cart_txt @@ -178,7 +179,6 @@ def log_variable(self, variable_name, DATAin, dim_array, longname_txt="", self._def_variable(variable_name, dim_array, longname_txt, units_txt,datatype) self.var_dict[variable_name].long_name = longname_txt - self.var_dict[variable_name].dim_name = str(dim_array) self.var_dict[variable_name].units = units_txt self.var_dict[variable_name][:] = DATAin @@ -228,6 +228,10 @@ def add_dim_with_content(self, dimension_name, DATAin, longname_txt="", self.var_dict[dimension_name].units = units_txt self.var_dict[dimension_name].cartesian_axis = cart_txt self.var_dict[dimension_name][:] = DATAin + + # Add standard attributes for vertical pressure coordinates + if dimension_name in ['pfull', 'phalf']: + self.var_dict[dimension_name].positive = "down" # .. note:: The attribute ``name`` was replaced by ``_name`` for # compatibility with MFDataset: @@ -245,8 +249,57 @@ def copy_Ncaxis_with_content(self, Ncdim_var): longname_txt = getattr(Ncdim_var, "long_name", Ncdim_var._name) units_txt = getattr(Ncdim_var, "units", "") cart_txt = getattr(Ncdim_var, "cartesian_axis", "") - self.add_dim_with_content(Ncdim_var._name, Ncdim_var[:], longname_txt, - units_txt, cart_txt) + + # Get _FillValue if it exists + fill_value = getattr(Ncdim_var, "_FillValue", None) + + # Add dimension with content + if Ncdim_var._name not in self.dim_dict.keys(): + self.add_dimension(Ncdim_var._name, len(Ncdim_var)) + + if Ncdim_var._name not in self.var_dict.keys(): + # Check if this is actually a 1D variable or multi-dimensional + if len(Ncdim_var.shape) == 1: + self._def_axis1D(Ncdim_var._name, Ncdim_var._name, longname_txt, + units_txt, cart_txt, fill_value=fill_value) + + else: + # For multi-dimensional coordinate variables (like areo + # with shape (time, scalar_axis)) + # we need to create them with the proper dimensions + self.log_variable(Ncdim_var._name, + Ncdim_var[:], # Get the data + Ncdim_var.dimensions, # Use original dimensions + longname_txt, + units_txt) + # Set cartesian_axis if it exists + if cart_txt: + self.var_dict[Ncdim_var._name].cartesian_axis = cart_txt + # Return early since log_variable already copied the data + if Ncdim_var._name in ['pfull', 'phalf', 'plev', 'level']: + self.var_dict[Ncdim_var._name].positive = "down" + for attr_name in Ncdim_var.ncattrs(): + if attr_name not in ['long_name', 'units', 'cartesian_axis', '_FillValue', 'positive']: + setattr(self.var_dict[Ncdim_var._name], attr_name, + getattr(Ncdim_var, attr_name)) + return + + self.var_dict[Ncdim_var._name].long_name = longname_txt + self.var_dict[Ncdim_var._name].units = units_txt + self.var_dict[Ncdim_var._name].cartesian_axis = cart_txt + + # Use chunked copying instead of loading all at once + self._copy_data_chunked(Ncdim_var, self.var_dict[Ncdim_var._name]) + + # Add standard attributes for vertical pressure coordinates + if Ncdim_var._name in ['pfull', 'phalf', 'plev', 'level']: + self.var_dict[Ncdim_var._name].positive = "down" + + # Copy any additional attributes (except _FillValue and the ones we already set) + for attr_name in Ncdim_var.ncattrs(): + if attr_name not in ['long_name', 'units', 'cartesian_axis', '_FillValue', 'positive']: + setattr(self.var_dict[Ncdim_var._name], attr_name, + getattr(Ncdim_var, attr_name)) def copy_Ncvar(self, Ncvar, swap_array=None): @@ -261,19 +314,74 @@ def copy_Ncvar(self, Ncvar, swap_array=None): dim_array = Ncvar.dimensions longname_txt = getattr(Ncvar, "long_name", Ncvar._name) units_txt = getattr(Ncvar, "units", "") + + # Get _FillValue if it exists (must be set at creation time) + fill_value = getattr(Ncvar, "_FillValue", None) + self._def_variable(Ncvar._name, Ncvar.dimensions, longname_txt, - units_txt,Ncvar.dtype) - if np.any(swap_array): - self.log_variable(Ncvar._name, swap_array[:], Ncvar.dimensions, - longname_txt, units_txt) + units_txt, Ncvar.dtype, fill_value=fill_value) + + # Copy ALL OTHER attributes (except _FillValue which was already set) + for attr_name in Ncvar.ncattrs(): + if attr_name != "_FillValue": + setattr(self.var_dict[Ncvar._name], attr_name, + getattr(Ncvar, attr_name)) + + # Set the data - use chunked copying for memory efficiency + if swap_array is not None: + self.var_dict[Ncvar._name][:] = swap_array else: - self.log_variable(Ncvar._name, Ncvar[:], Ncvar.dimensions, - longname_txt, units_txt) + # For large variables, copy in chunks along the first dimension + self._copy_data_chunked(Ncvar, self.var_dict[Ncvar._name]) else: print(f"***Warning***, '{Ncvar._name}' is already defined, " - f"skipping it") - + f"skipping it") + def _copy_data_chunked(self, src_var, dst_var, chunk_size=100): + """ + Copy data from source to destination variable in chunks to avoid + memory issues with large files. + + :param src_var: source netCDF variable + :param dst_var: destination netCDF variable + :param chunk_size: number of records to copy at once along first dimension + """ + shape = src_var.shape + print(f"Copying variable '{src_var._name}' of shape {shape} in chunks...") + + if len(shape) == 0: + # Scalar variable + dst_var[:] = src_var[:] + elif len(shape) == 1: + # 1D variable - copy in chunks + for i in range(0, shape[0], chunk_size): + end_i = min(i + chunk_size, shape[0]) + dst_var[i:end_i] = src_var[i:end_i] + elif len(shape) == 2: + # 2D variable - copy in chunks along first dimension + for i in range(0, shape[0], chunk_size): + end_i = min(i + chunk_size, shape[0]) + dst_var[i:end_i, :] = src_var[i:end_i, :] + elif len(shape) == 3: + # 3D variable - copy in chunks along first dimension + for i in range(0, shape[0], chunk_size): + end_i = min(i + chunk_size, shape[0]) + dst_var[i:end_i, :, :] = src_var[i:end_i, :, :] + elif len(shape) == 4: + # 4D variable - copy in chunks along first dimension + for i in range(0, shape[0], chunk_size): + end_i = min(i + chunk_size, shape[0]) + dst_var[i:end_i, :, :, :] = src_var[i:end_i, :, :, :] + elif len(shape) == 5: + print("5D variable - copy in chunks along first dimension") + for i in range(0, shape[0], chunk_size): + end_i = min(i + chunk_size, shape[0]) + dst_var[i:end_i, :, :, :, :] = src_var[i:end_i, :, :, :, :] + else: + # For higher dimensions (there shouldn't be any), fall back to full copy + print("***Warning***, variable has more than 5 dimensions, copying all at once") + dst_var[:] = src_var[:] + def copy_all_dims_from_Ncfile(self, Ncfile_in, exclude_dim=[], time_unlimited=True): """ diff --git a/amescap/Script_utils.py b/amescap/Script_utils.py index 5bb0e483..cd5247fd 100644 --- a/amescap/Script_utils.py +++ b/amescap/Script_utils.py @@ -27,6 +27,8 @@ import re from matplotlib.colors import ListedColormap, hex2color +from amescap.FV3_utils import area_weights_deg + # ====================================================================== # DEFINITIONS # ====================================================================== @@ -183,23 +185,37 @@ def print_varContent(fileNcdf, list_varfull, print_stat=False): try: slice = "[:]" if "[" in varfull: - varname, slice = varfull.strip().split("[") + cmd_txt, slice = varfull.strip().split("[") slice = f"[{slice}" else: - varname = varfull.strip() + cmd_txt = varfull.strip() - cmd_txt = f"f.variables['{varname}']{slice}" + varname = f"f.variables['{cmd_txt}']{slice}" f = Dataset(fileNcdf.name, "r") - var = eval(cmd_txt) + var = eval(varname) + + # Get the full latitude array (not sliced) + lat = f.variables['lat'][:] + if print_stat: Min = np.nanmin(var) Mean = np.nanmean(var) Max = np.nanmax(var) + + try: + weight = area_weights_deg(var.shape, lat) + Wmean = np.nanmean(var * weight) # print at end + last_wmean = Wmean # Store it + last_varfull = varfull # Store variable name + except: + # For non-spatial variables or if weighting fails, use regular mean + Wmean = Mean + print(f"{Cyan}{varfull:>26s}|{Min:>15g}|{Mean:>15g}|" f"{Max:>15g}|{Nclr}") - if varname == "areo": + if cmd_txt == "areo": # If variable is areo then print modulo print(f"{Cyan}{varfull:>17s}(mod 360)|" f"({(np.nanmin(var%360)):>13g})|" @@ -207,7 +223,7 @@ def print_varContent(fileNcdf, list_varfull, print_stat=False): f"({(np.nanmax(var%360)):>13g})|{Nclr}") else: - if varname != "areo": + if cmd_txt != "areo": print(f"{Cyan}{varfull}= {Nclr}") print(f"{Cyan}{var}{Nclr}") else: @@ -227,39 +243,18 @@ def print_varContent(fileNcdf, list_varfull, print_stat=False): ) else: print(f"{Red}{varfull}") + if print_stat: # Last line for the table - print(f"{Cyan}__________________________|_______________|_______________|_______________|{Nclr}") + print(f"{Cyan}__________________________|_______________|__" + f"_____________|_______________|") + # Only print weighted mean if we successfully calculated one + if last_wmean is not None: + print(f" Global area-weighted mean {last_varfull}: " + f"{last_wmean:.3f}{Nclr}\n") f.close() -def give_permission(filename): - """ - Sets group file permissions for the NAS system - - :param filename: full path to the netCDF file - :type filename: str - - :return: None - - :raises subprocess.CalledProcessError: if the setfacl command fails - :raises FileNotFoundError: if the file is not found - """ - - try: - # catch error and standard output - subprocess.check_call( - ["setfacl -v"], - shell=True, - stdout=open(os.devnull, "w"), - stderr=open(os.devnull, "w") - ) - cmd_txt = f"setfacl -R -m g:s0846:r {filename}" - subprocess.call(cmd_txt, shell=True) - except subprocess.CalledProcessError: - pass - - def check_file_tape(fileNcdf): """ Checks whether a file exists on the disk. diff --git a/amescap/Spectral_utils.py b/amescap/Spectral_utils.py index f609fc7c..2e4419d9 100644 --- a/amescap/Spectral_utils.py +++ b/amescap/Spectral_utils.py @@ -477,11 +477,19 @@ def zonal_construct(COEFFS_flat, VAR_shape, btype=None, low_highcut=None): VAR = np.zeros((nflatten, VAR_shape[-2], VAR_shape[-1])) if btype == "low": - kmax= int(low_highcut) + # Handle numpy arrays properly + kmax_val = low_highcut.item() if hasattr(low_highcut, 'item') else low_highcut + kmax = int(kmax_val) if btype == "high": - kmin= int(low_highcut) + # Handle numpy arrays properly + kmin_val = low_highcut.item() if hasattr(low_highcut, 'item') else low_highcut + kmin = int(kmin_val) if btype == "band": - kmin, kmax= int(low_highcut[0]), int(low_highcut[1]) + # Handle numpy arrays properly + kmin_val = low_highcut[0].item() if hasattr(low_highcut[0], 'item') else low_highcut[0] + kmax_val = low_highcut[1].item() if hasattr(low_highcut[1], 'item') else low_highcut[1] + kmin = int(kmin_val) + kmax = int(kmax_val) for ii in range(0, nflatten): # Filtering diff --git a/bin/MarsFiles.py b/bin/MarsFiles.py index 38e6bd95..1cc73dbc 100755 --- a/bin/MarsFiles.py +++ b/bin/MarsFiles.py @@ -26,7 +26,7 @@ * ``[-bps, --band_pass_spatial]`` Spatial filter: band-pass * ``[-tide, --tide_decomp]`` Extract diurnal tide and its harmonics * ``[-recon, --reconstruct]`` Reconstruct the first N harmonics - * ``[-norm, --normalize]`` Provide ``-tide`` result in % amplitude + * ``[-norm, --normalize]`` Provide ``-tide`` result in percent amplitude * ``[-prop, --prop_tides]`` Extract propagating tide harmonics * ``[-regrid, --regrid_XY_to_match]`` Regrid a target file to match a source file * ``[-zavg, --zonal_average]`` Zonally average all variables in a file @@ -499,16 +499,16 @@ def wrapper(*args, **kwargs): f"Capabilities' in the installation instructions at \n" f"https://amescap.readthedocs.io/en/latest/installation.html" f"{Nclr}\n" - f"Use fourier decomposition to break down the signal into N " - f"harmonics.\nOnly works with 'diurn' files.\nReturns the phase " - f"and amplitude of the variable.\n" + f"Use fourier decomposition to break down a variable into N " + f"harmonics to isolate the diurnal (N=1), semi-diurnal (N=2), etc. phase and " + f"amplitude of a surface or 3D field.\nOnly works with 'diurn' files.\nReturns the phase" + f"and amplitude in separate variables.\n" f"N = 1 diurnal tide, N = 2 semi-diurnal, etc.\n" f"Works on 'diurn' files only.\n" f"{Yellow}Generates a new file ending in ``_tide_decomp.nc``\n" f"{Green}Example:\n" - f"> MarsFiles 01336.atmos_diurn.nc -tide 2 -incl ps temp\n" - f"{Blue}(extracts semi-diurnal tide component of ps and\ntemp " - f"variables; 2 harmonics)" + f"> MarsFiles 01336.atmos_diurn.nc -tide 2 -incl ps\n" + f"{Blue}(extracts semi-diurnal tide component of ps; 2 harmonics)" f"{Nclr}\n\n" ) ) @@ -598,8 +598,8 @@ def wrapper(*args, **kwargs): parser=parser, nargs=0, help=( - f"For use with ``-tide``: Returns the result in percent " - f"amplitude.\n" + f"For use with ``-tide``: Returns the amplitude as a percent daily mean " + f"rather than an absolute value.\n" f"N = 1 diurnal tide, N = 2 semi-diurnal, etc.\n" f"Works on 'diurn' files only.\n" f"{Yellow}Generates a new file ending in ``_norm.nc``\n" @@ -1719,27 +1719,27 @@ def main(): if btype == "low": fnew.add_constant( "sol_max", - nsol, + float(nsol.item()) if hasattr(nsol, 'item') else float(nsol), "Low-pass filter cut-off period ", "sol" ) elif btype == "high": fnew.add_constant( "sol_min", - nsol, + float(nsol.item()) if hasattr(nsol, 'item') else float(nsol), "High-pass filter cut-off period ", "sol" ) elif btype == "band": fnew.add_constant( "sol_min", - nsol[0], + float(nsol[0]), "High-pass filter low cut-off period ", "sol" ) fnew.add_constant( "sol_max", - nsol[1], + float(nsol[1]), "High-pass filter high cut-off period ", "sol" ) @@ -1752,8 +1752,12 @@ def main(): if btype == "band": # Flip the sols so that the low frequency comes first low_highcut = 1/nsol[::-1] + # Ensure it's a proper array/list for butter function + low_highcut = [float(low_highcut[0]), float(low_highcut[1])] else: - low_highcut = 1./nsol + # Extract scalar value for single cutoff + low_highcut_val = nsol.item() if hasattr(nsol, 'item') else nsol + low_highcut = float(1./low_highcut_val) # Loop over all variables in the file for ivar in var_list: @@ -1879,27 +1883,27 @@ def main(): if btype == "low": fnew.add_constant( "kmax", - nk, + float(nk.item()) if hasattr(nk, 'item') else float(nk), "Low-pass filter zonal wavenumber ", "wavenumber" ) elif btype == "high": fnew.add_constant( "kmin", - nk, + float(nk.item()) if hasattr(nk, 'item') else float(nk), "High-pass filter zonal wavenumber ", "wavenumber" ) elif btype == "band": fnew.add_constant( "kmin", - nk[0], + float(nk[0]), "Band-pass filter low zonal wavenumber ", "wavenumber" ) fnew.add_constant( "kmax", - nk[1], + float(nk[1]), "Band-pass filter high zonal wavenumber ", "wavenumber" ) @@ -1970,7 +1974,7 @@ def main(): fdiurn = Dataset(input_file_name, "r", format="NETCDF4_CLASSIC") - var_list = filter_vars(fdiurn, args.include) + var_list = filter_vars(fdiurn, args.include + ["ps"]) # Find time_of_day variable name tod_name = find_tod_in_diurn(fdiurn) @@ -2126,7 +2130,7 @@ def main(): fdiurn = Dataset(input_file_name, "r", format="NETCDF4_CLASSIC") - var_list = filter_vars(fdiurn, args.include) + var_list = filter_vars(fdiurn, args.include + ["ps"]) # Find time_of_day variable name tod_name = find_tod_in_diurn(fdiurn) @@ -2214,10 +2218,11 @@ def main(): "hr" ) - elif ivar in ["pfull", "lat", "lon", "phalf", "pk", + elif ivar in ["pfull", "lat", "lon", "phalf", "pk", "bk", "pstd", "zstd", "zagl", "time"]: print(f"{Cyan}Copying axis: {ivar}...{Nclr}") fnew.copy_Ncaxis_with_content(fdiurn.variables[ivar]) + elif ivar in ["areo"]: print(f"{Cyan}Processing: {ivar}...{Nclr}") # Create areo variable reflecting the diff --git a/bin/MarsPlot.py b/bin/MarsPlot.py index c10e6864..699061cb 100755 --- a/bin/MarsPlot.py +++ b/bin/MarsPlot.py @@ -600,7 +600,6 @@ def main(): with open(output_pdf, "wb") as f: writer.write(f) - give_permission(output_pdf) print(f"{output_pdfq} was generated") except subprocess.CalledProcessError: # If gs command fails again, prompt user to try @@ -1811,33 +1810,9 @@ def make_template(): customFileIN.write(" \n") customFileIN.close() - # NAS system only: set group permissions & print confirmation - give_permission(newname) print(f"{newname} was created") -def give_permission(filename): - """ - Sets group permissions for files created on NAS. - - :param filename: name of the file - :type filename: str - :raises ValueError: If the input filename is not a valid type - for file name. - """ - - # NAS system only: set group permissions to file - try: - # Catch error and standard output - subprocess.check_call(["setfacl -v"], shell = True, - stdout = open(os.devnull, "w"), - stderr = open(os.devnull, "w")) - cmd_txt = f"setfacl -R -m g:s0846:r {filename}" - subprocess.call(cmd_txt, shell = True) - except subprocess.CalledProcessError: - pass - - def namelist_parser(Custom_file): """ Parse a ``Custom.in`` template. diff --git a/docs/source/autoapi/amescap/FV3_utils/index.rst b/docs/source/autoapi/amescap/FV3_utils/index.rst index 2bbe381e..b664cf14 100644 --- a/docs/source/autoapi/amescap/FV3_utils/index.rst +++ b/docs/source/autoapi/amescap/FV3_utils/index.rst @@ -5,10 +5,10 @@ .. autoapi-nested-parse:: - FV3_utils contains internal Functions for processing data in MGCM + FV3_utils contains internal Functions for processing data in MGCM output files such as vertical interpolation. - These functions can be used on their own outside of CAP if they are + These functions can be used on their own outside of CAP if they are imported as a module:: from /u/path/FV3_utils import fms_press_calc @@ -18,7 +18,6 @@ * ``numpy`` * ``warnings`` * ``scipy`` - @@ -36,7 +35,6 @@ Functions amescap.FV3_utils.UT_LTtxt amescap.FV3_utils.add_cyclic amescap.FV3_utils.alt_KM - amescap.FV3_utils.area_meridional_cells_deg amescap.FV3_utils.area_weights_deg amescap.FV3_utils.areo_avg amescap.FV3_utils.axis_interp @@ -72,7 +70,6 @@ Functions amescap.FV3_utils.regression_2D amescap.FV3_utils.robin2cart amescap.FV3_utils.second_hhmmss - amescap.FV3_utils.sfc_area_deg amescap.FV3_utils.shiftgrid_180_to_360 amescap.FV3_utils.shiftgrid_360_to_180 amescap.FV3_utils.sol2ls @@ -95,15 +92,12 @@ Functions the Martian water cycle as inferred from a general circulation model :param ls: solar longitude [°] - :type ls: array - + :type ls: array :param lat : latitude [°] - :type lat: array - + :type lat: array :return: top altitude for the dust [km] - .. py:function:: MGSzmax_ls_lat(ls, lat) Return the max altitude for the dust from "MGS scenario" from @@ -111,50 +105,41 @@ Functions the Martian water cycle as inferred from a general circulation model :param ls: solar longitude [°] - :type ls: array - + :type ls: array :param lat : latitude [°] - :type lat: array - + :type lat: array :return: top altitude for the dust [km] - .. py:function:: UT_LTtxt(UT_sol, lon_180=0.0, roundmin=None) Returns the time in HH:MM:SS at a certain longitude. :param time_sol: the time in sols - :type time_sol: float - + :type time_sol: float :param lon_180: the center longitude in -180/180 coordinates. Increments by 1hr every 15° - :type lon_180: float - + :type lon_180: float :param roundmin: round to the nearest X minute. Typical values are ``roundmin = 1, 15, 60`` - :type roundmin: int + :type roundmin: int - .. note:: + .. note:: If ``roundmin`` is requested, seconds are not shown - .. py:function:: add_cyclic(data, lon) Add a cyclic (overlapping) point to a 2D array. Useful for azimuth and orthographic projections. :param data: variable of size ``[nlat, nlon]`` - :type data: array - + :type data: array :param lon: longitudes - :type lon: array - + :type lon: array :return: a 2D array of size ``[nlat, nlon+1]`` with last column identical to the 1st; and a 1D array of longitudes size [nlon+1] where the last element is ``lon[-1] + dlon`` - .. py:function:: alt_KM(press, scale_height_KM=8.0, reference_press=610.0) @@ -162,71 +147,25 @@ Functions Gives the approximate altitude [km] for a given pressure :param press: the pressure [Pa] - :type press: 1D array - + :type press: 1D array :param scale_height_KM: scale height [km] (default is 8 km, an isothermal at 155K) - :type scale_height_KM: float - + :type scale_height_KM: float :param reference_press: reference surface pressure [Pa] (default is 610 Pa) - :type reference_press: float - + :type reference_press: float :return: ``z_KM`` the equivalent altitude for that pressure [km] - .. note:: + .. note:: Scale height is ``H = rT/g`` - -.. py:function:: area_meridional_cells_deg(lat_c, dlon, dlat, normalize=False, R=3390000.0) - - Return area of invidual cells for a meridional band of thickness - ``dlon`` where ``S = int[R^2 dlon cos(lat) dlat]`` with - ``sin(a)-sin(b) = 2 cos((a+b)/2)sin((a+b)/2)`` - so ``S = 2 R^2 dlon 2cos(lat)sin(dlat/2)``:: - - _________lat + dlat/2 - \ lat \ ^ - \lon + \ | dlat - \________\lat - dlat/2 v - lon - dlon/2 lon + dlon/2 - <------> - dlon - - :param lat_c: latitude of cell center [°] - :type lat_c: float - - :param dlon: cell angular width [°] - :type dlon: float - - :param dlat: cell angular height [°] - :type dlat: float - - :param R: planetary radius [m] - :type R: float - - :param normalize: if True, the sum of the output elements = 1 - :type normalize: bool - - :return: ``S`` areas of the cells, same size as ``lat_c`` in [m2] - or normalized by the total area - - - .. py:function:: area_weights_deg(var_shape, lat_c, axis=-2) - Return weights for averaging the variable. - - :param var_shape: variable shape - :type var_shape: tuple - - :param lat_c: latitude of cell centers [°] - :type lat_c: float - - :param axis: position of the latitude axis for 2D and higher - dimensional arrays. The default is the SECOND TO LAST dimension - :type axis: int + Returns weights scaled so that np.mean(var*W) gives an area-weighted + average. This works because grid cells near the poles have smaller + areas than those at the equator, so they should contribute less to + a global average. Expected dimensions are: @@ -241,13 +180,19 @@ Functions may be truncated on either end (e.g., ``lat = [-20 ..., 0... 50]``) but must be continuous. + :param var_shape: the shape/dimensions of your data array + :type var_shape: tuple + :param lat_c: latitude cell centers in degrees [°] + :type lat_c: float + :param axis: which dimension contains latitude, default: 2nd-to-last + :type axis: int :return: ``W`` weights for the variable ready for standard averaging as ``np.mean(var*W)`` [condensed form] or ``np.average(var, weights=W)`` [expanded form] .. note:: Given a variable var: - + ``var = [v1, v2, ...vn]`` The regular average is @@ -263,40 +208,19 @@ Functions ``W = [w1, w2, ... , wn]*N / (w1 + w2 + ...wn)`` Therfore taking a regular average of (``var*W``) with - ``np.mean(var*W)`` or ``np.average(var, weights=W)`` - + ``np.mean(var*W)`` or ``np.average(var, weights=W)`` + returns the weighted average of the variable. Use ``np.average(var, weights=W, axis = X)`` to average over a specific axis. - .. py:function:: areo_avg(VAR, areo, Ls_target, Ls_angle, symmetric=True) Return a value average over a central solar longitude - :param VAR: a variable with ``time`` in the 1st dimension - :type VAR: ND array - - :param areo: solar longitude of the input variable (0-720) - :type areo: 1D array - - :param Ls_target: central solar longitude of interest - :type Ls_target: float - - :param Ls_angle: requested window angle centered at ``Ls_target`` - :type Ls_angle: float - - :param symmetric: If ``True`` and the requested window is out of range, - ``Ls_angle`` is reduced. If False, the time average is performed - on the data available - :type symmetric: bool (defaults to True) - - :return: the variable averaged over solar longitudes - ``Ls_target-Ls_angle/2`` to ``Ls_target+Ls_angle/2`` - EX:: - + ``Ls_target = 90.`` ``Ls_angle = 10.`` @@ -309,51 +233,57 @@ Functions If ``symmetric = False`` and the input data range is Ls = 88-100° then ``88 < Ls_target < 95`` (7°, assymetric) + :param VAR: a variable with ``time`` in the 1st dimension + :type VAR: ND array + :param areo: solar longitude of the input variable (0-720) + :type areo: 1D array + :param Ls_target: central solar longitude of interest + :type Ls_target: float + :param Ls_angle: requested window angle centered at ``Ls_target`` + :type Ls_angle: float + :param symmetric: If ``True`` and the requested window is out of range, + ``Ls_angle`` is reduced. If False, the time average is performed + on the data available + :type symmetric: bool (defaults to True) + :return: the variable averaged over solar longitudes + ``Ls_target-Ls_angle/2`` to ``Ls_target+Ls_angle/2`` + .. note:: The routine can bin data from muliples Mars years - .. py:function:: axis_interp(var_IN, x, xi, axis, reverse_input=False, type_int='lin', modulo=None) One dimensional linear/logarithmic interpolation along one axis. :param var_IN: Variable on a REGULAR grid (e.g., ``[lev, lat, lon]`` or ``[time, lev, lat, lon]``) - :type var_IN: ND array - + :type var_IN: ND array :param x: Original position array (e.g., ``time``) - :type x: 1D array - + :type x: 1D array :param xi: Target array to interpolate the array on - :type xi: 1D array - + :type xi: 1D array :param axis: Position of the interpolation axis (e.g., ``0`` for a temporal interpolation on ``[time, lev, lat, lon]``) - :type axis: int - + :type axis: int :param reverse_input: Reverse input arrays (e.g., if ``zfull(0)``= 120 km, ``zfull(N)``= 0 km, which is typical) - :type reverse_input: bool - + :type reverse_input: bool :param type_int: "log" for logarithmic (typically pressure), "lin" for linear - :type type_int: str - + :type type_int: str :param modulo: For "lin" interpolation only, use cyclic input (e.g., when using ``modulo = 24`` for time of day, 23.5 and 00 am are considered 30 min apart, not 23.5 hr apart) - :type modulo: float - + :type modulo: float :return: ``VAR_OUT`` interpolated data on the requested axis - .. note:: + .. note:: This routine is similar but simpler than the vertical interpolation ``vinterp()`` as the interpolation axis is assumed to be fully defined by a 1D array such as ``time``, ``pstd`` or ``zstd`` rather than 3D arrays like ``pfull`` or ``zfull``. - For lon/lat interpolation, consider using ``interp_KDTree()``. Calculation:: @@ -361,7 +291,6 @@ Functions X_OUT = Xn*A + (1-A)*Xn + 1 with ``A = log(xi/xn + 1) / log(xn/xn + 1)`` in "log" mode or ``A = (xi-xn + 1)/(xn-xn + 1)`` in "lin" mode - .. py:function:: azimuth2cart(LAT, LON, lat0, lon0=0) @@ -370,36 +299,28 @@ Functions to cartesian coordinates. :param LAT: latitudes[°] size [nlat] - :type LAT: 1D or 2D array - + :type LAT: 1D or 2D array :param LON: longitudes [°] size [nlon] - :type LON: 1D or 2D array - + :type LON: 1D or 2D array :param lat0: latitude coordinate of the pole - :type lat0: float - + :type lat0: float :param lon0: longitude coordinate of the pole - :type lon0: float - + :type lon0: float :return: the cartesian coordinates for the latitudes and longitudes - .. py:function:: broadcast(var_1D, shape_out, axis) Broadcast a 1D array based on a variable's dimensions :param var_1D: variable (e.g., ``lat`` size = 36, or ``time`` size = 133) - :type var_1D: 1D array - + :type var_1D: 1D array :param shape_out: broadcasting shape (e.g., ``temp.shape = [133, lev, 36, lon]``) - :type shape_out: list - + :type shape_out: list :return: (ND array) broadcasted variables (e.g., size = ``[1,36,1,1]`` for ``lat`` or ``[133,1,1,1]`` for ``time``) - .. py:function:: cart_to_azimut_TR(u, v, mode='from') @@ -407,46 +328,36 @@ Functions Convert cartesian coordinates or wind vectors to radians using azimuth angle. :param x: the cartesian coordinate - :type x: 1D array - + :type x: 1D array :param y: the cartesian coordinate - :type y: 1D array - + :type y: 1D array :param mode: "to" for the direction that the vector is pointing, "from" for the direction from which the vector is coming - :type mode: str - + :type mode: str :return: ``Theta`` [°] and ``R`` the polar coordinates - .. py:function:: compute_uneven_sigma(num_levels, N_scale_heights, surf_res, exponent, zero_top) Construct an initial array of sigma based on the number of levels and an exponent :param num_levels: the number of levels - :type num_levels: float - + :type num_levels: float :param N_scale_heights: the number of scale heights to the top of the model (e.g., ``N_scale_heights`` = 12.5 ~102 km assuming an 8 km scale height) - :type N_scale_heights: float - + :type N_scale_heights: float :param surf_res: the resolution at the surface - :type surf_res: float - + :type surf_res: float :param exponent: an exponent to increase the thickness of the levels - :type exponent: float - + :type exponent: float :param zero_top: if True, force the top pressure boundary (in N = 0) to 0 Pa - :type zero_top: bool - + :type zero_top: bool :return: an array of sigma layers - .. py:function:: daily_to_average(varIN, dt_in, nday=5, trim=True) Bin a variable from an ``atmos_daily`` file format to the @@ -454,22 +365,18 @@ Functions :param varIN: variable with ``time`` dimension first (e.g., ``ts[time, lat, lon]``) - :type varIN: ND array - + :type varIN: ND array :param dt_in: delta of time betwen timesteps in sols (e.g., ``dt_in = time[1] - time[0]``) - :type dt_in: float - + :type dt_in: float :param nday: bining period in sols, default is 5 sols - :type nday: int - + :type nday: int :param trim: whether to discard any leftover data at the end of file before binning - :type trim: bool - + :type trim: bool :return: the variable bin over ``nday`` - .. note:: + .. note:: If ``varIN[time, lat, lon]`` from ``atmos_daily`` is ``[160, 48, 96]`` and has 4 timesteps per day (every 6 hours), then the resulting variable for ``nday = 5`` is @@ -481,7 +388,6 @@ Functions time is 133 and last 3 sols the are discarded. If ``trim = False``, the time is 134 and last bin contains only 3 sols of data. - .. py:function:: daily_to_diurn(varIN, time_in) @@ -491,12 +397,10 @@ Functions :param varIN: variable with time dimension first (e.g., ``[time, lat, lon]``) - :type varIN: ND array - + :type varIN: ND array :param time_in: time array in sols. Only the first N elements are actually required if saving memory is important - :type time_in: ND array - + :type time_in: ND array :return: the variable binned in the ``atmos_diurn`` format (``[time, time_of_day, lat, lon]``) and the time of day array [hr] @@ -512,7 +416,6 @@ Functions Since the time dimension is first, the output variables may be passed to the ``daily_to_average()`` function for further binning. - .. py:function:: dvar_dh(arr, h=None) @@ -520,20 +423,6 @@ Functions Differentiate an array ``A[dim1, dim2, dim3...]`` w.r.t ``h``. The differentiated dimension must be the first dimension. - If ``h`` is 1D, then ``h``and ``dim1`` must have the same length - - If ``h`` is 2D, 3D or 4D, then ``arr`` and ``h`` must have the - same shape - - :param arr: variable - :type arr: ND array - - :param h: the dimension (``Z``, ``P``, ``lat``, ``lon``) - :type h: str - - Returns: - d_arr: the array differentiated w.r.t ``h``, e.g., d(array)/dh - EX: Compute ``dT/dz`` where ``T[time, lev, lat, lon]`` is the temperature and ``Zkm`` is the array of level heights [km]. @@ -545,6 +434,17 @@ Functions dTdz = dvar_dh(t.transpose([1, 0, 2, 3]), Zkm).transpose([1, 0, 2, 3]) + If ``h`` is 1D, then ``h``and ``dim1`` must have the same length + + If ``h`` is 2D, 3D or 4D, then ``arr`` and ``h`` must have the + same shape + + :param arr: variable + :type arr: ND array + :param h: the dimension (``Z``, ``P``, ``lat``, ``lon``) + :type h: str + :return: d_arr: the array differentiated w.r.t ``h``, e.g., d(array)/dh + .. py:function:: expand_index(Nindex, VAR_shape_axis_FIRST, axis_list) @@ -552,18 +452,15 @@ Functions :param Nindex: Interpolation indices, size is (``n_axis``, ``Nfull = [time, lat, lon]``) - :type Nindex: idx - + :type Nindex: idx :param VAR_shape_axis_FIRST: Shape for the variable to interpolate with interpolation axis first (e.g., ``[tod, time, lev, lat, lon]``) - :type VAR_shape_axis_FIRST: tuple - + :type VAR_shape_axis_FIRST: tuple :param axis_list: Position or list of positions for axis to insert (e.g., ``2`` for ``lev`` in ``[tod, time, lev, lat, lon]``, ``[2, 4]`` for ``lev`` and ``lon``). The axis positions are those for the final shape (``VAR_shape_axis_FIRST``) and must be INCREASING - :type axis_list: int or list - + :type axis_list: int or list :return: ``LFULL`` a 2D array (size ``n_axis``, ``NfFULL = [time, lev, lat, lon]``) with the indices expanded along the ``lev`` dimension and flattened @@ -576,7 +473,6 @@ Functions This routine expands the indices from ``[tod, time, lat, lon]`` to ``[tod, time, lev, lat, lon]`` with ``Nfull = [time, lev, lat, lon]`` for use in interpolation. - .. py:function:: find_n(X_IN, X_OUT, reverse_input=False, modulo=None) @@ -585,21 +481,17 @@ Functions just below the input values. :param X_IN: Source level [Pa] or [m] - :type X_IN: float or 1D array - + :type X_IN: float or 1D array :param X_OUT: Desired pressure [Pa] or altitude [m] at layer midpoints. Level dimension is FIRST - :type X_OUT: array - + :type X_OUT: array :param reverse_input: If input array is decreasing (e.g., if z(0) = 120 km, z(N) = 0 km, which is typical, or if data is p(0) = 1000 Pa, p(N) = 0 Pa, which is uncommon) - :type reverse_input: bool - + :type reverse_input: bool :return: The index for the level(s) where the pressure < ``plev`` - .. py:function:: find_n0(Lfull_IN, Llev_OUT, reverse_input=False) Return the index for the level(s) just below ``Llev_OUT``. @@ -608,16 +500,13 @@ Functions :param Lfull_IN: Input pressure [Pa] or altitude [m] at layer midpoints. ``Level`` dimension is FIRST - :type Lfull_IN: array - + :type Lfull_IN: array :param Llev_OUT: Desired level type for interpolation [Pa] or [m] - :type Llev_OUT: float or 1D array - + :type Llev_OUT: float or 1D array :param reverse_input: Reverse array (e.g., if ``z(0) = 120 km``, ``z(N) = 0km`` -- which is typical -- or if input data is ``p(0) = 1000Pa``, ``p(N) = 0Pa``) - :type reverse_input: bool - + :type reverse_input: bool :return: ``n`` index for the level(s) where the pressure is just below ``plev`` @@ -629,38 +518,31 @@ Functions If ``Lfull_IN`` is ND ``[lev, time, lat, lon]`` and ``Llev_OUT`` is a 1D array of size ``klev`` then ``n`` is an array of size ``[klev, Ndim]`` with ``Ndim = [time, lat, lon]`` - .. py:function:: fms_Z_calc(psfc, ak, bk, T, topo=0.0, lev_type='full') Returns the 3D altitude field [m] AGL (or above aeroid). - :param psfc: The surface pressure [Pa] or array of surface + :param psfc: The surface pressure [Pa] or array of surface pressures (1D, 2D, or 3D) - :type psfc: array - + :type psfc: array :param ak: 1st vertical coordinate parameter - :type ak: array - + :type ak: array :param bk: 2nd vertical coordinate parameter - :type bk: array - - :param T: The air temperature profile. 1D array (for a single grid + :type bk: array + :param T: The air temperature profile. 1D array (for a single grid point), ND array with VERTICAL AXIS FIRST - :type T: 1D array or ND array - - :param topo: The surface elevation. Same dimension as ``psfc``. + :type T: 1D array or ND array + :param topo: The surface elevation. Same dimension as ``psfc``. If None is provided, AGL is returned - :type topo: array - - :param lev_type: "full" (layer midpoint) or "half" (layer + :type topo: array + :param lev_type: "full" (layer midpoint) or "half" (layer interfaces). Defaults to "full" - :type lev_type: str - - :return: The layer altitude at the full level ``Z_f(:, :, Nk-1)`` - or half-level ``Z_h(:, :, Nk)`` [m]. ``Z_f`` and ``Z_h`` are - AGL if ``topo = None``. ``Z_f`` and ``Z_h`` are above aeroid + :type lev_type: str + :return: The layer altitude at the full level ``Z_f(:, :, Nk-1)`` + or half-level ``Z_h(:, :, Nk)`` [m]. ``Z_f`` and ``Z_h`` are + AGL if ``topo = None``. ``Z_f`` and ``Z_h`` are above aeroid if topography is not None. Calculation:: @@ -674,14 +556,14 @@ Functions --- Nk --- SFC ======== z_half / / / / / - .. note:: + .. note:: Expands to the time dimension using:: topo = np.repeat(zsurf[np.newaxis, :], ps.shape[0], axis = 0) - Calculation is derived from + Calculation is derived from ``./atmos_cubed_sphere_mars/Mars_phys.F90``:: - + # (dp/dz = -rho g) => (dz = dp/(-rho g)) and # (rho = p/(r T)) => (dz = rT/g * (-dp/p)) @@ -697,12 +579,12 @@ Functions Z_half calculation:: - # Hydrostatic relation within the layer > (P(k+1)/P(k) = + # Hydrostatic relation within the layer > (P(k+1)/P(k) = # exp(-DZ(k)/H)) - + # layer thickness: DZ(k) = rT/g * -(du) - + # previous layer altitude + thickness of layer: Z_h k) = Z_h(k+1) +DZ_h(h) @@ -711,7 +593,7 @@ Functions # previous altitude + half the thickness of previous layer and # half of current layer Z_f(k) = Z_f(k+1) + (0.5 DZ(k) + 0.5 DZ(k+1)) - + # Add ``+0.5 DZ(k)-0.5 DZ(k)=0`` and re-organiz the equation Z_f(k) = Z_f(k+1) + DZ(k) + 0.5 (DZ(k+1) - DZ(k)) Z_f(k) = Z_h(k+1) + 0.5 (DZ(k+1) - DZ(k)) @@ -737,15 +619,14 @@ Functions Line 7: ``Z_full[k] = Z_half[k] + (T_half[k+1]-T_full[k])/Γ`` Pulling out ``Tfull`` from above equation and using ``Γ = (gγ)/R``:: - + Z_full[k] = (Z_half[k+1] + (R Tfull[k]) / (gγ)(T_half[k+1] / T_full[k] - 1)) Using the isentropic relation above:: - + Z_full = (Z_half[k+1] + (R Tfull[k]) / (gγ)(p_half[k+1] / p_full[k])**(R/Cp)-1)) - .. py:function:: fms_press_calc(psfc, ak, bk, lev_type='full') @@ -753,21 +634,17 @@ Functions Returns the 3D pressure field from the surface pressure and the ak/bk coefficients. - :param psfc: the surface pressure [Pa] or an array of surface + :param psfc: the surface pressure [Pa] or an array of surface pressures (1D, 2D, or 3D if time dimension) - :type psfc: array - + :type psfc: array :param ak: 1st vertical coordinate parameter - :type ak: array - + :type ak: array :param bk: 2nd vertical coordinate parameter - :type bk: array: - - :param lev_type: "full" (layer midpoints) or "half" + :type bk: array: + :param lev_type: "full" (layer midpoints) or "half" (layer interfaces). Defaults to "full." - :type lev_type: str - - :return: the 3D pressure field at the full levels + :type lev_type: str + :return: the 3D pressure field at the full levels ``PRESS_f(Nk-1:,:,:)`` or half-levels ``PRESS_h(Nk,:,:,)`` [Pa] Calculation:: @@ -784,7 +661,6 @@ Functions .. note:: Some literature uses pk (pressure) instead of ak with ``p3d = ps * bk + P_ref * ak`` instead of ``p3d = ps * bk + ak`` - .. py:function:: frontogenesis(U, V, theta, lon_deg, lat_deg, R=3400 * 1000.0, spacing='varying') @@ -799,43 +675,34 @@ Functions :param U: wind field with ``lat`` SECOND TO LAST and ``lon`` as last dimensions (e.g., ``[lat, lon]`` or ``[time, lev, lat, lon``] etc.) - :type U: array - + :type U: array :param V: wind field with ``lat`` SECOND TO LAST and ``lon`` as last dimensions (e.g., ``[lat, lon]`` or ``[time, lev, lat, lon``] etc.) - :type V: array - + :type V: array :param theta: potential temperature [K] - :type theta: array - + :type theta: array :param lon_deg: longitude [°] (2D if irregularly-spaced) - :type lon_deg: 1D array - + :type lon_deg: 1D array :param lat_deg: latitude [°] (2D if irregularly-spaced) - :type lat_deg: 1D array - + :type lat_deg: 1D array :param R: planetary radius [m] - :type R: float - + :type R: float :param spacing: when ``lon`` and ``lat`` are 1D arrays, using spacing = "varying" differentiates latitude and longitude. When spacing = "regular", ``dx = lon[1]-lon[0]``, `` dy=lat[1]-lat[0]`` and the ``numpy.gradient()`` method are used - :type spacing: str (defaults to "varying") - + :type spacing: str (defaults to "varying") :return: the frontogenesis field [m-1] - .. py:function:: gauss_profile(x, alpha, x0=0.0) Return Gaussian line shape at x. This can be used to generate a bell-shaped mountain. - .. py:function:: get_trend_2D(VAR, LON, LAT, type_trend='wmean') Extract spatial trends from the data. The output can be directly @@ -844,25 +711,20 @@ Functions :param VAR: Variable for decomposition. ``lat`` is SECOND TO LAST and ``lon`` is LAST (e.g., ``[time, lat, lon]`` or ``[time, lev, lat, lon]``) - :type VAR: ND array - + :type VAR: ND array :param LON: longitude coordinates - :type LON: 2D array - + :type LON: 2D array :param LAT: latitude coordinates - :type LAT: 2D array - + :type LAT: 2D array :param type_trend: type of averaging to perform: "mean" - use a constant average over all lat/lon "wmean" - use a area-weighted average over all lat/lon "zonal" - detrend over the zonal axis only "2D" - use a 2D planar regression (not area-weighted) - :type type_trend: str - + :type type_trend: str :return: the trend, same size as ``VAR`` - .. py:function:: interp_KDTree(var_IN, lat_IN, lon_IN, lat_OUT, lon_OUT, N_nearest=10) Inverse distance-weighted interpolation using nearest neighboor for @@ -871,33 +733,26 @@ Functions :param var_IN: ND variable to regrid (e.g., ``[lev, lat, lon]``, ``[time, lev, lat, lon]`` with ``[lat, lon]`` dimensions LAST [°]) - :type var_IN: ND array - + :type var_IN: ND array :param lat_IN: latitude [°] (``LAT[y, x]`` array for irregular grids) - :type lat_IN: 1D or 2D array - + :type lat_IN: 1D or 2D array :param lon_IN: latitude [°] (``LAT[y, x]`` array for irregular grids) - :type lon_IN: 1D or 2D array - + :type lon_IN: 1D or 2D array :param lat_OUT: latitude [°] for the TARGET grid structure (or ``LAT1[y,x]`` for irregular grids) - :type lat_OUT: 1D or 2D array - + :type lat_OUT: 1D or 2D array :param lon_OUT: longitude [°] for the TARGET grid structure (or ``LON1[y,x]`` for irregular grids) - :type lon_OUT: 1D or 2D array - + :type lon_OUT: 1D or 2D array :param N_nearest: number of nearest neighbours for the search - :type N_nearest: int - + :type N_nearest: int :return: ``VAR_OUT`` interpolated data on the target grid .. note:: This implementation is much FASTER than ``griddata`` and it supports unstructured grids like an MGCM tile. - The nearest neighbour interpolation is only done on the lon/lat axis (not level). Although this interpolation works well on the 3D field [x, y, z], this is typically not what is expected. In @@ -905,7 +760,6 @@ Functions on the target grid are 100's of km away while the closest points in the vertical are a few 10's -100's meter in the PBL. This would result in excessive weighting in the vertical. - .. py:function:: layers_mid_point_to_boundary(pfull, sfc_val) @@ -914,19 +768,17 @@ Functions p_half = ps*bk + pk - This routine converts the coordinate of the layer MIDPOINTS, - ``p_full`` or ``bk``, into the coordinate of the layer BOUNDARIES + This routine converts the coordinate of the layer MIDPOINTS, + ``p_full`` or ``bk``, into the coordinate of the layer BOUNDARIES ``p_half``. The surface value must be provided. :param p_full: Pressure/sigma values for the layer MIDPOINTS, INCREASING with ``N`` (e.g., [0.01 -> 720] or [0.001 -> 1]) - :type p_full: 1D array - + :type p_full: 1D array :param sfc_val: The surface value for the lowest layer's boundary ``p_half[N]`` (e.g., ``sfc_val`` = 720 Pa or ``sfc_val`` = 1 in sigma coordinates) - :type sfc_val: float - + :type sfc_val: float :return: ``p_half`` the pressure at the layer boundaries (size = ``N+1``) @@ -948,7 +800,7 @@ Functions = phalf[N] - pfull[N] log(phalf[N]) We want to solve for ``phalf[N-1] = X``:: - + v v v X - pfull[N] log(X) = B @@ -964,52 +816,50 @@ Functions ``100*(phalf - phalf_reconstruct)/phalf < 0.4%`` at the top. - .. py:function:: lin_interp(X_in, X_ref, Y_ref) Simple linear interpolation with no dependance on scipy :param X_in: input values - :type X_in: float or array - + :type X_in: float or array :param X_ref x values - :type X_ref: array - + :type X_ref: array :param Y_ref y values - :type Y_ref: array - + :type Y_ref: array :return: y value linearly interpolated at ``X_in`` - .. py:function:: lon180_to_360(lon) Transform a float or an array from the -180/180 coordinate system to 0-360 :param lon: longitudes in the -180/180 coordinate system - :type lon: float, 1D array, or 2D array - + :type lon: float, 1D array, or 2D array :return: the equivalent longitudes in the 0-360 coordinate system - .. py:function:: lon360_to_180(lon) + Transform a float or an array from the 0-360 coordinate system to + -180/180. + + :param lon: longitudes in the 0-360 coordinate system + :type lon: float, 1D array, or 2D array + :return: the equivalent longitudes in the -180/180 coordinate system + .. py:function:: ls2sol(Ls_in) Ls to sol converter. :param Ls_in: solar longitudes (0-360...720) - :type Ls_in: float or 1D array - + :type Ls_in: float or 1D array :return: the corresponding sol number .. note:: This function simply uses a numerical solver on the ``sol2ls()`` function. - .. py:function:: mass_stream(v_avg, lat, level, type='pstd', psfc=700, H=8000.0, factor=1e-08) @@ -1025,35 +875,28 @@ Functions :param v_avg: zonal wind [m/s] with ``lev`` dimension FIRST and ``lat`` dimension SECOND (e.g., ``[pstd, lat]``, ``[pstd, lat, lon]`` or ``[pstd, lat, lon, time]``) - :type v_avg: ND array - + :type v_avg: ND array :param lat: latitudes [°] - :type lat: 1D array - + :type lat: 1D array :param level: interpolated layers [Pa] or [m] - :type level: 1D array - + :type level: 1D array :param type: interpolation type (``pstd``, ``zstd`` or ``zagl``) - :type type: str - + :type type: str :param psfc: reference surface pressure [Pa] - :type psfc: float - + :type psfc: float :param H: reference scale height [m] when pressures are used - :type H: float - + :type H: float :param factor: normalize the mass stream function by a factor, use ``factor = 1`` for [kg/s] - :type factor: int - + :type factor: int :return: ``MSF`` the meridional mass stream function (in ``factor * [kg/s]``) - .. note:: + .. note:: This routine allows the time and zonal averages to be computed before OR after the MSF calculation. - .. note:: + .. note:: The expressions for MSF use log(pressure) Z coordinates, which integrate better numerically. @@ -1073,7 +916,6 @@ Functions ⌠ .g. ⌡ f(z)dz = (Zn-Zn-1){f(Zn) + f(Zn-1)}/2 n-1 - .. py:function:: mollweide2cart(LAT, LON) @@ -1082,41 +924,31 @@ Functions cartesian coordinates. :param LAT: latitudes[°] size [nlat] - :type LAT: 1D or 2D array - + :type LAT: 1D or 2D array :param LON: longitudes [°] size [nlon] - :type LON: 1D or 2D array - + :type LON: 1D or 2D array :param lat0: latitude coordinate of the pole - :type lat0: float - + :type lat0: float :param lon0: longitude coordinate of the pole - :type lon0: float - + :type lon0: float :return: the cartesian coordinates for the latitudes and longitudes - .. py:function:: ortho2cart(LAT, LON, lat0, lon0=0) Orthographic projection. Converts from latitude-longitude to cartesian coordinates. :param LAT: latitudes[°] size [nlat] - :type LAT: 1D or 2D array - + :type LAT: 1D or 2D array :param LON: longitudes [°] size [nlon] - :type LON: 1D or 2D array - + :type LON: 1D or 2D array :param lat0: latitude coordinate of the pole - :type lat0: float - + :type lat0: float :param lon0: longitude coordinate of the pole - :type lon0: float - + :type lon0: float :return: the cartesian coordinates for the latitudes and longitudes; and a mask (NaN array) that hides the back side of the planet - .. py:function:: polar2XYZ(lon, lat, alt, Re=3400 * 10**3) @@ -1124,23 +956,18 @@ Functions Spherical to cartesian coordinate transformation. :param lon: Longitude in radians - :type lon: ND array - + :type lon: ND array :param lat: Latitude in radians - :type lat: ND array - + :type lat: ND array :param alt: Altitude [m] - :type alt: ND array - + :type alt: ND array :param Re: Planetary radius [m], defaults to 3400*10^3 - :type Re: float - + :type Re: float :return: ``X``, ``Y``, ``Z`` in cartesian coordinates [m] .. note:: This is a classic polar coordinate system with ``colatitude = pi/2 - lat`` where ``cos(colat) = sin(lat)`` - .. py:function:: polar_warming(T, lat, outside_range=np.nan) @@ -1151,22 +978,18 @@ Functions :param T: temperature with the lat dimension FIRST (transpose as needed) - :type T: ND array - + :type T: ND array :param lat: latitude array - :type lat: 1D array - + :type lat: 1D array :param outside_range: values to set the polar warming to when outside pf the range. Default = NaN but 0 may be desirable. - :type outside_range: float - + :type outside_range: float :return: The polar warming [K] .. note:: ``polar_warming()`` concatenates the results from both hemispheres obtained from the nested function ``PW_half_hemisphere()`` - .. py:function:: press_pa(alt_KM, scale_height_KM=8.0, reference_press=610.0) @@ -1174,43 +997,36 @@ Functions Gives the approximate altitude [km] for a given pressure :param alt_KM: the altitude [km] - :type alt_KM: 1D array - + :type alt_KM: 1D array :param scale_height_KM: scale height [km] (default is 8 km, an isothermal at 155K) - :type scale_height_KM: float - + :type scale_height_KM: float :param reference_press: reference surface pressure [Pa] (default is 610 Pa) - :type reference_press: float - + :type reference_press: float :return: ``press_pa`` the equivalent pressure at that altitude [Pa] - .. note:: + .. note:: Scale height is ``H = rT/g`` - .. py:function:: press_to_alt_atmosphere_Mars(Pi) Return the altitude [m] as a function of pressure from the analytical calculation above. :param Pi: input pressure [Pa] (<= 610 Pa) - :type Pi: float or 1D array - + :type Pi: float or 1D array :return: the corresponding altitude [m] (float or 1D array) - .. py:function:: ref_atmosphere_Mars_PTD(Zi) Analytical atmospheric model for Martian pressure, temperature, and density. Alex Kling, June 2021 :param Zi: input altitude [m] (must be >= 0) - :type Zi: float or 1D array - + :type Zi: float or 1D array :return: tuple of corresponding pressure [Pa], temperature [K], and density [kg/m3] floats or arrays @@ -1233,7 +1049,6 @@ Functions altitudes, we provide a fit in the form of ``P = P0 exp(-az-bz^2)`` based on diurnal average of the MCD database at lat = 0, Ls = 150. - .. py:function:: regression_2D(X, Y, VAR, order=1) @@ -1241,18 +1056,15 @@ Functions Linear and quadratic regression on the plane. :param X: first coordinate - :type X: 2D array - + :type X: 2D array :param Y: second coordinate - :type Y: 2D array - + :type Y: 2D array :param VAR: variable of the same size as X - :type VAR: 2D array - + :type VAR: 2D array :param order: 1 (linear) or 2 (quadratic) - :type order: int + :type order: int - .. note:: + .. note:: When ``order = 1``, the equation is: ``aX + bY + C = Z``. When ``order = 2``, the equation is: ``aX^2 + 2bX*Y + cY^2 + 2dX + 2eY + f = Z`` @@ -1271,100 +1083,58 @@ Functions minimizes ``||b – A x||^2`` - .. py:function:: robin2cart(LAT, LON) Robinson projection. Converts from latitude-longitude to cartesian coordinates. :param LAT: latitudes[°] size [nlat] - :type LAT: 1D or 2D array - + :type LAT: 1D or 2D array :param LON: longitudes [°] size [nlon] - :type LON: 1D or 2D array - + :type LON: 1D or 2D array :param lat0: latitude coordinate of the pole - :type lat0: float - + :type lat0: float :param lon0: longitude coordinate of the pole - :type lon0: float - + :type lon0: float :return: the cartesian coordinates for the latitudes and longitudes - .. py:function:: second_hhmmss(seconds, lon_180=0.0) Given the time [sec], return local true solar time at a certain longitude. :param seconds: the time [sec] - :type seconds: float - + :type seconds: float :param lon_180: the longitude in -180/180 coordinate - :type lon_180: float - + :type lon_180: float :return: the local time [float] or a tuple (hours, minutes, seconds) - -.. py:function:: sfc_area_deg(lon1, lon2, lat1, lat2, R=3390000.0) - - Return the surface between two sets of latitudes/longitudes:: - - S = int[R^2 dlon cos(lat) dlat] _____lat2 - \ \____\lat1 - lon1 lon2 - :param lon1: longitude from set 1 [°] - :type lon1: float - - :param lon2: longitude from set 2 [°] - :type lon2: float - - :param lat1: latitude from set 1 [°] - :type lat1: float - - :param lat2: longitude from set 2 [°] - :type lat2: float - - :param R: planetary radius [m] - :type R: int - - .. note:: - qLon and Lat define the corners of the area not the grid cell center. - - - .. py:function:: shiftgrid_180_to_360(lon, data) This function shifts ND data from a -180/180 to a 0-360 grid. :param lon: longitudes in the 0-360 coordinate system - :type lon: 1D array - + :type lon: 1D array :param data: variable with ``lon`` in the last dimension - :type data: ND array - + :type data: ND array :return: shifted data - .. py:function:: shiftgrid_360_to_180(lon, data) This function shifts ND data from a 0-360 to a -180/180 grid. :param lon: longitudes in the 0-360 coordinate system - :type lon: 1D array - + :type lon: 1D array :param data: variable with ``lon`` in the last dimension - :type data: ND array - + :type data: ND array :return: shifted data - .. note:: + .. note:: Use ``np.ma.hstack`` instead of ``np.hstack`` to keep the masked array properties - .. py:function:: sol2ls(jld, cumulative=False) @@ -1373,31 +1143,25 @@ Functions Sol=0 is the spring equinox. :param jld: sol number after perihelion - :type jld: float or 1D array - + :type jld: float or 1D array :param cumulative: if True, result is cumulative (Ls=0-360, 360-720 etc..) - :type cumulative: bool - + :type cumulative: bool :return: the corresponding solar longitude - .. py:function:: sol_hhmmss(time_sol, lon_180=0.0) Given the time in days, return return local true solar time at a certain longitude. :param time_sol: the time in sols - :type seconds: float - + :type seconds: float :param lon_180: the longitude in -180/180 coordinate - :type lon_180: float - + :type lon_180: float :return: the local time [float] or a tuple (hours, minutes, seconds) - .. py:function:: spherical_curl(U, V, lon_deg, lat_deg, R=3400 * 1000.0, spacing='varying') Compute the vertical component of the relative vorticity using @@ -1409,33 +1173,26 @@ Functions :param U: wind field with ``lat`` SECOND TO LAST and ``lon`` as last dimensions (e.g., ``[lat, lon]`` or ``[time, lev, lat, lon``] etc.) - :type U: array - + :type U: array :param V: wind field with ``lat`` SECOND TO LAST and ``lon`` as last dimensions (e.g., ``[lat, lon]`` or ``[time, lev, lat, lon``] etc.) - :type V: array - + :type V: array :param lon_deg: longitude [°] (2D if irregularly-spaced) - :type lon_deg: 1D array - + :type lon_deg: 1D array :param lat_deg: latitude [°] (2D if irregularly-spaced) - :type lat_deg: 1D array - + :type lat_deg: 1D array :param R: planetary radius [m] - :type R: float - + :type R: float :param spacing: when ``lon`` and ``lat`` are 1D arrays, using spacing = "varying" differentiates latitude and longitude. When spacing = "regular", ``dx = lon[1]-lon[0]``, `` dy=lat[1]-lat[0]`` and the ``numpy.gradient()`` method are used - :type spacing: str (defaults to "varying") - + :type spacing: str (defaults to "varying") :return: the vorticity of the wind field [m-1] - .. py:function:: spherical_div(U, V, lon_deg, lat_deg, R=3400 * 1000.0, spacing='varying') Compute the divergence of the wind fields using finite difference:: @@ -1446,137 +1203,114 @@ Functions :param U: wind field with ``lat`` SECOND TO LAST and ``lon`` as last dimensions (e.g., ``[lat, lon]`` or ``[time, lev, lat, lon``] etc.) - :type U: array - + :type U: array :param V: wind field with ``lat`` SECOND TO LAST and ``lon`` as last dimensions (e.g., ``[lat, lon]`` or ``[time, lev, lat, lon``] etc.) - :type V: array - + :type V: array :param lon_deg: longitude [°] (2D if irregularly-spaced) - :type lon_deg: 1D array - + :type lon_deg: 1D array :param lat_deg: latitude [°] (2D if irregularly-spaced) - :type lat_deg: 1D array - + :type lat_deg: 1D array :param R: planetary radius [m] - :type R: float - + :type R: float :param spacing: when ``lon`` and ``lat`` are 1D arrays, using spacing = "varying" differentiates latitude and longitude. When spacing = "regular", ``dx = lon[1]-lon[0]``, `` dy=lat[1]-lat[0]`` and the ``numpy.gradient()`` method are used - :type spacing: str (defaults to "varying") - + :type spacing: str (defaults to "varying") :return: the horizonal divergence of the wind field [m-1] - .. py:function:: swinbank(plev, psfc, ptrans=1.0) Compute ``ak`` and ``bk`` values with a transition based on Swinbank :param plev: the pressure levels [Pa] - :type plev: 1D array - + :type plev: 1D array :param psfc: the surface pressure [Pa] - :type psfc: 1D array - + :type psfc: 1D array :param ptrans: the transition pressure [Pa] - :type ptrans: 1D array - + :type ptrans: 1D array :return: the coefficients for the new layers - -.. py:function:: time_shift_calc(array, lon, timeo, timex=None) +.. py:function:: time_shift_calc(var_in, lon, tod, target_times=None) Conversion to uniform local time. - :param array: variable to be shifted. Assume ``lon`` is the first - dimension and ``time_of_day`` is the last dimension - :type array: ND array + Mars rotates approx. 14.6° lon per Mars-hour (360° ÷ 24.6 hr) + Each 14.6° shift in lon represents a 1-hour shift in local time + This code uses the more precise calculation: lon_shift = 24.0 * lon / 360. + :param var_in: variable to be shifted. Assume ``lon`` is the first + dimension and ``tod`` is the last dimension + :type var_in: ND array :param lon: longitude - :type lon: 1D array - - :param timeo: ``time_of_day`` index from the input file - :type timeo: 1D array - - :param timex: local time(s) [hr] to shift to (e.g., ``"3. 15."``) - :type timex: float (optional) - + :type lon: 1D array + :param tod: ``time_of_day`` index from the input file + :type tod: 1D array + :param target_times: local time(s) [hr] to shift to (e.g., ``"3. 15."``) + :type target_times: float (optional) :return: the array shifted to uniform local time .. note:: - If ``timex`` is not specified, the file is interpolated - on the same ``time_of_day`` as the input - + If ``target_times`` is not specified, the file is interpolated + on the same ``tod`` as the input .. py:function:: transition(pfull, p_sigma=0.1, p_press=0.05) Return the transition factor to construct ``ak`` and ``bk`` - :param pfull: the pressure [Pa] - :type pfull: 1D array + In the MGCM code, the full pressures are computed from:: + del(phalf) + pfull = ----------------------------- + log(phalf(k+1/2)/phalf(k-1/2)) + + :param pfull: the pressure [Pa] + :type pfull: 1D array :param p_sigma: the pressure level where the vertical grid starts transitioning from sigma to pressure - :type p_sigma: float - + :type p_sigma: float :param p_press: the pressure level above which the vertical grid is pure (constant) pressure - :type p_press: float - + :type p_press: float :return: the transition factor. = 1 for pure sigma, = 0 for pure pressure and =0-1 for the transition - In the MGCM code, the full pressures are computed from:: - - del(phalf) - pfull = ----------------------------- - log(phalf(k+1/2)/phalf(k-1/2)) - - .. py:function:: vinterp(varIN, Lfull, Llev, type_int='log', reverse_input=False, masktop=True, index=None) Vertical linear or logarithmic interpolation for pressure or altitude. :param varIN: Variable to interpolate (VERTICAL AXIS FIRST) - :type varIN: ND array - + :type varIN: ND array :param Lfull: Pressure [Pa] or altitude [m] at full layers, same dimensions as ``varIN`` - :type Lfull: array - + :type Lfull: array :param Llev: Desired level for interpolation [Pa] or [m]. May be increasing or decreasing as the output levels are processed one at a time - :type Llev: 1D array - + :type Llev: 1D array :param type_int: "log" for logarithmic (typically pressure), "lin" for linear (typically altitude) - :type type_int: str - + :type type_int: str :param reverse_input: Reverse input arrays. e.g., if ``zfull[0]`` = 120 km then ``zfull[N]`` = 0km (typical) or if input data is ``pfull[0]``=1000 Pa, ``pfull[N]``=0 Pa - :type reverse_input: bool - + :type reverse_input: bool :param masktop: Set to NaN values if above the model top - :type masktop: bool - + :type masktop: bool :param index: Indices for the interpolation, already processed as ``[klev, Ndim]``. Indices calculated if not provided - :type index: None or array - + :type index: None or array :return: ``varOUT`` variable interpolated on the ``Llev`` pressure or altitude levels - .. note:: + .. note:: This interpolation assumes pressure decreases with height:: -- 0 -- TOP [0 Pa] : [120 km]| X_OUT = Xn*A + (1-A)*Xn + 1 @@ -1592,7 +1326,6 @@ Functions with ``A = log(Llev/pn + 1) / log(pn/pn + 1)`` in "log" mode or ``A = (zlev-zn + 1) / (zn-zn + 1)`` in "lin" mode - .. py:function:: vw_from_MSF(msf, lat, lev, ztype='pstd', norm=True, psfc=700, H=8000.0) @@ -1600,41 +1333,33 @@ Functions Return the V and W components of the circulation from the mass stream function. - :param msf: the mass stream function with ``lev`` SECOND TO - LAST and the ``lat`` dimension LAST (e.g., ``[lev, lat]``, - ``[time, lev, lat]``, ``[time, lon, lev, lat]``) - :type msf: ND array - - :param lat: latitude [°] - :type lat: 1D array - - :param lev: level [Pa] or [m] (``pstd``, ``zagl``, ``zstd``) - :type lev: 1D array - - :param ztype: Use ``pstd`` for pressure so vertical - differentation is done in log space. - :type ztype: str - - :param norm: if True, normalize ``lat`` and ``lev`` before - differentiating to avoid having to rescale manually the - vectors in quiver plots - :type norm: bool - - :param psfc: surface pressure for pseudo-height when - ``ztype = pstd`` - :type psfc: float - - :param H: scale height for pseudo-height when ``ztype = pstd`` - :type H: float - + :param msf: the mass stream function with ``lev`` SECOND TO + LAST and the ``lat`` dimension LAST (e.g., ``[lev, lat]``, + ``[time, lev, lat]``, ``[time, lon, lev, lat]``) + :type msf: ND array + :param lat: latitude [°] + :type lat: 1D array + :param lev: level [Pa] or [m] (``pstd``, ``zagl``, ``zstd``) + :type lev: 1D array + :param ztype: Use ``pstd`` for pressure so vertical + differentation is done in log space. + :type ztype: str + :param norm: if True, normalize ``lat`` and ``lev`` before + differentiating to avoid having to rescale manually the + vectors in quiver plots + :type norm: bool + :param psfc: surface pressure for pseudo-height when + ``ztype = pstd`` + :type psfc: float + :param H: scale height for pseudo-height when ``ztype = pstd`` + :type H: float :return: the meditional and altitude components of the mass stream function for plotting as a quiver or streamlines. - .. note:: + .. note:: The components are: ``[v]= g/(2 pi cos(lat)) dphi/dz`` ``[w]= -g/(2 pi cos(lat)) dphi/dlat`` - .. py:function:: zonal_detrend(VAR) @@ -1642,15 +1367,13 @@ Functions Substract the zonal average mean value from a field. :param VAR: variable with detrending dimension last - :type VAR: ND array - + :type VAR: ND array :return: detrented field (same size as input) - .. note:: + .. note:: ``RuntimeWarnings`` are expected if the slice contains only NaNs which is the case below the surface and above the model top in the interpolated files. This routine disables such warnings temporarily. - diff --git a/docs/source/autoapi/amescap/Ncdf_wrapper/index.rst b/docs/source/autoapi/amescap/Ncdf_wrapper/index.rst index 8f8aeef1..1078a321 100644 --- a/docs/source/autoapi/amescap/Ncdf_wrapper/index.rst +++ b/docs/source/autoapi/amescap/Ncdf_wrapper/index.rst @@ -16,7 +16,6 @@ * ``netCDF4`` * ``os`` * ``datetime`` - @@ -62,11 +61,9 @@ Classes :param object: _description_ :type object: _type_ - :return: _description_ :rtype: _type_ - .. py:class:: Fort_var(input_vals, name_txt, long_name_txt, units_txt, dimensions_tuple) @@ -87,11 +84,9 @@ Classes :param np.ndarray: _description_ :type np.ndarray: _type_ - :return: _description_ :rtype: _type_ - .. py:method:: __abs__() @@ -433,27 +428,23 @@ Classes Create average file (e.g., N-day averages [N=5 usually]) - .. py:method:: write_to_daily() Create daily file (continuous time series) - .. py:method:: write_to_diurn(day_average=5) Create diurn file (variables organized by time of day & binned (typically 5-day bins) - .. py:method:: write_to_fixed() Create ``fixed`` file (all static variables) - .. py:class:: Ncdf(filename=None, description_txt='', action='w', ncformat='NETCDF4_CLASSIC') @@ -487,10 +478,8 @@ Classes :param object: _description_ :type object: _type_ - :return: netCDF file - .. py:method:: add_constant(variable_name, value, longname_txt='', units_txt='') @@ -506,7 +495,6 @@ Classes Log.add_dim_with_content("lon", lon_array, "longitudes", "degree", "X") - .. py:method:: add_dimension(dimension_name, length) @@ -522,7 +510,6 @@ Classes yet, it will be created - .. py:method:: copy_Ncvar(Ncvar, swap_array=None) Copy a netCDF variable from another file (e.g., @@ -531,14 +518,12 @@ Classes swapped with this array. - .. py:method:: copy_all_dims_from_Ncfile(Ncfile_in, exclude_dim=[], time_unlimited=True) Copy all variables, dimensions, and attributes from another netCDF file - .. py:method:: copy_all_vars_from_Ncfile(Ncfile_in, exclude_var=[]) @@ -547,16 +532,14 @@ Classes EX:: Log.log_axis1D("areo", areo, "time", "degree", "T") - - .. py:method:: log_variable(variable_name, DATAin, dim_array, longname_txt='', units_txt='') + .. py:method:: log_variable(variable_name, DATAin, dim_array, longname_txt='', units_txt='', datatype='f4') EX:: Log.log_variable("sfcT", sfcT, ("time", "Nx"), "soil temperature", "K") - .. py:method:: merge_files_from_list(Ncfilename_list, exclude_var=[]) diff --git a/docs/source/autoapi/amescap/Script_utils/index.rst b/docs/source/autoapi/amescap/Script_utils/index.rst index 51c30fd8..6e2e6829 100644 --- a/docs/source/autoapi/amescap/Script_utils/index.rst +++ b/docs/source/autoapi/amescap/Script_utils/index.rst @@ -18,7 +18,9 @@ * ``re`` * ``os`` * ``subprocess`` - + * ``sys`` + * ``math`` + * ``matplotlib.colors`` @@ -35,23 +37,18 @@ Functions amescap.Script_utils.MY_func amescap.Script_utils.ak_bk_loader amescap.Script_utils.alt_FV3path + amescap.Script_utils.check_bounds amescap.Script_utils.check_file_tape amescap.Script_utils.dkass_dust_cmap amescap.Script_utils.dkass_temp_cmap + amescap.Script_utils.except_message amescap.Script_utils.extract_path_basename amescap.Script_utils.filter_vars amescap.Script_utils.find_fixedfile amescap.Script_utils.find_tod_in_diurn amescap.Script_utils.get_Ncdf_path amescap.Script_utils.get_longname_unit - amescap.Script_utils.give_permission amescap.Script_utils.hot_cold_cmap - amescap.Script_utils.prCyan - amescap.Script_utils.prGreen - amescap.Script_utils.prLightPurple - amescap.Script_utils.prPurple - amescap.Script_utils.prRed - amescap.Script_utils.prYellow amescap.Script_utils.pretty_print_to_fv_eta amescap.Script_utils.print_fileContent amescap.Script_utils.print_varContent @@ -83,305 +80,381 @@ Attributes .. py:function:: FV3_file_type(fNcdf) - Return the type of the netCDF file (i.e., ``fixed``, ``diurn``, - ``average``, ``daily``) and the format of the Ls array ``areo`` - (i.e., ``fixed``, ``continuous``, or ``diurn``). + Return the type of the netCDF file. - :param fNcdf: an open Netcdf file - :type fNcdf: Netcdf file object + Returns netCDF file type (i.e., ``fixed``, ``diurn``, ``average``, + ``daily``) and the format of the Ls array ``areo`` (i.e., ``fixed``, + ``continuous``, or ``diurn``). + :param fNcdf: an open Netcdf file + :type fNcdf: Netcdf file object :return: The Ls array type (string, ``fixed``, ``continuous``, or ``diurn``) and the netCDF file type (string ``fixed``, ``diurn``, ``average``, or ``daily``) - + :rtype: tuple + :raises ValueError: if the file is not a netCDF file + :raises FileNotFoundError: if the file does not exist + :raises TypeError: if the file is not a Dataset or MFDataset + :raises AttributeError: if the file does not have the _files or + filepath attribute .. py:function:: MY_func(Ls_cont) - Returns the Mars Year + Returns the Mars Year. :param Ls_cont: solar longitude (continuous) - :type Ls_cont: array + :type Ls_cont: array :return: the Mars year - :rtype: int + :rtype: int + :raises ValueError: if Ls_cont is not in the range [0, 360) .. py:function:: ak_bk_loader(fNcdf) - Return ``ak`` and ``bk`` arrays from the current netCDF file. If - these are not found in the current file, search the fixed file in - the same directory. If not there, then search the tiled fixed files. + Loads the ak and bk variables from a netCDF file. - :param fNcdf: an open netCDF file - :type fNcdf: a netCDF file object + This function will first check the current netCDF file for the + ``ak`` and ``bk`` variables. If they are not found, it will + search the fixed file in the same directory. If they are still + not found, it will search the tiled fixed files. The function + will return the ``ak`` and ``bk`` arrays. + :param fNcdf: an open netCDF file + :type fNcdf: a netCDF file object :return: the ``ak`` and ``bk`` arrays + :rtype: tuple + :raises ValueError: if the ``ak`` and ``bk`` variables are not + found in the netCDF file, the fixed file, or the tiled fixed files .. note:: This routine will look for both ``ak`` and ``bk``. There - are cases when it is convenient to load the ``ak``, ``bk`` once - when the files are first opened in ``MarsVars``, but the ``ak`` - and ``bk`` arrays may not be necessary for in the calculation - as is the case for ``MarsVars XXXXX.atmos_average_psd.nc - --add msf``, which operates on a pressure interpolated + are cases when it is convenient to load the ``ak``, ``bk`` once + when the files are first opened in ``MarsVars``, but the ``ak`` + and ``bk`` arrays may not be necessary for in the calculation + as is the case for ``MarsVars XXXXX.atmos_average_psd.nc + --add msf``, which operates on a pressure interpolated (``_pstd.nc``) file. - .. py:function:: alt_FV3path(fullpaths, alt, test_exist=True) - Returns the original or fixed file given an interpolated daily, - diurn or average file. + Returns the original or fixed file for a given path. :param fullpaths: full path to a file or a list of full paths to more than one file - :type fullpaths: str - + :type fullpaths: str :param alt: type of file to return (i.e., original or fixed) - :type alt: str - + :type alt: str :param test_exist: Whether file exists on the disk, defaults to True - :type test_exist: bool, optional + :type test_exist: bool, optional + :return: path to original or fixed file (e.g., + /u/path/00010.atmos_average.nc or /u/path/00010.fixed.nc) + :rtype: str + :raises ValueError: if the file is not a netCDF file + :raises FileNotFoundError: if the file does not exist + :raises TypeError: if the file is not a Dataset or MFDataset + :raises AttributeError: if the file does not have the _files or + filepath attribute + :raises OSError: if the file is not a valid path - :raises ValueError: _description_ - :return: path to original or fixed file - (e.g., "/u/path/00010.atmos_average.nc" or - "/u/path/00010.fixed.nc") - :rtype: str +.. py:function:: check_bounds(values, min_val, max_val, dx) + Checks the bounds of a variable in a netCDF file. + This function checks if the values in a netCDF file are within + the specified bounds. If any value is out of bounds, it will + print an error message and exit the program. + The function can handle both single values and arrays. -.. py:function:: check_file_tape(fileNcdf, abort=False) + Parameters: + :param values: Single value or array of values to check + :type values: array-like + :param min_val: Minimum allowed value + :type min_val: float + :param max_val: Maximum allowed value + :type max_val: float + :return values: The validated value(s) + :rtype: array or float + :raises ValueError: if values is out of bounds or if values is not + a number, array, or list - For use in the NAS environnment only. - Checks whether a file is exists on the disk by running the command - ``dmls -l`` on NAS. This prevents the program from stalling if the - files need to be migrated from the disk to the tape. - :param fileNcdf: full path to a netcdf file or a file object with a name attribute - :type fileNcdf: str or file object +.. py:function:: check_file_tape(fileNcdf) + + Checks whether a file exists on the disk. - :param abort: If True, exit the program. Defaults to False - :type abort: bool, optional + If on a NAS system, also checks if the file needs to be migrated + from tape. + + :param fileNcdf: full path to a netcdf file or a file object with a name attribute + :type fileNcdf: str or file object :return: None + :raises ValueError: if the file is not a netCDF file + :raises FileNotFoundError: if the file does not exist + :raises subprocess.CalledProcessError: if the dmls command fails .. py:function:: dkass_dust_cmap() - Returns a color map useful for dust cross-sections. - (yellow -> orange -> red -> purple) - Provided by Courtney Batterson. + Color map (From Courtney Batterson). + + Returns a color map useful for dust cross-sections that highlight + dust mixing ratios > 4 ppm. The color map goes from + white -> yellow -> orange -> red -> purple. + [Courtney Batterson, Nov 2022] + :return: color map + :rtype: array .. py:function:: dkass_temp_cmap() - Returns a color map that highlights the 200K temperatures. - (black -> purple -> blue -> green -> yellow -> orange -> red) - Provided by Courtney Batterson. + Color map (From Courtney Batterson). + Returns a color map useful for highlighting the 200 K temperature + level. The color map goes from + black -> purple -> blue -> green -> yellow -> orange -> red. + [Courtney Batterson, Nov 2022] + + :return: color map + :rtype: array + + +.. py:function:: except_message(debug, exception, varname, ifile, pre='', ext='') + + Prints an error message if a variable is not found. + + It also contains a special error in the case of a pre-existing + variable. + + :param debug: Flag for debug mode + :type debug: logical + :param exception: Exception from try statement + :type exception: class object + :param varname: Name of variable causing exception + :type varname: string + :param ifile: Name of input file + :type ifile: string + :param pre: Prefix to new variable + :type pre: string + :param ext: Extension to new variable + :type ext: string + :return: None + :rtype: None + :raises ValueError: if debug is True, exception is not a class + object or string, varname is not a string, ifile is not a + string, pre is not a string, or ext is not a string .. py:function:: extract_path_basename(filename) - Returns the path and basename of a file. If only the filename is - provided, assume it is in the current directory. + Returns the path and basename of a file. - :param filename: name of the netCDF file (may include full path) - :type filename: str + If only the filename is provided, assume it is in the current + directory. + :param filename: name of the netCDF file (may include full path) + :type filename: str :return: full file path & name of file + :rtype: tuple + :raises ValueError: if the filename is not a string + :raises FileNotFoundError: if the file does not exist + :raises TypeError: if the filename is not a string + :raises OSError: if the filename is not a valid path .. note:: - This routine does not confirm that the file exists. - It operates on the provided input string. - + This routine does not confirm that the file exists. It operates + on the provided input string. .. py:function:: filter_vars(fNcdf, include_list=None, giveExclude=False) - Filters the variable names in a netCDF file for processing. Returns - all dimensions (``lon``, ``lat``, etc.), the ``areo`` variable, and - any other variable listed in ``include_list``. + Filters the variable names in a netCDF file for processing. + + Returns all dimensions (``lon``, ``lat``, etc.), the ``areo`` + variable, and any other variable listed in ``include_list``. :param fNcdf: an open netCDF object for a diurn, daily, or average file - :type fNcdf: netCDF file object - + :type fNcdf: netCDF file object :param include_list:list of variables to include (e.g., [``ucomp``, ``vcomp``], defaults to None - :type include_list: list or None, optional - + :type include_list: list or None, optional :param giveExclude: if True, returns variables to be excluded from the file, defaults to False - :type giveExclude: bool, optional - + :type giveExclude: bool, optional :return: list of variable names to include in the processed file - + :rtype: list + :raises ValueError: if the file is not a netCDF file + :raises FileNotFoundError: if the file does not exist + :raises TypeError: if the file is not a Dataset or MFDataset + :raises AttributeError: if the file does not have the _files or + filepath attribute + :raises OSError: if the file is not a valid path + :raises KeyError: if the variable is not found in the file .. py:function:: find_fixedfile(filename) - Finds the relevant fixed file for a given average, daily, or diurn - file. - [Batterson, Updated by Alex Nov 29 2022] + Finds the relevant fixed file for an average, daily, or diurn file. - :param filename: an average, daily, or diurn netCDF file - :type filename: str + [Courtney Batterson, updated by Alex Nov 29 2022] + :param filename: an average, daily, or diurn netCDF file + :type filename: str :return: full path to the correspnding fixed file - :rtype: str + :rtype: str + :raises ValueError: if the file is not a netCDF file + :raises FileNotFoundError: if the file does not exist + :raises TypeError: if the file is not a Dataset or MFDataset + :raises AttributeError: if the file does not have the _files or + filepath attribute + :raises OSError: if the file is not a valid path Compatible file types:: - DDDDD.atmos_average.nc -> DDDDD.fixed.nc - atmos_average.tileX.nc -> fixed.tileX.nc - DDDDD.atmos_average_plevs.nc -> DDDDD.fixed.nc - DDDDD.atmos_average_plevs_custom.nc -> DDDDD.fixed.nc - atmos_average.tileX_plevs.nc -> fixed.tileX.nc - atmos_average.tileX_plevs_custom.nc -> fixed.tileX.nc - atmos_average_custom.tileX_plevs.nc -> fixed.tileX.nc - + DDDDD.atmos_average.nc -> DDDDD.fixed.nc + atmos_average.tileX.nc -> fixed.tileX.nc + DDDDD.atmos_average_plevs.nc -> DDDDD.fixed.nc + DDDDD.atmos_average_plevs_custom.nc -> DDDDD.fixed.nc + atmos_average.tileX_plevs.nc -> fixed.tileX.nc + atmos_average.tileX_plevs_custom.nc -> fixed.tileX.nc + atmos_average_custom.tileX_plevs.nc -> fixed.tileX.nc .. py:function:: find_tod_in_diurn(fNcdf) - Returns the variable for the local time axis in diurn files - (e.g., time_of_day_24). - Original implementation by Victoria H. + Returns the variable for the local time axis in diurn files. + + (e.g., time_of_day_24). Original implementation by Victoria H. :param fNcdf: a netCDF file - :type fNcdf: netCDF file object + :type fNcdf: netCDF file object :return: the name of the time of day dimension - :rtype: str - + :rtype: str + :raises ValueError: if the time of day variable is not found .. py:function:: get_Ncdf_path(fNcdf) Returns the full path for a netCDF file object. - .. note:: - ``Dataset`` and multi-file dataset (``MFDataset``) have - different attributes for the path, hence the need for this - function. + ``Dataset`` and multi-file dataset (``MFDataset``) have different + attributes for the path, hence the need for this function. :param fNcdf: Dataset or MFDataset object - :type fNcdf: netCDF file object + :type fNcdf: netCDF file object :return: string list for the Dataset (MFDataset) - :rtype: str(list) + :rtype: str(list) + :raises TypeError: if the file is not a Dataset or MFDataset + :raises AttributeError: if the file does not have the _files or + filepath attribute + :raises ValueError: if the file is not a netCDF file + :raises FileNotFoundError: if the file does not exist .. py:function:: get_longname_unit(fNcdf, varname) - Returns the longname and unit attributes of a variable in a netCDF - file. If the attributes are unavailable, returns blank strings to - avoid an error. + Returns the longname and unit of a variable. - :param fNcdf: an open netCDF file - :type fNcdf: netCDF file object + If the attributes are unavailable, returns blank strings to avoid + an error. + :param fNcdf: an open netCDF file + :type fNcdf: netCDF file object :param varname: variable to extract attribute from - :type varname: str - + :type varname: str :return: longname and unit attributes - :rtype: str + :rtype: str + :raises ValueError: if the file is not a netCDF file + :raises FileNotFoundError: if the file does not exist + :raises TypeError: if the file is not a Dataset or MFDataset + :raises AttributeError: if the file does not have the _files or + filepath attribute + :raises OSError: if the file is not a valid path + :raises KeyError: if the variable is not found in the file .. note:: Some functions in MarsVars edit the units - (e.g., [kg] -> [kg/m]), therefore the empty string is 4 - characters in length (" " instead of "") to allow for + (e.g., [kg] -> [kg/m]), therefore the empty string is 4 + characters in length (" " instead of "") to allow for editing by ``editing units_txt[:-2]``, for example. - -.. py:function:: give_permission(filename) - - Sets group file permissions for the NAS system - - - .. py:function:: hot_cold_cmap() - Returns Dark blue > light blue>white>yellow>red colormap - Based on Matlab's bipolar colormap - - - -.. py:function:: prCyan(skk) - - -.. py:function:: prGreen(skk) - - -.. py:function:: prLightPurple(skk) - + Returns a color map (From Alex Kling, based on bipolar cmap). -.. py:function:: prPurple(skk) + Color map goes from dark blue -> light blue -> white -> yellow -> red. + Based on Matlab's bipolar colormap. + [Alex Kling, Nov 2022] - -.. py:function:: prRed(skk) - - -.. py:function:: prYellow(skk) + :return: color map + :rtype: array .. py:function:: pretty_print_to_fv_eta(var, varname, nperline=6) - Print the ``ak`` or ``bk`` coefficients for copying to - ``fv_eta.f90``. + Print the ``ak`` or ``bk`` coefficients for copying to ``fv_eta.f90``. - :param var: ak or bk data - :type var: array + The ``ak`` and ``bk`` coefficients are used to calculate the + vertical coordinate transformation. + :param var: ak or bk data + :type var: array :param varname: the variable name ("a" or "b") - :type varname: str - + :type varname: str :param nperline: the number of elements per line, defaults to 6 - :type nperline: int, optional - + :type nperline: int, optional :return: a print statement for copying into ``fv_eta.f90`` - + :rtype: None + :raises ValueError: if varname is not "a" or "b" + :raises ValueError: if nperline is not a positive integer + :raises ValueError: if var is not a 1D array of length NLAY+1 .. py:function:: print_fileContent(fileNcdf) - Prints the contents of a netCDF file to the screen. Variables sorted - by dimension. + Prints the contents of a netCDF file to the screen. + + Variables sorted by dimension. :param fileNcdf: full path to the netCDF file - :type fileNcdf: str + :type fileNcdf: str :return: None - .. py:function:: print_varContent(fileNcdf, list_varfull, print_stat=False) - Print variable contents from a variable in a netCDF file. Requires - a XXXXX.fixed.nc file in the current directory. + Print variable contents from a variable in a netCDF file. - :param fileNcdf: full path to a netcdf file - :type fileNcdf: str + Requires a XXXXX.fixed.nc file in the current directory. + :param fileNcdf: full path to a netcdf file + :type fileNcdf: str :param list_varfull: list of variable names and optional slices (e.g., ``["lon", "ps[:, 10, 20]"]``) - :type list_varfull: list - + :type list_varfull: list :param print_stat: If True, print min, mean, and max. If False, print values. Defaults to False - :type print_stat: bool, optional + :type print_stat: bool, optional :return: None + :raises ValueError: if the variable is not found in the file + :raises FileNotFoundError: if the file is not found + :raises Exception: if the variable is not found in the file + :raises Exception: if the file is not found .. py:function:: progress(k, Nmax) @@ -389,150 +462,164 @@ Attributes Displays a progress bar to monitor heavy calculations. :param k: current iteration of the outer loop - :type k: int - + :type k: int :param Nmax: max iteration of the outer loop - :type Nmax: int - + :type Nmax: int + :return: None + :raises ValueError: if k or Nmax are not integers, k > Nmax, or k < 0 .. py:function:: read_variable_dict_amescap_profile(f_Ncdf=None) - Inspect a Netcdf file and return the name of the variables and - dimensions based on the content of ~/.amescap_profile. + Reads a variable dictionary from the ``amescap_profile`` file. - Calling this function allows to remove hard-coded calls in CAP. - For example, to f.variables['ucomp'] is replaced by - f.variables["ucomp"], with "ucomp" taking the values of'ucomp', 'U' + This function will read the variable dictionary from the + ``amescap_profile`` file and return a dictionary with the + variable names and dimensions. The function will also check + the opened netCDF file for the variable names and dimensions. :param f_Ncdf: An opened Netcdf file object - :type f_Ncdf: File object - - :return: Model, a dictionary with the dimensions and variables, - e.g. "ucomp"='U' or "dim_lat"='latitudes' - - .. note:: - The defaut names for variables are defined in () - parenthesis in ~/.amescap_profile:: - - 'X direction wind [m/s] (ucomp)>' - - The defaut names for dimensions are defined in {} parenthesis in - ~/.amescap_profile:: - - Ncdf Y latitude dimension [integer] {lat}>lats - - The dimensions (lon, lat, pfull, pstd) are loaded in the dictionary - as "dim_lon", "dim_lat" - + :type f_Ncdf: File object + :return: MOD, a class object with the variable names and dimensions + (e.g., ``MOD.ucomp`` = 'U' or ``MOD.dim_lat`` = 'latitudes') + :rtype MOD: class object + :raises ValueError: if the ``amescap_profile`` file is not found + or if the variable dictionary is not found + :raises ValueError: if the variable or dimension name is not found + in the netCDF file + :raises ValueError: if the variable or dimension name is not found + in the``amescap_profile`` .. py:function:: regrid_Ncfile(VAR_Ncdf, file_Nc_in, file_Nc_target) Regrid a netCDF variable from one file structure to another. + Requires a file with the desired file structure to mimic. [Alex Kling, May 2021] :param VAR_Ncdf: a netCDF variable object to regrid (e.g., ``f_in.variable["temp"]``) - :type VAR_Ncdf: netCDF file variable - + :type VAR_Ncdf: netCDF file variable :param file_Nc_in: an open netCDF file to source for the variable (e.g., ``f_in = Dataset("filename", "r")``) - :type file_Nc_in: netCDF file object - + :type file_Nc_in: netCDF file object :param file_Nc_target: an open netCDF file with the desired file structure (e.g., ``f_out = Dataset("filename", "r")``) - :type file_Nc_target: netCDF file object - + :type file_Nc_target: netCDF file object :return: the values of the variable interpolated to the target file grid. - :rtype: array + :rtype: array + :raises ValueError: if the file is not a netCDF file + :raises FileNotFoundError: if the file does not exist + :raises TypeError: if the file is not a Dataset or MFDataset + :raises AttributeError: if the file does not have the _files or + filepath attribute + :raises OSError: if the file is not a valid path .. note:: While the KDTree interpolation can handle a 3D dataset - (lon/lat/lev instead of just 2D lon/lat), the grid points in - the vertical are just a few (10--100s) meters in the PBL vs a - few (10-100s) kilometers in the horizontal. This results in - excessive weighting in the vertical, which is why the vertical + (lon/lat/lev instead of just 2D lon/lat), the grid points in + the vertical are just a few (10--100s) meters in the PBL vs a + few (10-100s) kilometers in the horizontal. This results in + excessive weighting in the vertical, which is why the vertical dimension is handled separately. - .. py:function:: replace_dims(Ncvar_dim, vert_dim_name=None) + Replaces the dimensions of a variable in a netCDF file. + Updates the name of the variable dimension to match the format of the new NASA Ames Mars GCM output files. :param Ncvar_dim: netCDF variable dimensions (e.g., ``f_Ncdf.variables["temp"].dimensions``) - :type Ncvar_dim: str - + :type Ncvar_dim: str :param vert_dim_name: the vertical dimension if it is ambiguous (``pstd``, ``zstd``, or ``zagl``). Defaults to None - :type vert_dim_name: str, optional - + :type vert_dim_name: str, optional :return: updated dimensions - :rtype: str - + :rtype: str + :raises ValueError: if Ncvar_dim is not a list + :raises ValueError: if vert_dim_name is not a string + :raises ValueError: if vert_dim_name is not in the list of + recognized vertical dimensions .. py:function:: reset_FV3_names(MOD) - This function reset the model dictionary to the native FV3's - variables, e.g.:: - + Resets the FV3 variable names in a netCDF file. + + This function resets the model dictionary to the native FV3 + variable names, e.g.:: + model.dim_lat = 'latitude' > model.dim_lat = 'lat' model.ucomp = 'U' > model.ucomp = 'ucomp' :param MOD: Generated with read_variable_dict_amescap_profile() - :type MOD: class object - - :return: same object with updated names for the dimensions and + :type MOD: class object + :return: same object with updated names for the dimensions and variables - + :rtype: class object + :raises ValueError: if the MOD object is not a class object or does + not contain the expected attributes .. py:function:: rjw_cmap() - Returns John Wilson's preferred color map - (red -> jade -> wisteria) + Returns a color map (From R. John Wilson). + Color map goes from red -> jade -> wisteria. + [R. John Wilson, Nov 2022] + + :return: color map + :rtype: array .. py:function:: section_content_amescap_profile(section_ID) - Executes first code section in ``~/.amescap_profile`` to read in - user-defined plot & interpolation settings. + Executes first code section in ``~/.amescap_profile``. + + Reads in user-defined plot & interpolation settings. + [Alex Kling, Nov 2022] :param section_ID: the section to load (e.g., Pressure definitions for pstd) - :type section_ID: str - + :type section_ID: str :return: the relevant line with Python syntax - + :rtype: str + :raises FileNotFoundError: if the file is not found .. py:function:: smart_reader(fNcdf, var_list, suppress_warning=False) + Reads a variable from a netCDF file. + + If the variable is not found in the file, it checks for the variable + in the original file (e.g., atmos_average.nc) or the fixed file + (e.g., 00010.fixed.nc). + Alternative to ``var = fNcdf.variables["var"][:]`` for handling - *processed* files that also checks for a matching average or daily - and XXXXX.fixed.nc file. + *processed* files. :param fNcdf: an open netCDF file - :type fNcdf: netCDF file object - + :type fNcdf: netCDF file object :param var_list: a variable or list of variables (e.g., ``areo`` or [``pk``, ``bk``, ``areo``]) - :type var_list: _type_ - + :type var_list: _type_ :param suppress_warning: suppress debug statement. Useful if a variable is not expected to be in the file anyway. Defaults to False - :type suppress_warning: bool, optional - + :type suppress_warning: bool, optional :return: variable content (single or values to unpack) - :rtype: list + :rtype: list + :raises ValueError: if the file is not a netCDF file + :raises FileNotFoundError: if the file does not exist + :raises TypeError: if the file is not a Dataset or MFDataset + :raises AttributeError: if the file does not have the _files or + filepath attribute + :raises OSError: if the file is not a valid path Example:: @@ -554,12 +641,15 @@ Attributes Only the variable content is returned, not attributes - .. py:function:: wbr_cmap() - Returns a color map that goes from - white -> blue -> green -> yellow -> red + Returns a color map (From R. John Wilson). + + Color map goes from white -> blue -> green -> yellow -> red + [R. John Wilson, Nov 2022] + :return: color map + :rtype: array .. py:data:: Blue diff --git a/docs/source/autoapi/amescap/Spectral_utils/index.rst b/docs/source/autoapi/amescap/Spectral_utils/index.rst index 77d68b5d..b372fb37 100644 --- a/docs/source/autoapi/amescap/Spectral_utils/index.rst +++ b/docs/source/autoapi/amescap/Spectral_utils/index.rst @@ -15,7 +15,7 @@ * ``scipy.signal`` * ``ImportError`` * ``Exception`` - + * ``pyshtools`` (optional) @@ -29,9 +29,13 @@ Functions .. autoapisummary:: amescap.Spectral_utils.diurn_extract + amescap.Spectral_utils.extract_diurnal_harmonics + amescap.Spectral_utils.import_pyshtools amescap.Spectral_utils.reconstruct_diurn - amescap.Spectral_utils.space_time amescap.Spectral_utils.zeroPhi_filter + amescap.Spectral_utils.zonal_construct + amescap.Spectral_utils.zonal_decomposition + .. py:function:: diurn_extract(VAR, N, tod, lon) @@ -42,88 +46,50 @@ Functions :param VAR: field to process. Time of day dimension must be first (e.g., ``[tod, time, lat, lon]`` or ``[tod]`` :type VAR: 1D or ND array - :param N: number of harmonics to extract (``N=1`` for diurnal, ``N=2`` for diurnal AND semi diurnal, etc.) :type N: int - :param tod: universal time of day in sols (``0-1.``) If provided in hours (``0-24``), it will be normalized. :type tod: 1D array - :param lon: longitudes ``0-360`` :type lon: 1D array or float - :return: the amplitudes & phases for the Nth first harmonics, (e.g., size ``[Nh, time, lat, lon]``) :rtype: ND arrays - -.. py:function:: reconstruct_diurn(amp, phas, tod, lon, sumList=[]) - - Reconstructs a field wave based on its diurnal harmonics - - :param amp: amplitude of the signal. Harmonics dimension FIRST - (e.g., ``[N, time, lat, lon]``) - :type amp: array - - :param phas: phase of the signal [hr]. Harmonics dimension FIRST - :type phas: array - - :param tod: time of day in universal time [hr] - :type tod: 1D array - - :param lon: longitude for converting universal -> local time - :type lon: 1D array or float - - :param sumList: the harmonics to include when reconstructing the - wave (e.g., ``sumN=[1, 2, 4]``), defaults to ``[]`` - :type sumList: list, optional - - :return: a variable with reconstructed harmonics with N dimension - FIRST and time of day SECOND (``[N, tod, time, lat, lon]``). If - sumList is provided, the wave output harmonics will be - aggregated (i.e., size = ``[tod, time, lat, lon]``) - :rtype: _type_ - - - -.. py:function:: space_time(lon, timex, varIN, kmx, tmx) +.. py:function:: extract_diurnal_harmonics(kmx, tmx, varIN, tod, lon) Obtain west and east propagating waves. This is a Python implementation of John Wilson's ``space_time`` routine. - Alex Kling 2019. + Richard Urata 2025. :param lon: longitude [°] (0-360) :type lon: 1D array - - :param timex: time [sol] (e.g., 1.5 days sampled every hour is - ``[0/24, 1/24, 2/24,.. 1,.. 1.5]``) - :type timex: 1D array - - :param varIN: variable for the Fourier analysis. First axis must be - ``lon`` and last axis must be ``time`` (e.g., - ``varIN[lon, time]``, ``varIN[lon, lat, time]``, or - ``varIN[lon, lev, lat, time]``) + :param tod: time_of_day [hour] (e.g., a 24 hour sol sampled every hour could be + ``[0.5, 1.5, 2.5, ..., 23.5]``) + :type tod: 1D array + :param varIN: variable for the Fourier analysis, expected to be scaled by zonal mean. + First axis must be + ``time`` and last axis must be ``lon`` (e.g., + ``varIN[time, tod, lon]``, ``varIN[time, tod, lat, lon]``, or + ``varIN[time, tod, lat, lon, lev]``) :type varIN: array - :param kmx: the number of longitudinal wavenumbers to extract (max = ``nlon/2``) :type kmx: int - :param tmx: the number of tidal harmonics to extract (max = ``nsamples/2``) :type tmx: int - :return: (ampe) East propagating wave amplitude [same unit as varIN]; (ampw) West propagating wave amplitude [same unit as varIN]; (phasee) East propagating phase [°]; (phasew) West propagating phase [°] - .. NOTE:: 1. ``ampe``, ``ampw``, ``phasee``, and ``phasew`` have + .. note:: 1. ``ampe``, ``ampw``, ``phasee``, and ``phasew`` have dimensions ``[kmx, tmx]`` or ``[kmx, tmx, lat]`` or - ``[kmx, tmx, lev, lat]`` etc. + ``[kmx, tmx, lat, lev]`` etc. 2. The x and y axes may be constructed as follows, which will display the eastern and western modes:: @@ -133,7 +99,39 @@ Functions KTIME, KLON = np.meshgrid(ktime, klon) amplitude = np.concatenate((ampw[:, ::-1], ampe), axis = 1) phase = np.concatenate((phasew[:, ::-1], phasee), axis = 1) - + + +.. py:function:: import_pyshtools() + + Attempt to import pyshtools and set the global variable + PYSHTOOLS_AVAILABLE accordingly. + Check if pyshtools is available and return its availability status + + :return: True if pyshtools is available, False otherwise + :rtype: bool + + +.. py:function:: reconstruct_diurn(amp, phas, tod, lon, sumList=[]) + + Reconstructs a field wave based on its diurnal harmonics + + :param amp: amplitude of the signal. Harmonics dimension FIRST + (e.g., ``[N, time, lat, lon]``) + :type amp: array + :param phas: phase of the signal [hr]. Harmonics dimension FIRST + :type phas: array + :param tod: time of day in universal time [hr] + :type tod: 1D array + :param lon: longitude for converting universal -> local time + :type lon: 1D array or float + :param sumList: the harmonics to include when reconstructing the + wave (e.g., ``sumN=[1, 2, 4]``), defaults to ``[]`` + :type sumList: list, optional + :return: a variable with reconstructed harmonics with N dimension + FIRST and time of day SECOND (``[N, tod, time, lat, lon]``). If + sumList is provided, the wave output harmonics will be + aggregated (i.e., size = ``[tod, time, lat, lon]``) + :rtype: _type_ .. py:function:: zeroPhi_filter(VAR, btype, low_highcut, fs, axis=0, order=4, add_trend=False) @@ -144,30 +142,69 @@ Functions :param VAR: values for filtering 1D or ND array. Filtered dimension must be FIRST. Adjusts axis as necessary. :type VAR: array - :param btype: filter type (i.e., "low", "high" or "band") :type btype: str - :param low_high_cut: low, high, or [low, high] cutoff frequency depending on the filter [Hz or m-1] :type low_high_cut: int - :param fs: sampling frequency [Hz or m-1] :type fs: int - :param axis: if data is an ND array, this identifies the filtering dimension :type axis: int - :param order: order for the filter :type order: int - :param add_trend: if True, return the filtered output. If false, return the trend and filtered output. :type add_trend: bool - :return: the filtered data - .. NOTE:: ``Wn=[low, high]`` are expressed as a function of the + .. note:: ``Wn=[low, high]`` are expressed as a function of the Nyquist frequency - + + +.. py:function:: zonal_construct(COEFFS_flat, VAR_shape, btype=None, low_highcut=None) + + Recomposition into spherical harmonics + + :param COEFFS_flat: Spherical harmonic coefficients as a flattened + array, (e.g., ``[time, lat, lon]`` or + ``[time x lev, 2, lat, lon]``) + :type COEFFS_flat: array + :param VAR_shape: shape of the original variable + :type VAR_shape: tuple + :param btype: filter type: "low", "high", or "band". If None, + returns reconstructed array using all zonal wavenumbers + :type btype: str or None + :param low_high_cut: low, high or [low, high] zonal wavenumber + :type low_high_cut: int or list + :return: detrended variable reconstructed to original size + (e.g., [time, lev, lat, lon]) + + .. note:: The minimum and maximum wavelenghts in [km] are computed:: + dx = 2*np.pi * 3400 + L_min = (1./kmax) * dx + L_max = 1./max(kmin, 1.e-20) * dx + if L_max > 1.e20: + L_max = np.inf + print("(kmin,kmax) = ({kmin}, {kmax}) + -> dx min = {L_min} km, dx max = {L_max} km") + + +.. py:function:: zonal_decomposition(VAR) + + Decomposition into spherical harmonics. [A. Kling, 2020] + + :param VAR: Detrend variable for decomposition. Lat is SECOND to + LAST dimension and lon is LAST (e.g., ``[time,lat,lon]`` or + ``[time,lev,lat,lon]``) + :return: (COEFFS) coefficient for harmonic decomposion, shape is + flattened (e.g., ``[time, 2, lat/2, lat/2]`` + ``[time x lev, 2, lat/2, lat/2]``); + (power_per_l) power spectral density, shape is re-organized + (e.g., [time, lat/2] or [time, lev, lat/2]) + + .. note:: Output size is (``[...,lat/2, lat/2]``) as lat is the + smallest dimension. This matches the Nyquist frequency. + + diff --git a/docs/source/autoapi/amescap/pdf2image/index.rst b/docs/source/autoapi/amescap/pdf2image/index.rst index 1e2853b9..032a4bdb 100644 --- a/docs/source/autoapi/amescap/pdf2image/index.rst +++ b/docs/source/autoapi/amescap/pdf2image/index.rst @@ -18,7 +18,8 @@ * ``os`` * ``subprocess`` * ``PIL`` - + * ``uuid`` + * ``Pillow`` @@ -38,38 +39,30 @@ Functions .. py:function:: convert_from_bytes(pdf_file, dpi=200, output_folder=None, first_page=None, last_page=None, fmt='ppm', thread_count=1, userpw=None, use_cropbox=False) - Convert PDF to Image will throw an error whenever one of the condition is reached + Convert PDF to Image will throw an error whenever one of the + condition is reached :param pdf_file: Bytes representing the PDF file :type pdf_file: float - :param dpi: image quality in DPI (default 200) :type dpi: int - :param output_folder: folder to write the images to (instead of directly in memory) :type output_folder: str - :param first_page: first page to process :type first_page: int - :param last_page: last page to process before stopping :type last_page: int - :param fmt: output image format :type fmt: str - :param thread_count: how many threads to spawn for processing :type thread_count: int - :param userpw: PDF password :type userpw: str - :param use_cropbox: use cropbox instead of mediabox :type use_cropbox: bool - .. py:function:: convert_from_path(pdf_path, dpi=200, output_folder=None, first_page=None, last_page=None, fmt='ppm', thread_count=1, userpw=None, use_cropbox=False) Convert PDF to Image will throw an error whenever one of the @@ -77,31 +70,22 @@ Functions :param pdf_path: path to the PDF that you want to convert :type pdf_path: str - :param dpi: image quality in DPI (default 200) :type dpi: int - :param output_folder: folder to write the images to (instead of directly in memory) :type output_folder: str - :param first_page: first page to process :type first_page: int - :param last_page: last page to process before stopping :type last_page: int - :param fmt: output image format :type fmt: str - :param thread_count: how many threads to spawn for processing :type thread_count: int - :param userpw: PDF password :type userpw: str - :param use_cropbox: use cropbox instead of mediabox :type use_cropbox: bool - diff --git a/docs/source/autoapi/bin/MarsCalendar/index.rst b/docs/source/autoapi/bin/MarsCalendar/index.rst index c480121d..11797a4c 100644 --- a/docs/source/autoapi/bin/MarsCalendar/index.rst +++ b/docs/source/autoapi/bin/MarsCalendar/index.rst @@ -10,18 +10,22 @@ The executable requires 1 of the following arguments: - * ``[-sol --sol]`` The sol to convert to Ls, OR - * ``[-ls --ls]`` The Ls to convert to sol + * ``[-sol --sol]`` The sol to convert to Ls, OR + * ``[-ls --ls]`` The Ls to convert to sol and optionally accepts: - * ``[-my --marsyear]`` The Mars Year of the simulation to compute sol or Ls from, AND/OR - * ``[-c --continuous]`` Returns Ls in continuous form + * ``[-my --marsyear]`` The Mars Year of the simulation to compute sol or Ls from, AND/OR + * ``[-c --continuous]`` Returns Ls in continuous form Third-party Requirements: * ``numpy`` * ``argparse`` + * ``functools`` + * ``traceback`` + * ``sys`` + * ``amescap`` @@ -34,6 +38,7 @@ Functions .. autoapisummary:: + bin.MarsCalendar.debug_wrapper bin.MarsCalendar.main bin.MarsCalendar.parse_array @@ -45,13 +50,81 @@ Attributes .. autoapisummary:: bin.MarsCalendar.args + bin.MarsCalendar.debug bin.MarsCalendar.exclusive_group + bin.MarsCalendar.exit_code bin.MarsCalendar.group bin.MarsCalendar.parser +.. py:function:: debug_wrapper(func) + + A decorator that wraps a function with error handling + based on the --debug flag. + + If the --debug flag is set, it prints the full traceback + of any exception that occurs. Otherwise, it prints a + simplified error message. + + :param func: The function to wrap. + :type func: function + :return: The wrapped function. + :rtype: function + :raises Exception: If an error occurs during the function call. + :raises TypeError: If the function is not callable. + :raises ValueError: If the function is not found. + :raises NameError: If the function is not defined. + :raises AttributeError: If the function does not have the + specified attribute. + :raises ImportError: If the function cannot be imported. + :raises RuntimeError: If the function cannot be run. + :raises KeyError: If the function does not have the + specified key. + :raises IndexError: If the function does not have the + specified index. + :raises IOError: If the function cannot be opened. + :raises OSError: If the function cannot be accessed. + :raises EOFError: If the function cannot be read. + :raises MemoryError: If the function cannot be allocated. + :raises OverflowError: If the function cannot be overflowed. + :raises ZeroDivisionError: If the function cannot be divided by zero. + :raises StopIteration: If the function cannot be stopped. + :raises KeyboardInterrupt: If the function cannot be interrupted. + :raises SystemExit: If the function cannot be exited. + :raises AssertionError: If the function cannot be asserted. + + .. py:function:: main() + Main function for MarsCalendar command-line tool. + + This function processes user-specified arguments to convert between + Mars solar longitude (Ls) and sol (Martian day) values for a given + Mars year. It supports both continuous and discrete sol + calculations, and can handle input in either Ls or sol, returning + the corresponding converted values. The results are printed in a + formatted table. + + Arguments are expected to be provided via the `args` namespace: + + - args.marsyear: Mars year (default is 0) + - args.continuous: If set, enables continuous sol calculation + - args.ls: List of Ls values to convert to sol + - args.sol: List of sol values to convert to Ls + + :param args: Command-line arguments parsed using argparse. + :type args: argparse.Namespace + :raises ValueError: If the input is not a valid number or + range. + :returns: 0 if the program runs successfully, 1 if an error occurs. + :rtype: int + :raises RuntimeError: If the input is not a valid runtime. + :raises TypeError: If the input is not a valid type. + :raises IndexError: If the input is not a valid index. + :raises KeyError: If the input is not a valid key. + :raises AttributeError: If the input is not a valid attribute. + :raises ImportError: If the input is not a valid import. + .. py:function:: parse_array(len_input) @@ -63,20 +136,31 @@ Attributes :param len_input: The input Ls or sol to convert. Can either be one number (e.g., 300) or start stop step (e.g., 300 310 2). - :type len_input: float - + :type len_input: float or list of floats :raises: Error if neither ``[-ls --ls]`` or ``[-sol --sol]`` are provided. - :return: ``input_as_arr`` An array formatted for input into ``ls2sol`` or ``sol2ls``. If ``len_input = 300``, then ``input_as_arr=[300]``. If ``len_input = 300 310 2``, then ``input_as_arr = [300, 302, 304, 306, 308]``. + :rtype: list of floats + :raises ValueError: If the input is not a valid number or + range. + :raises TypeError: If the input is not a valid type. + :raises IndexError: If the input is not a valid index. + :raises KeyError: If the input is not a valid key. + :raises AttributeError: If the input is not a valid attribute. + :raises ImportError: If the input is not a valid import. + :raises RuntimeError: If the input is not a valid runtime. + :raises OverflowError: If the input is not a valid overflow. + :raises MemoryError: If the input is not a valid memory. +.. py:data:: args + -.. py:data:: args +.. py:data:: debug @@ -84,6 +168,10 @@ Attributes +.. py:data:: exit_code + + + .. py:data:: group diff --git a/docs/source/autoapi/bin/MarsFiles/index.rst b/docs/source/autoapi/bin/MarsFiles/index.rst index e0f916ab..276e9688 100644 --- a/docs/source/autoapi/bin/MarsFiles/index.rst +++ b/docs/source/autoapi/bin/MarsFiles/index.rst @@ -18,6 +18,7 @@ * ``[-bin, --bin_files]`` Produce MGCM 'fixed', 'diurn', 'average' and 'daily' files from Legacy output * ``[-c, --concatenate]`` Combine sequential files of the same type into one file + * ``[-split, --split]`` Split file along a specified dimension or extracts slice at one point along the dim * ``[-t, --time_shift]`` Apply a time-shift to 'diurn' files * ``[-ba, --bin_average]`` Bin MGCM 'daily' files like 'average' files * ``[-bd, --bin_diurn]`` Bin MGCM 'daily' files like 'diurn' files @@ -30,7 +31,8 @@ * ``[-bps, --band_pass_spatial]`` Spatial filter: band-pass * ``[-tide, --tide_decomp]`` Extract diurnal tide and its harmonics * ``[-recon, --reconstruct]`` Reconstruct the first N harmonics - * ``[-norm, --normalize]`` Provide ``-tide`` result in % amplitude + * ``[-norm, --normalize]`` Provide ``-tide`` result in percent amplitude + * ``[-prop, --prop_tides]`` Extract propagating tide harmonics * ``[-regrid, --regrid_XY_to_match]`` Regrid a target file to match a source file * ``[-zavg, --zonal_average]`` Zonally average all variables in a file * ``[-incl, --include]`` Only include specific variables in a calculation @@ -38,13 +40,16 @@ Third-party Requirements: - * ``numpy`` - * ``netCDF4`` * ``sys`` * ``argparse`` * ``os`` - * ``subprocess`` * ``warnings`` + * ``re`` + * ``numpy`` + * ``netCDF4`` + * ``shutil`` + * ``functools`` + * ``traceback`` @@ -68,6 +73,7 @@ Functions bin.MarsFiles.change_vname_longname_unit bin.MarsFiles.concatenate_files + bin.MarsFiles.debug_wrapper bin.MarsFiles.do_avg_vars bin.MarsFiles.ls2sol_1year bin.MarsFiles.main @@ -86,6 +92,8 @@ Attributes bin.MarsFiles.all_args bin.MarsFiles.args + bin.MarsFiles.debug + bin.MarsFiles.exit_code bin.MarsFiles.out_ext bin.MarsFiles.out_ext bin.MarsFiles.parser @@ -96,54 +104,35 @@ Attributes Bases: :py:obj:`argparse.Action` - Information about how to convert command line strings to Python objects. - - Action objects are used by an ArgumentParser to represent the information - needed to parse a single argument from one or more strings from the - command line. The keyword arguments to the Action constructor are also - all attributes of Action instances. - - Keyword Arguments: - - - option_strings -- A list of command-line option strings which - should be associated with this action. - - - dest -- The name of the attribute to hold the created object(s) - - - nargs -- The number of command-line arguments that should be - consumed. By default, one argument will be consumed and a single - value will be produced. Other values include: - - N (an integer) consumes N arguments (and produces a list) - - '?' consumes zero or one arguments - - '*' consumes zero or more arguments (and produces a list) - - '+' consumes one or more arguments (and produces a list) - Note that the difference between the default and nargs=1 is that - with the default, a single value will be produced, while with - nargs=1, a list containing a single value will be produced. - - - const -- The value to be produced if the option is specified and the - option uses an action that takes no values. - - - default -- The value to be produced if the option is not specified. - - - type -- A callable that accepts a single string argument, and - returns the converted value. The standard Python types str, int, - float, and complex are useful examples of such callables. If None, - str is used. - - - choices -- A container of values that should be allowed. If not None, - after a command-line argument has been converted to the appropriate - type, an exception will be raised if it is not a member of this - collection. - - - required -- True if the action must always be specified at the - command line. This is only meaningful for optional command-line - arguments. - - - help -- The help string describing the argument. - - - metavar -- The name to be used for the option's argument with the - help string. If None, the 'dest' value will be used as the name. + Custom action for argparse to handle file extensions. + + This action is used to add an extension to the output file. + + :param ext_content: The content to be added to the file name. + :type ext_content: str + :param parser: The parser instance to which this action belongs. + :type parser: argparse.ArgumentParser + :param args: Additional positional arguments. + :type args: tuple + :param kwargs: Additional keyword arguments. + :type kwargs: dict + :param ext_content: The content to be added to the file name + :type ext_content: str + :return: None + :rtype: None + :raises ValueError: If ext_content is not provided. + :raises TypeError: If parser is not an instance of argparse.ArgumentParser. + :raises Exception: If an error occurs during the action. + :raises AttributeError: If the parser does not have the specified attribute. + :raises ImportError: If the parser cannot be imported. + :raises RuntimeError: If the parser cannot be run. + :raises KeyError: If the parser does not have the specified key. + :raises IndexError: If the parser does not have the specified index. + :raises IOError: If the parser cannot be opened. + :raises OSError: If the parser cannot be accessed. + :raises EOFError: If the parser cannot be read. + :raises MemoryError: If the parser cannot be allocated. + :raises OverflowError: If the parser cannot be overflowed. .. py:method:: __call__(parser, namespace, values, option_string=None) @@ -162,25 +151,34 @@ Attributes Bases: :py:obj:`argparse.ArgumentParser` - Object for parsing command line strings into Python objects. - - Keyword Arguments: - - prog -- The name of the program (default: - ``os.path.basename(sys.argv[0])``) - - usage -- A usage message (default: auto-generated from arguments) - - description -- A description of what the program does - - epilog -- Text following the argument descriptions - - parents -- Parsers whose arguments should be copied into this one - - formatter_class -- HelpFormatter class for printing help messages - - prefix_chars -- Characters that prefix optional arguments - - fromfile_prefix_chars -- Characters that prefix files containing - additional arguments - - argument_default -- The default value for all arguments - - conflict_handler -- String indicating how to handle conflicts - - add_help -- Add a -h/-help option - - allow_abbrev -- Allow long options to be abbreviated unambiguously - - exit_on_error -- Determines whether or not ArgumentParser exits with - error info when an error occurs + Custom ArgumentParser that handles file extensions for output files. + + This class extends the argparse.ArgumentParser to add functionality + for handling file extensions in the command-line arguments. + + :param args: Additional positional arguments. + :type args: tuple + :param kwargs: Additional keyword arguments. + :type kwargs: dict + :param ext_content: The content to be added to the file name. + :type ext_content: str + :param parser: The parser instance to which this action belongs. + :type parser: argparse.ArgumentParser + :return: None + :rtype: None + :raises ValueError: If ext_content is not provided. + :raises TypeError: If parser is not an instance of argparse.ArgumentParser. + :raises Exception: If an error occurs during the action. + :raises AttributeError: If the parser does not have the specified attribute. + :raises ImportError: If the parser cannot be imported. + :raises RuntimeError: If the parser cannot be run. + :raises KeyError: If the parser does not have the specified key. + :raises IndexError: If the parser does not have the specified index. + :raises IOError: If the parser cannot be opened. + :raises OSError: If the parser cannot be accessed. + :raises EOFError: If the parser cannot be read. + :raises MemoryError: If the parser cannot be allocated. + :raises OverflowError: If the parser cannot be overflowed. .. py:method:: __repr__() @@ -259,16 +257,22 @@ Attributes designed to work specifically with LegacyCGM.nc files. :param vname: variable name - :type vname: str - + :type vname: str :param longname_txt: variable description - :type longname_txt: str - + :type longname_txt: str :param units_txt: variable units - :type units_txt: str - + :type units_txt: str :return: variable name and corresponding description and unit + (e.g. ``vname = "ps"``) + :rtype: tuple + :raises KeyError: if the required variables are not found + :raises ValueError: if the required dimensions are not found + :raises AttributeError: if the required attributes are not found + :note:: + The ``diurn`` file is created by binning the Legacy files. + The ``average`` and ``daily`` files are created by + averaging over the ``diurn`` file. .. py:function:: concatenate_files(file_list, full_file_list) @@ -276,11 +280,55 @@ Attributes Concatenates sequential output files in chronological order. :param file_list: list of file names - :type file_list: list - + :type file_list: list :param full_file_list: list of file names and full paths - :type full_file_list: list - + :type full_file_list: list + :returns: None + :rtype: None + :raises OSError: If the file cannot be removed. + :raises IOError: If the file cannot be moved. + :raises Exception: If the file cannot be opened. + :raises ValueError: If the file cannot be accessed. + :raises TypeError: If the file is not of the correct type. + :raises IndexError: If the file does not have the correct index. + :raises KeyError: If the file does not have the correct key. + + +.. py:function:: debug_wrapper(func) + + A decorator that wraps a function with error handling + based on the --debug flag. + + If the --debug flag is set, it prints the full traceback + of any exception that occurs. Otherwise, it prints a + simplified error message. + + :param func: The function to wrap. + :type func: function + :return: The wrapped function. + :rtype: function + :raises Exception: If an error occurs during the function call. + :raises TypeError: If the function is not callable. + :raises ValueError: If the function is not found. + :raises NameError: If the function is not defined. + :raises AttributeError: If the function does not have the + specified attribute. + :raises ImportError: If the function cannot be imported. + :raises RuntimeError: If the function cannot be run. + :raises KeyError: If the function does not have the + specified key. + :raises IndexError: If the function does not have the + specified index. + :raises IOError: If the function cannot be opened. + :raises OSError: If the function cannot be accessed. + :raises EOFError: If the function cannot be read. + :raises MemoryError: If the function cannot be allocated. + :raises OverflowError: If the function cannot be overflowed. + :raises ZeroDivisionError: If the function cannot be divided by zero. + :raises StopIteration: If the function cannot be stopped. + :raises KeyboardInterrupt: If the function cannot be interrupted. + :raises SystemExit: If the function cannot be exited. + :raises AssertionError: If the function cannot be asserted. .. py:function:: do_avg_vars(histfile, newf, avgtime, avgtod, bin_period=5) @@ -288,100 +336,157 @@ Attributes Performs a time average over all fields in a file. :param histfile: file to perform time average on - :type histfile: str - + :type histfile: str :param newf: path to target file - :type newf: str - + :type newf: str :param avgtime: whether ``histfile`` has averaged fields (e.g., ``atmos_average``) - :type avgtime: bool - + :type avgtime: bool :param avgtod: whether ``histfile`` has a diurnal time dimenion (e.g., ``atmos_diurn``) - :type avgtod: bool - + :type avgtod: bool :param bin_period: the time binning period if `histfile` has averaged fields (i.e., if ``avgtime==True``), defaults to 5 - :type bin_period: int, optional - + :type bin_period: int, optional :return: a time-averaged file + :rtype: None + :raises KeyError: if the required variables are not found + :raises ValueError: if the required dimensions are not found + :raises AttributeError: if the required attributes are not found + :note:: + The ``diurn`` file is created by binning the Legacy files. + The ``average`` and ``daily`` files are created by + averaging over the ``diurn`` file. .. py:function:: ls2sol_1year(Ls_deg, offset=True, round10=True) Returns a sol number from the solar longitude. - :param Ls_deg: solar longitude [°] - :type Ls_deg: float + This is consistent with the MGCM model. The Ls is the solar + longitude in degrees. The sol number is the number of sols since + the perihelion (Ls = 250.99 degrees). + :param Ls_deg: solar longitude [°] + :type Ls_deg: float :param offset: if True, force year to start at Ls 0 - :type offset: bool - + :type offset: bool :param round10: if True, round to the nearest 10 sols - :type round10: bool - + :type round10: bool :returns: ``Ds`` the sol number + :rtype: float + :raises ValueError: if the required variables are not found + :raises KeyError: if the required variables are not found + :raises AttributeError: if the required attributes are not found ..note:: - For the moment, this is consistent with 0 <= Ls <= 359.99, but - not for monotically increasing Ls. - + This is consistent with 0 <= Ls <= 359.99, but not for + monotically increasing Ls. .. py:function:: main() + Main entry point for MarsFiles data processing utility. + + This function processes input NetCDF or legacy MGCM files according + to command-line arguments. It supports a variety of operations, + including: + + - Conversion of legacy MGCM files to FV3 format. + - Concatenation and splitting of NetCDF files along specified + dimensions. + - Temporal binning of daily files into average or diurnal files. + - Temporal filtering (high-pass, low-pass, band-pass) using spectral + methods. + - Spatial (zonal) filtering and decomposition using spherical + harmonics. + - Tidal analysis and harmonic decomposition of diurnal files. + - Regridding of data to match a target NetCDF file's grid. + - Zonal averaging of variables over longitude. + + The function handles file path resolution, argument validation, and + orchestrates the appropriate processing routines based on + user-specified options. Output files are written in NetCDF format, + with new variables and dimensions created as needed. + + Global Variables: + data_dir (str): The working directory for input/output files. + Arguments: + None directly. Uses global 'args' parsed from command-line. + Returns: + None. Outputs are written to disk. + Raises: + SystemExit: For invalid argument combinations or processing + errors. + .. py:function:: make_FV3_files(fpath, typelistfv3, renameFV3=True) Make MGCM-like ``average``, ``daily``, and ``diurn`` files. - Used if call to [``-bin --bin_files``] is made AND Legacy files are in - netCDFformat (not fort.11). - :param fpath: Full path to the Legacy netcdf files - :type fpath: str + Used if call to [``-bin --bin_files``] is made AND Legacy files are + in netCDF format (not fort.11). + :param fpath: Full path to the Legacy netcdf files + :type fpath: str :param typelistfv3: MGCM-like file type: ``average``, ``daily``, or ``diurn`` - :type typelistfv3: list - + :type typelistfv3: list :param renameFV3: Rename the files from Legacy_LsXXX_LsYYY.nc to ``XXXXX.atmos_average.nc`` following MGCM output conventions - :type renameFV3: bool - + :type renameFV3: bool :return: The MGCM-like files: ``XXXXX.atmos_average.nc``, ``XXXXX.atmos_daily.nc``, ``XXXXX.atmos_diurn.nc``. + :rtype: None + :note:: + The ``average`` and ``daily`` files are created by + averaging over the ``diurn`` file. The ``diurn`` file is + created by binning the Legacy files. + + :note:: + The ``diurn`` file is created by binning the Legacy files. .. py:function:: process_time_shift(file_list) - This function converts the data in diurn files with a time_of_day_XX - dimension to universal local time. + Converts diurn files to local time. - :param file_list: list of file names - :type file_list: list + This function is used to convert diurn files to local time. + :param file_list: list of file names + :type file_list: list + :returns: None + :rtype: None + :raises OSError: If the file cannot be removed. + :raises IOError: If the file cannot be moved. + :raises Exception: If the file cannot be opened. + :raises ValueError: If the file cannot be accessed. + :raises TypeError: If the file is not of the correct type. + :raises IndexError: If the file does not have the correct index. .. py:function:: replace_at_index(tuple_dims, idx, new_name) - Updates variable dimensions. + Replaces the dimension at the given index with a new name. + + If ``new_name`` is None, the dimension is removed. + This is designed to work specifically with LegacyCGM.nc files. :param tuple_dims: the dimensions as tuples e.g. (``pfull``, ``nlat``, ``nlon``) - :type tuple_dims: tuple - + :type tuple_dims: tuple :param idx: index indicating axis with the dimensions to update (e.g. ``idx = 1`` for ``nlat``) - :type idx: int - + :type idx: int :param new_name: new dimension name (e.g. ``latitude``) - :type new_name: str - + :type new_name: str :return: updated dimensions - + :rtype: tuple + :raises KeyError: if the required variables are not found + :raises ValueError: if the required dimensions are not found + :raises AttributeError: if the required attributes are not found .. py:function:: replace_dims(dims, todflag) @@ -390,26 +495,51 @@ Attributes This is designed to work specifically with LegacyCGM.nc files. :param dims: dimensions of the variable - :type dims: str - + :type dims: str :param todflag: indicates whether there exists a ``time_of_day`` dimension - :type todflag: bool - + :type todflag: bool :return: new dimension names for the variable - + :rtype: tuple + :raises KeyError: if the required variables are not found + :raises ValueError: if the required dimensions are not found + :raises AttributeError: if the required attributes are not found .. py:function:: split_files(file_list, split_dim) - Extracts variables in the file along the time dimension, unless - other dimension is specified (lev, lat, or lon). + Extracts variables in the file along the specified dimension. - :param file_list: list of file names - :type split_dim: dimension along which to perform extraction + The function creates a new file with the same name as the input + file, but with the specified dimension sliced to the requested + value or range of values. The new file is saved in the same + directory as the input file. - :returns: new file with sliced dimensions + The function also checks the bounds of the requested dimension + against the actual values in the file. If the requested value + is outside the bounds of the dimension, an error message is + printed and the function exits. + + The function also checks if the requested dimension is valid + (i.e., time, lev, lat, or lon). If the requested dimension is + invalid, an error message is printed and the function exits. + The function also checks if the requested dimension is a + single dimension (i.e., areo). If the requested dimension is + a single dimension, the function removes all single dimensions + from the areo dimension (i.e., scalar_axis) before performing + the extraction. + + :param file_list: list of file names + :type split_dim: dimension along which to perform extraction + :returns: new file with sliced dimensions + :rtype: None + :raises OSError: If the file cannot be removed. + :raises IOError: If the file cannot be moved. + :raises Exception: If the file cannot be opened. + :raises ValueError: If the file cannot be accessed. + :raises TypeError: If the file is not of the correct type. + :raises IndexError: If the file does not have the correct index. .. py:data:: all_args @@ -420,6 +550,14 @@ Attributes +.. py:data:: debug + + + +.. py:data:: exit_code + + + .. py:data:: out_ext diff --git a/docs/source/autoapi/bin/MarsFormat/index.rst b/docs/source/autoapi/bin/MarsFormat/index.rst index 4e12402a..790ca56d 100644 --- a/docs/source/autoapi/bin/MarsFormat/index.rst +++ b/docs/source/autoapi/bin/MarsFormat/index.rst @@ -6,12 +6,12 @@ .. autoapi-nested-parse:: The MarsFormat executable is for converting non-MGCM data, such as that - from EMARS, OpenMARS, PCM, and MarsWRF, into MGCM-like netCDF data + from EMARS, OpenMARS, PCM, and MarsWRF, into MGCM-like netCDF data products. The MGCM is the NASA Ames Mars Global Climate Model developed - and maintained by the Mars Climate Modeling Center (MCMC). The MGCM + and maintained by the Mars Climate Modeling Center (MCMC). The MGCM data repository is available at data.nas.nasa.gov/mcmc. - The executable requires 2 arguments: + The executable requires two arguments: * ``[input_file]`` The file to be transformed * ``[-gcm --gcm_name]`` The GCM from which the file originates @@ -22,17 +22,18 @@ * ``[-ba, --bin_average]`` Bin non-MGCM files like 'average' files * ``[-bd, --bin_diurn]`` Bin non-MGCM files like 'diurn' files - Third-party Requirements: + Third-party requirements: * ``numpy`` * ``netCDF4`` * ``sys`` * ``argparse`` * ``os`` - - List of Functions: - - * download - Queries the requested file from the NAS Data Portal. + * ``re`` + * ``functools`` + * ``traceback`` + * ``xarray`` + * ``amescap`` @@ -45,6 +46,8 @@ Functions .. autoapisummary:: + bin.MarsFormat.debug_wrapper + bin.MarsFormat.get_time_dimension_name bin.MarsFormat.main @@ -55,18 +58,134 @@ Attributes .. autoapisummary:: bin.MarsFormat.args + bin.MarsFormat.debug + bin.MarsFormat.exit_code bin.MarsFormat.parser bin.MarsFormat.path2data bin.MarsFormat.ref_press +.. py:function:: debug_wrapper(func) + + A decorator that wraps a function with error handling + based on the --debug flag. + If the --debug flag is set, it prints the full traceback + of any exception that occurs. Otherwise, it prints a + simplified error message. + + :param func: The function to wrap. + :type func: function + :return: The wrapped function. + :rtype: function + :raises Exception: If an error occurs during the function call. + :raises TypeError: If the function is not callable. + :raises ValueError: If the function is not found. + :raises NameError: If the function is not defined. + :raises AttributeError: If the function does not have the + specified attribute. + :raises ImportError: If the function cannot be imported. + :raises RuntimeError: If the function cannot be run. + :raises KeyError: If the function does not have the + specified key. + :raises IndexError: If the function does not have the + specified index. + :raises IOError: If the function cannot be opened. + :raises OSError: If the function cannot be accessed. + :raises EOFError: If the function cannot be read. + :raises MemoryError: If the function cannot be allocated. + :raises OverflowError: If the function cannot be overflowed. + :raises ZeroDivisionError: If the function cannot be divided by zero. + :raises StopIteration: If the function cannot be stopped. + :raises KeyboardInterrupt: If the function cannot be interrupted. + :raises SystemExit: If the function cannot be exited. + :raises AssertionError: If the function cannot be asserted. + + +.. py:function:: get_time_dimension_name(DS, model) + + Find the time dimension name in the dataset. + + Updates the model object with the correct dimension name. + + :param DS: The xarray Dataset + :type DS: xarray.Dataset + :param model: Model object with dimension information + :type model: object + :return: The actual time dimension name found + :rtype: str + :raises KeyError: If no time dimension is found + :raises ValueError: If the model object is not defined + :raises TypeError: If the dataset is not an xarray Dataset + :raises AttributeError: If the model object does not have the + specified attribute + :raises ImportError: If the xarray module cannot be imported + + .. py:function:: main() + Main processing function for MarsFormat. + + This function processes NetCDF files from various Mars General + Circulation Models (GCMs) + including MarsWRF, OpenMars, PCM, and EMARS, and reformats them for + use in the AmesCAP + framework. + + It performs the following operations: + - Validates the selected GCM type and input files. + - Loads NetCDF files and reads model-specific variable and + dimension mappings. + - Applies model-specific post-processing, including: + - Unstaggering variables (for MarsWRF and EMARS). + - Creating and orienting pressure coordinates (pfull, phalf, + ak, bk). + - Standardizing variable and dimension names. + - Converting longitude ranges to 0-360 degrees east. + - Adding scalar axes where required. + - Handling vertical dimension orientation, especially for + PCM files. + - Optionally performs time binning: + - Daily, average (over N sols), or diurnal binning. + - Ensures correct time units and bin sizes. + - Preserves or corrects vertical orientation after binning. + - Writes processed datasets to new NetCDF files with appropriate + naming conventions. + + Args: + None. Uses global `args` for configuration and file selection. + + Raises: + KeyError: If required dimensions or variables are missing in + the input files. + ValueError: If dimension swapping fails for PCM files. + SystemExit: If no valid GCM type is specified. + + Outputs: + Writes processed NetCDF files to disk, with suffixes indicating + the type of processing + (e.g., _daily, _average, _diurn, _nat). + + Note: + This function assumes the presence of several helper functions + and global variables, + such as `read_variable_dict_amescap_profile`, + `get_time_dimension_name`, `reset_FV3_names`, and color + constants for printing. + + .. py:data:: args +.. py:data:: debug + + + +.. py:data:: exit_code + + + .. py:data:: parser diff --git a/docs/source/autoapi/bin/MarsInterp/index.rst b/docs/source/autoapi/bin/MarsInterp/index.rst index aecc534e..3d338660 100644 --- a/docs/source/autoapi/bin/MarsInterp/index.rst +++ b/docs/source/autoapi/bin/MarsInterp/index.rst @@ -22,7 +22,6 @@ * ``[-ext --extension]`` Custom extension for the new file * ``[-print --print_grid]`` Print the vertical grid to the screen - Third-party Requirements: * ``numpy`` @@ -31,6 +30,11 @@ * ``os`` * ``time`` * ``matplotlib`` + * ``re`` + * ``functools`` + * ``traceback`` + * ``sys`` + * ``amescap`` @@ -43,6 +47,7 @@ Functions .. autoapisummary:: + bin.MarsInterp.debug_wrapper bin.MarsInterp.main @@ -56,6 +61,8 @@ Attributes bin.MarsInterp.M_co2 bin.MarsInterp.R bin.MarsInterp.args + bin.MarsInterp.debug + bin.MarsInterp.exit_code bin.MarsInterp.filepath bin.MarsInterp.fill_value bin.MarsInterp.g @@ -63,8 +70,90 @@ Attributes bin.MarsInterp.rgas +.. py:function:: debug_wrapper(func) + + A decorator that wraps a function with error handling + based on the --debug flag. + If the --debug flag is set, it prints the full traceback + of any exception that occurs. Otherwise, it prints a + simplified error message. + + :param func: The function to wrap. + :type func: function + :return: The wrapped function. + :rtype: function + :raises Exception: If an error occurs during the function call. + :raises TypeError: If the function is not callable. + :raises ValueError: If the function is not found. + :raises NameError: If the function is not defined. + :raises AttributeError: If the function does not have the + specified attribute. + :raises ImportError: If the function cannot be imported. + :raises RuntimeError: If the function cannot be run. + :raises KeyError: If the function does not have the + specified key. + :raises IndexError: If the function does not have the + specified index. + :raises IOError: If the function cannot be opened. + :raises OSError: If the function cannot be accessed. + :raises EOFError: If the function cannot be read. + :raises MemoryError: If the function cannot be allocated. + :raises OverflowError: If the function cannot be overflowed. + :raises ZeroDivisionError: If the function cannot be divided by zero. + :raises StopIteration: If the function cannot be stopped. + :raises KeyboardInterrupt: If the function cannot be interrupted. + :raises SystemExit: If the function cannot be exited. + :raises AssertionError: If the function cannot be asserted. + + .. py:function:: main() + Main function for performing vertical interpolation on Mars + atmospheric model NetCDF files. + + This function processes one or more input NetCDF files, + interpolating variables from their native vertical coordinate + (e.g., model pressure levels) to a user-specified standard vertical + grid (pressure, altitude, or altitude above ground level). + The interpolation type and grid can be customized via command-line + arguments. + + Workflow: + 1. Parses command-line arguments for input files, interpolation + type, custom vertical grid, and other options. + 2. Loads standard vertical grid definitions (pressure, altitude, + or altitude above ground level) or uses a custom grid. + 3. Optionally prints the vertical grid and exits if requested. + 4. For each input file: + - Checks file existence. + - Loads necessary variables (e.g., pk, bk, ps, temperature). + - Computes the 3D vertical coordinate field for + interpolation. + - Creates a new NetCDF output file with updated vertical + dimension. + - Interpolates selected variables to the new vertical grid. + - Copies or interpolates other variables as appropriate. + 5. Handles both regular and diurnal-cycle files, as well as + FV3-tiled and lat/lon grids. + + Command-line arguments (via `args`): + - input_file: List of input NetCDF files to process. + - interp_type: Type of vertical interpolation ('pstd', 'zstd', + or 'zagl'). + - vertical_grid: Custom vertical grid definition (optional). + - print_grid: If True, prints the vertical grid and exits. + - extension: Optional string to append to output filenames. + - include: List of variable names to include in interpolation. + - debug: Enable debug output. + + Notes: + - Requires several helper functions and classes (e.g., + section_content_amescap_profile, find_fixedfile, Dataset, + Ncdf, vinterp). + - Handles both FV3-tiled and regular lat/lon NetCDF files. + - Exits with an error message if required files or variables are + missing. + .. py:data:: Cp :value: 735.0 @@ -85,6 +174,14 @@ Attributes +.. py:data:: debug + + + +.. py:data:: exit_code + + + .. py:data:: filepath diff --git a/docs/source/autoapi/bin/MarsNest/index.rst b/docs/source/autoapi/bin/MarsNest/index.rst new file mode 100644 index 00000000..fc9be0ee --- /dev/null +++ b/docs/source/autoapi/bin/MarsNest/index.rst @@ -0,0 +1,219 @@ +:py:mod:`bin.MarsNest` +====================== + +.. py:module:: bin.MarsNest + +.. autoapi-nested-parse:: + + The MarsNest executable does not accept any input. + + The executable requires the presence of grid_spec files in the current directory. + Test grid_spec files can be generated by running the AmesGCM with 0 days as the input. + + + Third-party Requirements: + + * ``numpy`` + * ``argparse`` + * ``netCDF4`` + * ``matplotlib`` + * ``os`` + + + +Module Contents +--------------- + + +Functions +~~~~~~~~~ + +.. autoapisummary:: + + bin.MarsNest.debug_wrapper + bin.MarsNest.extract_grid + bin.MarsNest.extract_path_basename + bin.MarsNest.find_tile_number + bin.MarsNest.get_corners + bin.MarsNest.get_poly_path + bin.MarsNest.is_child + bin.MarsNest.list_full_paths + bin.MarsNest.main + + + +Attributes +~~~~~~~~~~ + +.. autoapisummary:: + + bin.MarsNest.args + bin.MarsNest.debug + bin.MarsNest.f + bin.MarsNest.lat + bin.MarsNest.lon + bin.MarsNest.map_dir + bin.MarsNest.nskip + bin.MarsNest.parser + bin.MarsNest.parser + bin.MarsNest.thick_line + bin.MarsNest.topo + + +.. py:function:: debug_wrapper(func) + + A decorator that wraps a function with error handling + based on the --debug flag. + + If the --debug flag is set, it prints the full traceback + of any exception that occurs. Otherwise, it prints a + simplified error message. + + :param func: The function to wrap. + :type func: function + :return: The wrapped function. + :rtype: function + :raises Exception: If an error occurs during the function call. + :raises TypeError: If the function is not callable. + :raises ValueError: If the function is not found. + :raises NameError: If the function is not defined. + :raises AttributeError: If the function does not have the + specified attribute. + :raises ImportError: If the function cannot be imported. + :raises RuntimeError: If the function cannot be run. + :raises KeyError: If the function does not have the + specified key. + :raises IndexError: If the function does not have the + specified index. + :raises IOError: If the function cannot be opened. + :raises OSError: If the function cannot be accessed. + :raises EOFError: If the function cannot be read. + :raises MemoryError: If the function cannot be allocated. + :raises OverflowError: If the function cannot be overflowed. + :raises ZeroDivisionError: If the function cannot be divided by zero. + :raises StopIteration: If the function cannot be stopped. + :raises KeyboardInterrupt: If the function cannot be interrupted. + :raises SystemExit: If the function cannot be exited. + :raises AssertionError: If the function cannot be asserted. + + +.. py:function:: extract_grid(lon, lat, varIN, lon_target, lat_target, window_lon, window_lat) + + This function extract a subsample of a larger grid variable. It deals with wrapping around the longitudes on the edges of the grid. + Args: + lon: 1D array of longitude -180/+180 + lat: 1D array of latitudes -90/+90 + varIN: 2D array var2IN(lat,lon) + lon_target: float, center longitude of subgrid (e.g =-136.) + lat_target: float, center latitude of subgrid (e.g =40.) + window_lon: width of window in degree (e.g =10.) + window_lat: heigth of window in degree (e.g. =5.) + Returns: + lon_window: 1D array of truncated longitude (may not be increasing) + lat_window: 1D array of truncated latitudes (may not be increasing) + varIN_out: 2D array of subsample of varIN + + +.. py:function:: extract_path_basename(filename) + + Return the path and basename of a file. If only the filename is provided, assume it is the current directory + Args: + filename: e.g. 'XXXXX.fixed.nc', './history/XXXXX.fixed.nc' or '/home/user/history/XXXXX.fixed.nc' + Returns: + filepath : '/home/user/history/XXXXX.fixed.nc' in all the cases above + basename: XXXXX.fixed.nc in all the cases above + + ***NOTE*** + This routine does not check for file existence and only operates on the provided input string. + + +.. py:function:: find_tile_number(filename) + + Return the line number given a file name. + (Spit grid_spec.nest03.tile8.nc as "grid_spec.nest03" ".tile" "8.nc" and remove ".nc" + + +.. py:function:: get_corners(poly_path) + + Return the North West and South East corner of a gridscspec file. + Args: + gridspec_filename: full path to a gridscpec file + Return: + NWc,SEc =(lat,lon) of the two corners + + +.. py:function:: get_poly_path(gridspec_filename) + + Return the North West and South East corner of a gridscspec file. + Args: + gridspec_filename: full path to a gridscpec file + Return: + poly = a numpy array size (N,(lat,lon) )with the corners + + +.. py:function:: is_child(poly_path_child, poly_path_parent) + + Test if a grid is nested in another grid. + Args: + poly_path_child: a polygon object wit the 4 corners of the child tile + poly_path_parent: a polygon object wit the 4 corners of the parent tile + Returns: + bool: True if all the corners of the child grids are contained within the path defined by the paretn grid + + + +.. py:function:: list_full_paths(input_path) + + List all files in directory. + + +.. py:function:: main() + + +.. py:data:: args + + + +.. py:data:: debug + + + +.. py:data:: f + + + +.. py:data:: lat + + + +.. py:data:: lon + + + +.. py:data:: map_dir + :value: '/nobackup/rurata/FMS_MARS_data/mars_topo.nc' + + + +.. py:data:: nskip + :value: 2 + + + +.. py:data:: parser + + + +.. py:data:: parser + + + +.. py:data:: thick_line + :value: 0.5 + + + +.. py:data:: topo + + + diff --git a/docs/source/autoapi/bin/MarsPlot/index.rst b/docs/source/autoapi/bin/MarsPlot/index.rst index fe84d363..e524d7ae 100644 --- a/docs/source/autoapi/bin/MarsPlot/index.rst +++ b/docs/source/autoapi/bin/MarsPlot/index.rst @@ -24,6 +24,7 @@ * ``warnings`` * ``subprocess`` * ``matplotlib`` + * ``pypdf`` @@ -56,6 +57,7 @@ Functions bin.MarsPlot.clean_comma_whitespace bin.MarsPlot.create_exec bin.MarsPlot.create_name + bin.MarsPlot.debug_wrapper bin.MarsPlot.fig_layout bin.MarsPlot.filter_input bin.MarsPlot.format_lon_lat @@ -69,7 +71,6 @@ Functions bin.MarsPlot.get_overwrite_dim_2D bin.MarsPlot.get_time_index bin.MarsPlot.get_tod_index - bin.MarsPlot.give_permission bin.MarsPlot.main bin.MarsPlot.make_template bin.MarsPlot.mean_func @@ -93,7 +94,9 @@ Attributes bin.MarsPlot.add_sol_time_axis bin.MarsPlot.args bin.MarsPlot.current_version + bin.MarsPlot.debug bin.MarsPlot.degr + bin.MarsPlot.exit_code bin.MarsPlot.include_NaNs bin.MarsPlot.lon_coord_type bin.MarsPlot.namespace @@ -189,6 +192,91 @@ Attributes Bases: :py:obj:`object` + Fig_1D is a parent class for generating and handling 1D plots of + Mars atmospheric data. + + Attributes: + title : str + Title of the plot. + legend : str + Legend label for the plot. + varfull : str + Full variable specification, including file and variable + name. + t : str or float + Time axis or identifier for the varying dimension. + lat : float or str + Latitude value or identifier. + lon : float or str + Longitude value or identifier. + lev : float or str + Vertical level value or identifier. + ftod : float or str + Time of day requested. + hour : float or str + Hour of day, used for diurnal plots. + doPlot : bool + Whether to generate the plot. + plot_type : str + Type of 1D plot (e.g., "1D_time", "1D_lat"). + sol_array : str + Sol array extracted from varfull. + filetype : str + File type extracted from varfull. + var : str + Variable name extracted from varfull. + simuID : str + Simulation ID extracted from varfull. + nPan : int + Number of panels in the plot. + subID : int + Subplot ID. + addLine : bool + Whether to add a line to an existing plot. + layout : list or None + Page layout for multipanel plots. + fdim_txt : str + Annotation for free dimensions. + success : bool + Indicates if the plot was successfully created. + vert_unit : str + Vertical unit, either "m" or "Pa". + Dlim : list or None + Dimension limits for the axis. + Vlim : list or None + Variable limits for the axis. + axis_opt1 : str + Line style or axis option. + axis_opt2 : str + Additional axis option (optional). + + Methods: + make_template(): + Writes a template for the plot configuration to a file. + read_template(): + Reads plot configuration from a template file. + get_plot_type(): + Determines the type of 1D plot to create based on which + dimension is set to "AXIS" or -88888. + data_loader_1D(varfull, plot_type): + Loads 1D data for plotting, handling variable expressions + and dimension overwrites. + read_NCDF_1D(var_name, file_type, simuID, sol_array, plot_type, + t_req, lat_req, lon_req, lev_req, ftod_req): + Reads and processes 1D data from a NetCDF file for the + specified variable and dimensions. + exception_handler(e, ax): + Handles exceptions during plotting, displaying an error + message on the plot. + fig_init(): + Initializes the figure and subplot for plotting. + fig_save(): + Saves the generated figure to disk. + do_plot(): + Main method to generate the 1D plot, handling all plotting + logic and exceptions. + + .. py:method:: data_loader_1D(varfull, plot_type) @@ -213,7 +301,6 @@ Attributes :return: type of 1D plot to create (1D_time, 1D_lat, etc.) - .. py:method:: make_template() @@ -224,43 +311,32 @@ Attributes plot. :param var_name: variable name (e.g., ``temp``) - :type var_name: str - + :type var_name: str :param file_type: MGCM output file type. Must be ``fixed`` or ``average`` - :type file_type: str - + :type file_type: str :param simuID: number identifier for netCDF file directory - :type simuID: str - + :type simuID: str :param sol_array: sol if different from default (e.g., ``02400``) - :type sol_array: str - + :type sol_array: str :param plot_type: ``1D_lon``, ``1D_lat``, ``1D_lev``, or ``1D_time`` - :type plot_type: str - + :type plot_type: str :param t_req: Ls requested - :type t_req: str - + :type t_req: str :param lat_req: lat requested - :type lat_req: str - + :type lat_req: str :param lon_req: lon requested - :type lon_req: str - + :type lon_req: str :param lev_req: level [Pa/m] requested - :type lev_req: str - + :type lev_req: str :param ftod_req: time of day requested - :type ftod_req: str - + :type ftod_req: str :return: (dim_array) the axis (e.g., an array of longitudes), (var_array) the variable extracted - .. py:method:: read_template() @@ -270,6 +346,32 @@ Attributes Bases: :py:obj:`object` + Base class for 2D figures. This class is not intended to be + instantiated directly. Instead, it is used as a base class for + specific 2D figure classes (e.g., ``Fig_2D_lon_lat``, ``Fig_2D_time_lat``, + ``Fig_2D_lat_lev``, etc.). It provides common attributes and methods + for all 2D figures, such as the variable name, file type, simulation + ID, and plotting options. The class also includes methods for + creating a template for the figure, reading the template from a + file, and loading data for 2D plots. The class is designed to be + extended by subclasses that implement specific plotting + functionality for different types of 2D figures. + + :param varfull: full variable name (e.g., ``fileYYY.XXX``) + :type varfull: str + :param doPlot: whether to plot the figure (default: ``False``) + :type doPlot: bool + :param varfull2: second variable name (default: ``None``) + :type varfull2: str + :return: None + :rtype: None + :raises ValueError: If the input varfull is not a valid type + for variable name. + :raises TypeError: If the input doPlot is not a valid type + for plotting. + :raises Exception: If the input varfull2 is not a valid type + for variable name. + .. py:method:: data_loader_2D(varfull, plot_type) @@ -315,11 +417,65 @@ Attributes Bases: :py:obj:`Fig_2D` + A subclass of Fig_2D for generating 2D plots with latitude and + vertical level (pressure or altitude) axes. + + This class customizes the plotting template and plotting logic for + visualizing data as a function of latitude and vertical level. + It supports filled contour plots for a primary variable, and + optionally overlays solid contour lines for a secondary variable. + + Methods: + make_template(): + Sets up the plot template with appropriate titles and axis + labels for latitude vs. level plots. + + do_plot(): + Loads data, creates a filled contour plot of the primary + variable, optionally overlays contours of a secondary + variable, configures axis scaling and formatting (including + logarithmic pressure axis if needed), sets axis limits and + tick formatting, handles exceptions, and saves the resulting + figure. + + Attributes (inherited and/or used): + varfull : str + Name of the primary variable to plot. + varfull2 : str or None + Name of the secondary variable to overlay as contours, if + any. + plot_type : str + Type of plot or data selection. + vert_unit : str + Unit for the vertical axis ("Pa" for pressure, otherwise + altitude in meters). + Xlim : tuple or None + Limits for the x-axis (latitude). + Ylim : tuple or None + Limits for the y-axis (level). + contour2 : list or None + Contour levels for the secondary variable. + nPan : int + Number of panels in the plot (affects tick label size). + success : bool + Indicates if the plot was successfully created. + .. py:method:: data_loader_2D(varfull, plot_type) .. py:method:: do_plot() + Generates a 2D latitude-level plot for the specified variable(s). + This method initializes the figure, loads the required data, and creates a filled contour plot + of the primary variable. If a secondary variable is specified, it overlays solid contours for + that variable. The y-axis is set to logarithmic scale and inverted if the vertical unit is pressure. + Axis limits, labels, and tick formatting are applied as specified by the instance attributes. + The plot title is generated based on the variable information. Handles exceptions during plotting + and saves the resulting figure. + + Raises: + Exception: Any exception encountered during plotting is handled and logged. + .. py:method:: exception_handler(e, ax) @@ -338,6 +494,17 @@ Attributes .. py:method:: make_template() + Creates and configures a plot template for a 2D latitude versus level plot. + This method calls the parent class's `make_template` method with predefined + titles and axis labels suitable for a plot displaying latitude against atmospheric + level data. + The plot is labeled as "Plot 2D lat X lev" with the following axis labels: + - X-axis: "Ls 0-360 " + - Y-axis: "Lon +/-180" + - Additional axes: "Lat", "Level[Pa/m]" + Returns: + None + .. py:method:: make_title(var_info, xlabel, ylabel) @@ -363,11 +530,91 @@ Attributes Bases: :py:obj:`Fig_2D` + Fig_2D_lon_lat is a class for creating 2D longitude-latitude plots. + + Fig_2D_lon_lat is a subclass of Fig_2D designed for generating 2D + plots of longitude versus latitude, primarily for visualizing Mars + climate data. It provides methods for figure creation, data loading, + plotting, and overlaying topography contours, with support for + various map projections and customization options. + + Attributes: + varfull (str): Full variable name (e.g., "fileYYY.XXX") to plot. + doPlot (bool): Whether to plot the figure (default: False). + varfull2 (str, optional): Second variable name for overlaying + contours (default: None). + plot_type (str): Type of plot (default: "2D_lon_lat"). + fdim1 (str, optional): First free dimension (default: None). + fdim2 (str, optional): Second free dimension (default: None). + ftod (str, optional): Time of day (default: None). + axis_opt1 (str, optional): First axis option, e.g., colormap + (default: None). + axis_opt2 (str, optional): Second axis option (default: None). + axis_opt3 (str, optional): Projection type (e.g., "cart", + "robin", "moll", "Npole", "Spole", "ortho"). + Xlim (tuple, optional): Longitude axis limits. + Ylim (tuple, optional): Latitude axis limits. + range (bool, optional): Whether to use a specified range for + color levels. + contour2 (float or list, optional): Contour levels for the + second variable. + title (str, optional): Custom plot title. + nPan (int): Number of panels (for multi-panel plots). + fdim_txt (str): Text describing free dimensions. + success (bool): Status flag indicating if plotting succeeded. + + Methods: + make_template(): + Sets up the plot template with appropriate axis labels and + titles. + + get_topo_2D(varfull, plot_type): + Loads and returns topography data (zsurf) for overlaying as + contours, matching the simulation and file type of the main variable. + + do_plot(): + Main plotting routine. Loads data, applies projection, + overlays topography and optional second variable contours, + customizes axes, and saves the figure. Handles both + standard and special map projections (cartesian, Robinson, + Mollweide, polar, orthographic). + + Usage: + This class is intended to be used within the MarsPlot software + for visualizing Mars climate model outputs as longitude-latitude + maps, with optional overlays and advanced projection support. + .. py:method:: data_loader_2D(varfull, plot_type) .. py:method:: do_plot() + Generate a 2D longitude-latitude plot with various projection options and optional overlays. + + This method creates a 2D plot of a variable (and optionally a second variable as contours) + on a longitude-latitude grid. It supports multiple map projections, including cartesian, + Robinson, Mollweide, and azimuthal (north pole, south pole, orthographic) projections. + Topography contours can be added if available. The method handles axis formatting, + colorbars, titles, and annotation of meridians and parallels. + + The plotting behavior is controlled by instance attributes such as: + - self.varfull: Main variable to plot. + - self.varfull2: Optional second variable for contour overlay. + - self.plot_type: Type of plot to generate. + - self.axis_opt1: Colormap or colormap option. + - self.axis_opt3: Projection type. + - self.contour2: Contour levels for the second variable. + - self.Xlim, self.Ylim: Axis limits for cartesian projection. + - self.range: Whether to use a specific range for color levels. + - self.title: Custom plot title. + - self.fdim_txt: Additional dimension text for the title. + - self.nPan: Panel index for multi-panel plots. + + The method handles exceptions and saves the figure upon completion. + + Raises: + Exception: Any error encountered during plotting is handled and reported. + .. py:method:: exception_handler(e, ax) @@ -386,29 +633,35 @@ Attributes This function returns the longitude, latitude, and topography to overlay as contours in a ``2D_lon_lat`` plot. Because the main variable requested may be complex - (e.g., ``[00668.atmos_average_psdt2.temp]/1000.``), we will - ensure to load the matching topography (here ``00668.fixed.nc`` + (e.g., ``[01336.atmos_average_psdt2.temp]/1000.``), we will + ensure to load the matching topography (here ``01336.fixed.nc`` from the 2nd simulation). This function essentially does a simple task in a complicated way. Note that a great deal of the code is borrowed from the ``data_loader_2D()`` function. :param varfull: variable input to main_variable in Custom.in (e.g., ``03340.atmos_average.ucomp``) - :type varfull: str - + :type varfull: str :param plot_type: plot type (e.g., ``Plot 2D lon X time``) - :type plot_type: str - + :type plot_type: str :return: topography or ``None`` if no matching ``fixed`` file - .. py:method:: make_colorbar(levs) .. py:method:: make_template() + Creates and configures a plot template for 2D longitude vs latitude data. + This method calls the parent class's `make_template` method with predefined + parameters to set up the plot title and axis labels specific to a 2D longitude-latitude plot. + The template includes: + - Title: "Plot 2D lon X lat" + - X-axis label: "Ls 0-360" + - Y-axis label: "Level Pa/m" + - Additional axis labels: "Lon" (longitude), "Lat" (latitude) + .. py:method:: make_title(var_info, xlabel, ylabel) @@ -434,13 +687,48 @@ Attributes Bases: :py:obj:`Fig_2D` + A subclass of Fig_2D for generating 2D plots with longitude and + vertical level (pressure or altitude) axes. + + This class customizes the template and plotting routines to + visualize data as a function of longitude and vertical level. + It supports plotting filled contours for a primary variable and + optional solid contours for a secondary variable. + The vertical axis can be displayed in pressure (Pa, logarithmic + scale) or altitude (m). + + Methods: + make_template(): + Sets up the plot template with appropriate titles and axis + labels for longitude vs. level plots. + + do_plot(): + Loads data, applies longitude shifting, creates filled and + optional solid contour plots, + configures axis scales and labels, and handles exceptions + during plotting. + .. py:method:: data_loader_2D(varfull, plot_type) .. py:method:: do_plot() - Create figure + Generates a 2D plot of a variable as a function of longitude + and vertical level (pressure or altitude). + + This method initializes the figure, loads the required data, + applies longitude shifting, and creates filled and/or solid + contour plots. + It handles plotting of a secondary variable if specified, sets + axis scales and labels based on the vertical coordinate unit, + applies axis limits if provided, customizes tick formatting and + font sizes, and manages exceptions during plotting. + The resulting figure is saved to file. + + Raises: + Exception: If any error occurs during the plotting process, + it is handled and logged by the exception handler. .. py:method:: exception_handler(e, ax) @@ -460,8 +748,17 @@ Attributes .. py:method:: make_template() - Calls method from parent class + Creates and configures a plot template for 2D lon x lev data. + + This method sets up the plot with predefined titles and axis + labels: + - Title: "Plot 2D lon X lev" + - X-axis: "Ls 0-360" + - Y-axis: "Latitude" + - Additional labels: "Lon +/-180" and "Level[Pa/m]" + Overrides the base class method to provide specific + configuration for this plot type. .. py:method:: make_title(var_info, xlabel, ylabel) @@ -488,6 +785,25 @@ Attributes Bases: :py:obj:`Fig_2D` + A specialized 2D plotting class for visualizing data as a function + of longitude and time (Ls). + + Inherits from: Fig_2D + + Methods: + make_template(): + Sets up the plot template with appropriate titles and axis + labels for longitude vs. time plots. + + do_plot(): + Generates a 2D plot with longitude on the x-axis and solar + longitude (Ls) on the y-axis. + Loads and processes data, applies shifting if necessary, + and creates filled and/or solid contours. + Handles axis formatting, tick labeling (including optional + sol time annotation), and plot saving. + Catches and handles exceptions during plotting. + .. py:method:: data_loader_2D(varfull, plot_type) @@ -536,11 +852,58 @@ Attributes Bases: :py:obj:`Fig_2D` + A 2D plotting class for visualizing data as a function of time (Ls) + and latitude. Inherits from: Fig_2D + + Methods: + make_template(): + Sets up the plot template with appropriate titles and axis + labels for a 2D time vs latitude plot. + do_plot(): + Loads 2D data (time and latitude), creates a filled contour + plot of the primary variable, and optionally overlays a + solid contour of a secondary variable. + Formats axes, customizes tick labels to show both Ls and + sol time (if enabled), and applies axis limits if specified. + Handles exceptions during plotting and saves the resulting + figure. + + Attributes (inherited and used): + varfull : str + Name of the primary variable to plot. + varfull2 : str or None + Name of the secondary variable to overlay as contours + (optional). + plot_type : str + Type of plot/data to load. + Xlim : tuple or None + Limits for the x-axis (sol time). + Ylim : tuple or None + Limits for the y-axis (latitude). + contour2 : list or None + Contour levels for the secondary variable. + nPan : int + Number of panels (used for label sizing). + success : bool + Indicates if the plot was successfully created. + .. py:method:: data_loader_2D(varfull, plot_type) .. py:method:: do_plot() + Generates a 2D time-latitude plot for the specified variable(s). + This method initializes the figure, loads the required 2D data arrays (time and latitude), + and creates a filled contour plot of the primary variable. If a secondary variable is specified, + it overlays solid contours for that variable. The method also formats the axes, including + custom tick labels for solar longitude (Ls) and optionally sol time, and applies axis limits + if specified. Additional plot formatting such as tick intervals and font sizes are set. + The plot is saved at the end of the method. Any exceptions encountered during plotting + are handled and reported. + + Raises: + Exception: If any error occurs during the plotting process, it is handled and reported. + .. py:method:: exception_handler(e, ax) @@ -559,6 +922,14 @@ Attributes .. py:method:: make_template() + Creates and configures a plot template for a 2D time versus latitude figure. + This method calls the superclass's `make_template` method with predefined + titles and axis labels suitable for a plot displaying data across longitude, + level, solar longitude (Ls), and latitude. + + Returns: + None + .. py:method:: make_title(var_info, xlabel, ylabel) @@ -584,11 +955,77 @@ Attributes Bases: :py:obj:`Fig_2D` + A specialized 2D plotting class for visualizing data as a function + of time (Ls) and vertical level (pressure or altitude). + + Inherits from: Fig_2D + + Methods: + make_template(): + Sets up the plot template with appropriate axis labels and + titles for 2D time vs. level plots. + + do_plot(): + Loads data and generates a filled contour plot of the + primary variable as a function of solar longitude (Ls) and + vertical level. Optionally overlays a solid contour of a + secondary variable. + Handles axis formatting, tick labeling (including optional + sol time axis), and y-axis scaling (logarithmic for + pressure). Sets plot titles and saves the figure. Catches + and handles exceptions during plotting. + + Attributes (inherited and/or used): + varfull : str + Name of the primary variable to plot. + varfull2 : str or None + Name of the secondary variable to overlay as contours + (optional). + plot_type : str + Type of plot/data selection. + Xlim : tuple or None + Limits for the x-axis (solar day). + Ylim : tuple or None + Limits for the y-axis (vertical level). + vert_unit : str + Unit for the vertical axis ("Pa" for pressure or other for + altitude). + nPan : int + Number of panels/subplots (affects label size). + contour2 : list or None + Contour levels for the secondary variable. + success : bool + Indicates if the plot was successfully generated. + + .. py:method:: data_loader_2D(varfull, plot_type) .. py:method:: do_plot() + Generates a 2D time-level plot for Mars atmospheric data. + + This method initializes the figure, loads the required data, and creates a filled contour plot + of the primary variable over solar longitude (Ls) and pressure or altitude. If a secondary variable + is specified, it overlays solid contours for that variable. The method also formats axes, applies + custom tick labels (optionally including sol time), and adjusts axis scales and labels based on + the vertical unit (pressure or altitude). The plot is titled and saved to file. + Handles exceptions by invoking a custom exception handler and always attempts to save the figure. + + Attributes used: + varfull (str): Name of the primary variable to plot. + plot_type (str): Type of plot/data to load. + varfull2 (str, optional): Name of the secondary variable for contour overlay. + contour2 (list, optional): Contour levels for the secondary variable. + Xlim (tuple, optional): Limits for the x-axis (solar day). + Ylim (tuple, optional): Limits for the y-axis (pressure or altitude). + vert_unit (str): Vertical axis unit, either "Pa" for pressure or other for altitude. + nPan (int): Number of panels (affects label size). + success (bool): Set to True if plotting succeeds. + + Raises: + Handles all exceptions internally and logs them via a custom handler. + .. py:method:: exception_handler(e, ax) @@ -607,6 +1044,14 @@ Attributes .. py:method:: make_template() + Creates and configures a plot template for 2D time versus level visualization. + This method calls the superclass's `make_template` method with predefined + titles and axis labels suitable for plotting data with latitude, longitude, + solar longitude (Ls), and atmospheric level (in Pa/m). + + Returns: + None + .. py:method:: make_title(var_info, xlabel, ylabel) @@ -629,14 +1074,14 @@ Attributes .. py:function:: MY_func(Ls_cont) - Returns the Mars Year + Returns the Mars Year. :param Ls_cont: solar longitude (``areo``; continuous) - :type Ls_cont: array [areo] - + :type Ls_cont: array [areo] :return: the Mars year - :rtype: int - + :rtype: int + :raises ValueError: If the input Ls_cont is not a valid type for + year calculation. .. py:function:: clean_comma_whitespace(raw_input) @@ -645,12 +1090,10 @@ Attributes :param raw_input: dimensions specified by user input to Variable (e.g., ``lat=3. , lon=2 , lev = 10.``) - :type raw_input: str - + :type raw_input: str :return: raw_input without whitespaces (e.g., ``lat=3.,lon=2,lev=10.``) - :rtype: str - + :rtype: str .. py:function:: create_exec(raw_input, varfull_list) @@ -662,12 +1105,52 @@ Attributes :param root_name: path + default name for the file type (e.g., ``/path/custom.in`` or ``/path/figure.png``) - :type root_name: str - + :type root_name: str :return: the modified name if the file already exists (e.g., ``/path/custom_01.in`` or ``/path/figure_01.png``) - :rtype: str - + :rtype: str + :raises ValueError: If the input root_name is not a valid type + for file name. + :raises TypeError: If the input root_name is not a valid type + for file name. + :raises Exception: If the file name creation fails for any + reason. + + +.. py:function:: debug_wrapper(func) + + A decorator that wraps a function with error handling + based on the --debug flag. + If the --debug flag is set, it prints the full traceback + of any exception that occurs. Otherwise, it prints a + simplified error message. + + :param func: The function to wrap. + :type func: function + :return: The wrapped function. + :rtype: function + :raises Exception: If an error occurs during the function call. + :raises TypeError: If the function is not callable. + :raises ValueError: If the function is not found. + :raises NameError: If the function is not defined. + :raises AttributeError: If the function does not have the + specified attribute. + :raises ImportError: If the function cannot be imported. + :raises RuntimeError: If the function cannot be run. + :raises KeyError: If the function does not have the + specified key. + :raises IndexError: If the function does not have the + specified index. + :raises IOError: If the function cannot be opened. + :raises OSError: If the function cannot be accessed. + :raises EOFError: If the function cannot be read. + :raises MemoryError: If the function cannot be allocated. + :raises OverflowError: If the function cannot be overflowed. + :raises ZeroDivisionError: If the function cannot be divided by zero. + :raises StopIteration: If the function cannot be stopped. + :raises KeyboardInterrupt: If the function cannot be interrupted. + :raises SystemExit: If the function cannot be exited. + :raises AssertionError: If the function cannot be asserted. .. py:function:: fig_layout(subID, nPan, vertical_page=False) @@ -675,36 +1158,44 @@ Attributes Return figure layout. :param subID: current subplot number - :type subID: int - + :type subID: int :param nPan: number of panels desired on page (max = 64, 8x8) - :type nPan: int - + :type nPan: int :param vertical_page: reverse the tuple for portrait format if ``True`` - :type vertical_page: bool - + :type vertical_page: bool :return: plot layout (e.g., ``plt.subplot(nrows = out[0], ncols = out[1], plot_number = out[2])``) - :rtype: tuple - + :rtype: tuple + :raises ValueError: If the input subID is not a valid type for + subplot number. + :raises TypeError: If the input nPan is not a valid type for + subplot number. + :raises Exception: If the input vertical_page is not a valid type + for subplot number. + :raises Exception: If the figure layout calculation fails for any + reason. .. py:function:: filter_input(txt, typeIn='char') - Read template for the type of data expected + Read template for the type of data expected. + + Returns value to ``rT()``. :param txt: text input into ``Custom.in`` to the right of an equal sign - :type txt: str - + :type txt: str :param typeIn: type of data expected: ``char``, ``float``, ``int``, ``bool``, defaults to ``char`` - :type typeIn: str, optional - + :type typeIn: str, optional :return: text input reformatted to ``[val1, val2]`` - :rtype: float or array - + :rtype: float or array + :raises ValueError: If the input txt is not a valid type for + filtering. + :raises TypeError: If the input typeIn is not a valid type for + filtering. + :raises Exception: If the filtering operation fails for any reason. .. py:function:: format_lon_lat(lon_lat, type) @@ -713,14 +1204,16 @@ Attributes 45°E) :param lon_lat: latitude or longitude (+180/-180) - :type lon_lat: float - + :type lon_lat: float :param type: ``lat`` or ``lon`` - :type type: str - + :type type: str :return: formatted label - :rtype: str - + :rtype: str + :raises ValueError: If the input lon_lat is not a valid type for + latitude or longitude. + :raises TypeError: If the input type is not a valid type for + latitude or longitude. + :raises Exception: If the formatting fails for any reason. .. py:function:: get_Ncdf_num() @@ -729,8 +1222,9 @@ Attributes Requires at least one ``fixed`` file in the directory. :return: a sorted array of sols - :rtype: array - + :rtype: array + :raises ValueError: If the input input_paths is not a valid type + for file name. .. py:function:: get_figure_header(line_txt) @@ -739,135 +1233,154 @@ Attributes :param line_txt: template header from Custom.in (e.g., ``<<<<<<<<<| Plot 2D lon X lat = True |>>>>>>>>``) - :type line_txt: str - + :type line_txt: str :return: (figtype) figure type (e.g., ``Plot 2D lon X lat``) - :rtype: str - + :rtype: str :return: (boolPlot) whether to plot (``True``) or skip (``False``) figure - :rtype: bool - + :rtype: bool + :raises ValueError: If the input line_txt is not a valid type for + figure header. + :raises TypeError: If the input line_txt is not a valid type for + figure header. + :raises Exception: If the figure header parsing fails for any + reason. .. py:function:: get_lat_index(lat_query, lats) - Returns the indices that will extract data from the netCDF file - according to a range of *latitudes*. + Returns the indices for a range of latitudes in a file. :param lat_query: requested latitudes (-90/+90) - :type lat_query: list - + :type lat_query: list :param lats: latitude - :type lats: array [lat] - + :type lats: array [lat] :return: 1d array of file indices - :rtype: text descriptor for the extracted longitudes + :rtype: text descriptor for the extracted longitudes + :rtype: str + :raises ValueError: If the input lat_query is not a valid type for + latitude calculation. .. note::T The keyword ``all`` passed as ``-99999`` by the ``rt()`` function - .. py:function:: get_level_index(level_query, levs) - Returns the indices that will extract data from the netCDF file - according to a range of *pressures* (resp. depth for ``zgrid``). + Returns the indices for a range of pressures in a file. :param level_query: requested pressure [Pa] (depth [m]) - :type level_query: float - + :type level_query: float :param levs: levels (in the native coordinates) - :type levs: array [lev] - + :type levs: array [lev] :return: file indices - :rtype: array - + :rtype: array :return: descriptor for the extracted pressure (depth) - :rtype: str + :rtype: str + :raises ValueError: If the input level_query is not a valid type for + level calculation. .. note:: The keyword ``all`` is passed as ``-99999`` by the ``rT()`` functions - .. py:function:: get_list_varfull(raw_input) Return requested variable from a complex ``varfull`` object with ``[]``. :param raw_input: complex user input to Variable (e.g., ``2*[atmos_average.temp]+[atmos_average2.ucomp]*1000``) - :type raw_input: str - + :type raw_input: str :return: list required variables (e.g., [``atmos_average.temp``, ``atmos_average2.ucomp``]) - :rtype: str - + :rtype: str + :raises ValueError: If the input raw_input is not a valid type for + variable extraction. .. py:function:: get_lon_index(lon_query_180, lons) - Returns the indices that will extract data from the netCDF file - according to a range of *longitudes*. + Returns the indices for a range of longitudes in a file. :param lon_query_180: longitudes in -180/180: value, ``[min, max]``, or `None` - :type lon_query_180: list - + :type lon_query_180: list :param lons: longitude in 0-360 - :type lons: array [lon] - + :type lons: array [lon] :return: 1D array of file indices - :rtype: array - + :rtype: array :return: text descriptor for the extracted longitudes - :rtype: str + :rtype: str + :raises ValueError: If the input lon_query_180 is not a valid type + for longitude calculation. .. note:: The keyword ``all`` passed as ``-99999`` by the rT() functions - .. py:function:: get_overwrite_dim_1D(varfull_bracket, t_in, lat_in, lon_in, lev_in, ftod_in) - Return new dimensions that will overwrite default dimensions for a - varfull object with ``{}`` for a 1D plot. + 1D plot: overwrite dimensions in ``varfull`` object with ``{}``. + + (e.g., ``atmos_average.temp{lev=10;ls=350;lon=155;lat=25}``) + This function is used to overwrite the default dimensions in a + ``varfull`` object with ``{}`` (e.g., ``atmos_average.temp{lev=10; + ls=350;lon=155;lat=25}``) for a 1D plot. The function will return + the new dimensions that will overwrite the default dimensions for + the ``varfull`` object. The function will also return the required + file and variable (e.g., ``atmos_average.temp``) and the X and Y + axis dimensions for the plot. :param varfull_bracket: a ``varfull`` object with ``{}`` (e.g., ``atmos_average.temp{lev=10;ls=350;lon=155;lat=25}``) - :type varfull_bracket: str - + :type varfull_bracket: str :param t_in: self.t variable - :type t_in: array [time] - + :type t_in: array [time] :param lat_in: self.lat variable - :type lat_in: array [lat] - + :type lat_in: array [lat] :param lon_in: self.lon variable - :type lon_in: array [lon] - + :type lon_in: array [lon] :param lev_in: self.lev variable - :type lev_in: array [lev] - + :type lev_in: array [lev] :param ftod_in: self.ftod variable - :type ftod_in: array [tod] - + :type ftod_in: array [tod] :return: ``varfull`` object without brackets (e.g., ``atmos_average.temp``); - :return: (t_out) dimension to update; - :return: (lat_out) dimension to update; - :return: (lon_out) dimension to update; - :return: (lev_out) dimension to update; - :return: (ftod_out) dimension to update; - + :return: (t_out) dimension to update; + :return: (lat_out) dimension to update; + :return: (lon_out) dimension to update; + :return: (lev_out) dimension to update; + :return: (ftod_out) dimension to update; + :rtype: str + :raises ValueError: If the input varfull_bracket is not a valid + type for variable extraction. + :raises TypeError: If the input t_in, lat_in, lon_in, lev_in, + ftod_in are not valid types for variable extraction. + :raises Exception: If the variable extraction fails for any reason. + + .. note:: This function is used for 1D plots only. The function + will return the new dimensions that will overwrite the default + dimensions for the ``varfull`` object. The function will also + return the required file and variable (e.g., + ``atmos_average.temp``) and the X and Y axis dimensions for the + plot. .. py:function:: get_overwrite_dim_2D(varfull_bracket, plot_type, fdim1, fdim2, ftod) - Return new dimensions that will overwrite default dimensions for a - varfull object with ``{}`` on a 2D plot. + 2D plot: overwrite dimensions in ``varfull`` object with ``{}``. + + (e.g., ``atmos_average.temp{lev=10;ls=350;lon=155;lat=25}``) + + This function is used to overwrite the default dimensions in a + ``varfull`` object with ``{}`` (e.g., ``atmos_average.temp{lev=10; + ls=350;lon=155;lat=25}``) for a 2D plot. The function will return + the new dimensions that will overwrite the default dimensions for + the ``varfull`` object. The function will also return the required + file and variable (e.g., ``atmos_average.temp``) and the X and Y + axis dimensions for the plot. ``2D_lon_lat: fdim1 = ls, fdim2 = lev`` ``2D_lat_lev: fdim1 = ls, fdim2 = lon`` @@ -878,85 +1391,95 @@ Attributes :param varfull_bracket: a ``varfull`` object with ``{}`` (e.g., ``atmos_average.temp{lev=10;ls=350;lon=155;lat=25}``) - :type varfull_bracket: str - + :type varfull_bracket: str :param plot_type: the type of the plot template - :type plot_type: str - + :type plot_type: str :param fdim1: X axis dimension for plot - :type fdim1: str - + :type fdim1: str :param fdim2: Y axis dimension for plot - :type fdim2: str - + :type fdim2: str :return: (varfull) required file and variable (e.g., ``atmos_average.temp``); (fdim_out1) X axis dimension for plot; (fdim_out2) Y axis dimension for plot; (ftod_out) if X or Y axis dimension is time of day - + :rtype: str + :raises ValueError: If the input varfull_bracket is not a valid + type for variable extraction. + :raises TypeError: If the input plot_type is not a valid type for + variable extraction. + :raises Exception: If the variable extraction fails for any reason. .. py:function:: get_time_index(Ls_query_360, LsDay) - Returns the indices that will extract data from the netCDF file - according to a range of solar longitudes [0-360]. + Returns the indices for a range of solar longitudes in a file. First try the Mars Year of the last timestep, then try the year before that. Use whichever Ls period is closest to the requested date. :param Ls_query_360: requested solar longitudes - :type Ls_query_360: list - + :type Ls_query_360: list :param LsDay: continuous solar longitudes - :type LsDay: array [areo] - + :type LsDay: array [areo] :return: file indices - :rtype: array - + :rtype: array :return: descriptor for the extracted solar longitudes - :rtype: str + :rtype: str + :raises ValueError: If the input Ls_query_360 is not a valid type + for solar longitude calculation. + :raises TypeError: If the input LsDay is not a valid type for + solar longitude calculation. + :raises Exception: If the time index calculation fails for any + reason. .. note:: The keyword ``all`` is passed as ``-99999`` by the ``rT()`` function - .. py:function:: get_tod_index(tod_query, tods) - Returns the indices that will extract data from the netCDF file - according to a range of *times of day*. + Returns the indices for a range of times of day in a file. :param tod_query: requested time of day (0-24) - :type tod_query: list - + :type tod_query: list :param tods: times of day - :type tods: array [tod] - + :type tods: array [tod] :return: file indices - :rtype: array [tod] - + :rtype: array [tod] :return: descriptor for the extracted time of day - :rtype: str + :rtype: str + :raises ValueError: If the input tod_query is not a valid type for + time of day calculation. .. note:: The keyword ``all`` is passed as ``-99999`` by the ``rT()`` function +.. py:function:: main() -.. py:function:: give_permission(filename) - - Sets group permissions for files created on NAS. - - :param filename: name of the file - :type filename: str + Main entry point for the MarsPlot script. + Handles argument parsing, global variable setup, figure object + initialization, and execution of the main plotting workflow. + Depending on the provided arguments, this function can: + - Inspect the contents of a NetCDF file and print variable + information or statistics. + - Generate a template configuration file. + - Parse a provided template file, select data based on optional + date bounds, and generate + diagnostic plots as individual files or as a merged + multipage PDF. + - Manage output directories and file naming conventions. + - Display progress and handle debug output. -.. py:function:: main() + Global variables are set for configuration and figure formatting. + The function also manages error handling and user feedback for + invalid arguments or file operations. .. py:function:: make_template() @@ -964,23 +1487,34 @@ Attributes Generate the ``Custom.in`` template file. :return: Custom.in blank template - + :rtype: file + :raises ValueError: If the input customFileIN is not a valid type + for template generation. + :raises TypeError: If the input customFileIN is not a valid type + for template generation. + :raises Exception: If the template generation fails for any + reason. .. py:function:: mean_func(arr, axis) + Calculate the mean of an array along a specified axis. + This function calculates a mean over the selected axis, ignoring or including NaN values as specified by ``show_NaN_in_slice`` in ``amescap_profile``. :param arr: the array to be averaged - :type arr: array - + :type arr: array :param axis: the axis over which to average the array - :type axis: int - + :type axis: int :return: the mean over the time axis - + :rtype: array + :raises ValueError: If the array is empty or the axis is out of bounds. + :raises RuntimeWarning: If the mean calculation encounters NaN values. + :raises TypeError: If the input array is not a valid type for mean + calculation. + :raises Exception: If the mean calculation fails for any reason. .. py:function:: namelist_parser(Custom_file) @@ -988,37 +1522,42 @@ Attributes Parse a ``Custom.in`` template. :param Custom_file: full path to ``Custom.in`` file - :type Custom_file: str - + :type Custom_file: str :return: updated global variables, ``FigLayout``, ``objectList`` - + ``panelList``, ``subplotList``, ``addLineList``, ``layoutList`` + :rtype: list + :raises ValueError: If the input Custom_file is not a valid type + for file name. .. py:function:: prep_file(var_name, file_type, simuID, sol_array) Open the file as a Dataset or MFDataset object depending on its - status on Lou. Note that the input arguments are typically - extracted from a ``varfull`` object (e.g., - ``03340.atmos_average.ucomp``) and not from a file whose disk - status is known beforehand. + status on Lou. Note that the input arguments are typically + extracted from a ``varfull`` object (e.g., + ``03340.atmos_average.ucomp``) and not from a file whose disk status + is known beforehand. :param var_name: variable to extract (e.g., ``ucomp``) - :type var_name: str - + :type var_name: str :param file_type: MGCM output file type (e.g., ``average``) - :type file_name: str - + :type file_name: str :param simuID: simulation ID number (e.g., 2 for 2nd simulation) - :type simuID: int - + :type simuID: int :param sol_array: date in file name (e.g., [3340,4008]) - :type sol_array: list - + :type sol_array: list :return: Dataset or MFDataset object; - (var_info) longname and units; - (dim_info) dimensions e.g., (``time``, ``lat``,``lon``); - (dims) shape of the array e.g., [133,48,96] - + :return: (var_info) longname and units; + :return: (dim_info) dimensions e.g., (``time``, ``lat``,``lon``); + :return: (dims) shape of the array e.g., [133,48,96] + :rtype: Dataset or MFDataset object, str, tuple, list + :raises ValueError: If the input var_name is not a valid type + for variable name. + :raises TypeError: If the input file_type is not a valid type + for file type. + :raises Exception: If the file preparation fails for any + reason. + :raises IOError: If the file is not found or cannot be opened. .. py:function:: progress(k, Nmax, txt='', success=True) @@ -1026,27 +1565,36 @@ Attributes Display a progress bar when performing heavy calculations. :param k: current iteration of the outer loop - :type k: float - + :type k: float :param Nmax: max iteration of the outer loop - :type Nmax: float - + :type Nmax: float :return: progress bar (EX: ``Running... [#---------] 10.64 %``) - + :rtype: str + :raises ValueError: If the input k is not a valid type for + progress bar. + :raises TypeError: If the input Nmax is not a valid type for + progress bar. + :raises Exception: If the progress bar creation fails for any + reason. .. py:function:: rT(typeIn='char') - Read template for the type of data expected. Returns value to + Read template for the type of data expected. + + Returns value to ``filter_input()``. :param typeIn: type of data expected: ``char``, ``float``, ``int``, ``bool``, defaults to ``char`` - :type typeIn: str, optional - + :type typeIn: str, optional :return: text input reformatted to ``[val1, val2]`` - :rtype: float or array - + :rtype: float or array + :raises ValueError: If the input typeIn is not a valid type for + filtering. + :raises TypeError: If the input typeIn is not a valid type for + filtering. + :raises Exception: If the filtering operation fails for any reason. .. py:function:: read_axis_options(axis_options_txt) @@ -1055,24 +1603,20 @@ Attributes :param axis_options_txt: a copy of the last line ``Axis Options`` in ``Custom.in`` templates - :type axis_options_txt: str - + :type axis_options_txt: str :return: X-axis bounds as a numpy array or ``None`` if undedefined - :rtype: array or None - + :rtype: array or None :return: Y-axis bounds as a numpy array or ``None`` if undedefined - :rtype: array or None - + :rtype: array or None :return: colormap (e.g., ``jet``, ``nipy_spectral``) or line options (e.g., ``--r`` for dashed red) - :rtype: str - + :rtype: str :return: linear (``lin``) or logarithmic (``log``) color scale - :rtype: str - + :rtype: str :return: projection (e.g., ``ortho -125,45``) - :rtype: str - + :rtype: str + :raises ValueError: If the input axis_options_txt is not a valid + type for axis options. .. py:function:: remove_whitespace(raw_input) @@ -1084,12 +1628,12 @@ Attributes :param raw_input: user input for variable, (e.g., ``[atmos_average.temp] + 2)`` - :type raw_input: str - + :type raw_input: str :return: raw_input without whitespaces (e.g., ``[atmos_average.temp]+2)`` - :rtype: str - + :rtype: str + :raises ValueError: If the input raw_input is not a valid type for + whitespace removal. .. py:function:: select_range(Ncdf_num, bound) @@ -1098,14 +1642,16 @@ Attributes within the user-defined range. :param Ncdf_num: a sorted array of sols - :type Ncdf_num: array - + :type Ncdf_num: array :param bound: a sol (e.g., 0350) or range of sols ``[min max]`` - :type bound: int or array - + :type bound: int or array :return: a sorted array of sols within the bounds - :rtype: array - + :rtype: array + :raises ValueError: If the input Ncdf_num is not a valid type for + file name. + :raises TypeError: If the input bound is not a valid type for + file name. + :raises Exception: If the range selection fails for any reason. .. py:function:: shift_data(lon, data) @@ -1113,45 +1659,41 @@ Attributes Shifts the longitude data from 0-360 to -180/180 and vice versa. :param lon: 1D array of longitude - :type lon: array [lon] - + :type lon: array [lon] :param data: 2D array with last dimension = longitude - :type data: array [1,lon] + :type data: array [1,lon] + :return: 1D array of longitude in -180/180 or 0-360 + :rtype: array [lon] + :return: 2D array with last dimension = longitude + :rtype: array [1,lon] + :raises ValueError: If the longitude coordinate type is invalid. + :raises TypeError: If the input data is not a valid type for + shifting. + :raises Exception: If the shifting operation fails for any reason. - :raises ValueError: Longitude coordinate type is not recognized. - - :return: longitude (-180/180) - :rtype: array [lon] - - :return: shifted data - :rtype: array [1,lon] .. note:: Use ``np.ma.hstack`` instead of ``np.hstack`` to keep the masked array properties - .. py:function:: split_varfull(varfull) Split ``varfull`` object into its component parts :param varfull: a ``varfull`` object (e.g, ``atmos_average@2.zsurf``, ``02400.atmos_average@2.zsurf``) - :type varfull: str - + :type varfull: str :return: (sol_array) a sol number or ``None`` (if none provided) - :rtype: int or None - + :rtype: int or None :return: (filetype) file type (e.g, ``atmos_average``) - :rtype: str - + :rtype: str :return: (var) variable of interest (e.g, ``zsurf``) - :rtype: str - + :rtype: str :return: (``simuID``) simulation ID (Python indexing starts at 0) - :rtype: int - + :rtype: int + :raises ValueError: If the input varfull is not a valid type for + splitting. .. py:data:: add_sol_time_axis @@ -1167,11 +1709,19 @@ Attributes +.. py:data:: debug + + + .. py:data:: degr :value: '°' +.. py:data:: exit_code + + + .. py:data:: include_NaNs diff --git a/docs/source/autoapi/bin/MarsPull/index.rst b/docs/source/autoapi/bin/MarsPull/index.rst index c56e5078..4fdecb4a 100644 --- a/docs/source/autoapi/bin/MarsPull/index.rst +++ b/docs/source/autoapi/bin/MarsPull/index.rst @@ -5,7 +5,9 @@ .. autoapi-nested-parse:: - The MarsPull executable is for querying data from the Mars Climate Modeling Center (MCMC) Mars Global Climate Model (MGCM) repository on the NASA NAS Data Portal at data.nas.nasa.gov/mcmc. + The MarsPull executable is for querying data from the Mars Climate + Modeling Center (MCMC) Mars Global Climate Model (MGCM) repository on + the NASA NAS Data Portal at data.nas.nasa.gov/mcmc. The executable requires 2 arguments: @@ -15,14 +17,15 @@ Third-party Requirements: - * ``numpy`` + * ``sys`` * ``argparse`` + * ``os`` + * ``re`` + * ``numpy`` + * ``functools`` + * ``traceback`` * ``requests`` - List of Functions: - - * download - Queries the requested file from the NAS Data Portal. - Module Contents @@ -34,9 +37,10 @@ Functions .. autoapisummary:: + bin.MarsPull.debug_wrapper bin.MarsPull.download - bin.MarsPull.file_list bin.MarsPull.main + bin.MarsPull.print_file_list @@ -48,30 +52,107 @@ Attributes bin.MarsPull.Ls_end bin.MarsPull.Ls_ini bin.MarsPull.args + bin.MarsPull.debug + bin.MarsPull.exit_code bin.MarsPull.parser - bin.MarsPull.saveDir + bin.MarsPull.save_dir -.. py:function:: download(url, filename) +.. py:function:: debug_wrapper(func) - Downloads a file from the NAS Data Portal (data.nas.nasa.gov). + A decorator that wraps a function with error handling + based on the --debug flag. + If the --debug flag is set, it prints the full traceback + of any exception that occurs. Otherwise, it prints a + simplified error message. - :param url: The url to download, e.g 'https://data.nas.nasa.gov/legacygcm/download_data.php?file=/legacygcmdata/LegacyGCM_Ls000_Ls004.nc' - :type url: str + :param func: The function to wrap. + :type func: function + :return: The wrapped function. + :rtype: function + :raises Exception: If an error occurs during the function call. + :raises TypeError: If the function is not callable. + :raises ValueError: If the function is not found. + :raises AttributeError: If the function does not have the + specified attribute. + :raises IndexError: If the function does not have the + specified index. - :param filename: The local filename e.g '/lou/la4/akling/Data/LegacyGCM_Ls000_Ls004.nc' - :type filename: str - :return: The requested file(s), downloaded and saved to the current directory. +.. py:function:: download(url, file_name) + Downloads a file from the NAS Data Portal (data.nas.nasa.gov). + The function takes a URL and a file name as input, and downloads the + file from the URL, saving it to the specified file name. It also + provides a progress bar to show the download progress if the file + size is known. If the file size is unknown, it simply downloads the + file without showing a progress bar. + The function handles errors during the download process and prints + appropriate messages to the console. + + :param url: The url to download from, e.g., + 'https://data.nas.nasa.gov/legacygcm/fv3betaout1data/03340.fixed.nc' + :type url: str + :param file_name: The local file_name e.g., + '/lou/la4/akling/Data/LegacyGCM_Ls000_Ls004.nc' + :type file_name: str + :return: The requested file(s), downloaded and saved to the current + directory. + :rtype: None :raises FileNotFoundError: A file-not-found error. + :raises PermissionError: A permission error. + :raises OSError: An operating system error. + :raises ValueError: A value error. + :raises TypeError: A type error. + :raises requests.exceptions.RequestException: A request error. + :raises requests.exceptions.HTTPError: An HTTP error. + :raises requests.exceptions.ConnectionError: A connection error. + :raises requests.exceptions.Timeout: A timeout error. + :raises requests.exceptions.TooManyRedirects: A too many redirects + error. + :raises requests.exceptions.URLRequired: A URL required error. + :raises requests.exceptions.InvalidURL: An invalid URL error. + :raises requests.exceptions.InvalidSchema: An invalid schema error. + :raises requests.exceptions.MissingSchema: A missing schema error. + :raises requests.exceptions.InvalidHeader: An invalid header error. + :raises requests.exceptions.InvalidProxyURL: An invalid proxy URL + error. + :raises requests.exceptions.InvalidRequest: An invalid request error. + :raises requests.exceptions.InvalidResponse: An invalid response + error. + +.. py:function:: main() + The main function that handles the command-line arguments -.. py:function:: file_list(list_of_files) + Handles the command-line arguments and coordinates the download + process. It checks for the presence of the required arguments, + validates the input, and calls the appropriate functions to download + the requested files. It also handles the logic for listing available + directories and files, as well as downloading files based on + specified solar longitudes (Ls) or file names. + :return: 0 if successful, 1 if an error occurred. + :rtype: int + :raises SystemExit: If an error occurs during the execution of the + program, the program will exit with a non-zero status code. -.. py:function:: main() + +.. py:function:: print_file_list(list_of_files) + + Prints a list of files. + + :param list_of_files: The list of files to print. + :type list_of_files: list + :return: None + :rtype: None + :raises TypeError: If list_of_files is not a list. + :raises ValueError: If list_of_files is empty. + :raises IndexError: If list_of_files is out of range. + :raises KeyError: If list_of_files is not found. + :raises OSError: If list_of_files is not accessible. + :raises IOError: If list_of_files is not open. .. py:data:: Ls_end @@ -86,11 +167,19 @@ Attributes +.. py:data:: debug + + + +.. py:data:: exit_code + + + .. py:data:: parser -.. py:data:: saveDir +.. py:data:: save_dir diff --git a/docs/source/autoapi/bin/MarsVars/index.rst b/docs/source/autoapi/bin/MarsVars/index.rst index 944cda27..2c597523 100644 --- a/docs/source/autoapi/bin/MarsVars/index.rst +++ b/docs/source/autoapi/bin/MarsVars/index.rst @@ -31,12 +31,21 @@ Third-party Requirements: - * ``numpy`` - * ``netCDF4`` + * ``sys`` * ``argparse`` * ``os`` - * ``subprocess`` + * ``warnings`` + * ``re`` + * ``numpy`` + * ``netCDF4`` + * ``shutil`` + * ``functools`` + * ``traceback`` * ``matplotlib`` + * ``time`` + * ``io`` + * ``locale`` + * ``amescap`` @@ -50,6 +59,8 @@ Functions .. autoapisummary:: bin.MarsVars.add_help + bin.MarsVars.check_dependencies + bin.MarsVars.check_variable_exists bin.MarsVars.compute_DP_3D bin.MarsVars.compute_DZ_3D bin.MarsVars.compute_DZ_full_pstd @@ -60,6 +71,8 @@ Functions bin.MarsVars.compute_Tco2 bin.MarsVars.compute_Vg_sed bin.MarsVars.compute_WMFF + bin.MarsVars.compute_dustref_per_pa + bin.MarsVars.compute_dustref_per_z bin.MarsVars.compute_mmr bin.MarsVars.compute_p_3D bin.MarsVars.compute_rho @@ -70,7 +83,17 @@ Functions bin.MarsVars.compute_xzTau bin.MarsVars.compute_zfull bin.MarsVars.compute_zhalf + bin.MarsVars.debug_wrapper + bin.MarsVars.ensure_file_closed + bin.MarsVars.force_close_netcdf_files + bin.MarsVars.get_existing_var_name bin.MarsVars.main + bin.MarsVars.patched_print_message + bin.MarsVars.process_add_variables + bin.MarsVars.safe_copy_replace + bin.MarsVars.safe_move_file + bin.MarsVars.safe_print + bin.MarsVars.safe_remove_file @@ -99,12 +122,15 @@ Attributes bin.MarsVars.amu_co2 bin.MarsVars.args bin.MarsVars.cap_str + bin.MarsVars.debug + bin.MarsVars.exit_code bin.MarsVars.filepath bin.MarsVars.fill_value bin.MarsVars.g bin.MarsVars.mass_co2 bin.MarsVars.master_list bin.MarsVars.n0 + bin.MarsVars.original_print_message bin.MarsVars.parser bin.MarsVars.psrf bin.MarsVars.rgas @@ -116,30 +142,64 @@ Attributes .. py:function:: add_help(var_list) + Create a help string for the add_variable argument. + + :param var_list: Dictionary of variables and their attributes + :type var_list: dict + :return: Formatted help string + :rtype: str + + +.. py:function:: check_dependencies(f, var, master_list, add_missing=True, dependency_chain=None) + + Check for variable dependencies in a file, add missing dependencies. + + :param f: NetCDF file object + :param var: Variable to check deps. for + :param master_list: Dict of supported vars and their deps. + :param add_missing: Whether to try adding missing deps. (default: True) + :param dependency_chain: List of vars in the current dep. chain (for detecting cycles) + :return: True if all deps. are present or successfully added, False otherwise + :raises RuntimeError: If the variable is not in the master list + :raises Exception: If any other error occurs + + +.. py:function:: check_variable_exists(var_name, file_vars) + + Check if a variable exists in a file. + + Considers alternative naming conventions. + + :param var_name: Variable name to check + :type var_name: str + :param file_vars: Set of variable names in the file + :type file_vars: set + :return: True if the variable exists, False otherwise + :rtype: bool + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs + .. py:function:: compute_DP_3D(ps, ak, bk, shape_out) Calculate the thickness of a layer in pressure units. :param ps: Surface pressure (Pa) - :type ps: array [time, lat, lon] - + :type ps: array [time, lat, lon] :param ak: Vertical coordinate pressure value (Pa) - :type ak: array [phalf] - + :type ak: array [phalf] :param bk: Vertical coordinate sigma value (None) - :type bk: array [phalf] - + :type bk: array [phalf] :param shape_out: Determines how to handle the dimensions of DP_3D. If len(time) = 1 (one timestep), DP_3D is returned as [1, lev, lat, lon] as opposed to [lev, lat, lon] - :type shape_out: float - - :raises: - + :type shape_out: float :return: ``DP`` Layer thickness in pressure units (Pa) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs .. py:function:: compute_DZ_3D(ps, ak, bk, temp, shape_out) @@ -147,30 +207,27 @@ Attributes Calculate the thickness of a layer in altitude units. :param ps: Surface pressure (Pa) - :type ps: array [time, lat, lon] - + :type ps: array [time, lat, lon] :param ak: Vertical coordinate pressure value (Pa) - :type ak: array [phalf] - + :type ak: array [phalf] :param bk: Vertical coordinate sigma value (None) - :type bk: array [phalf] - + :type bk: array [phalf] :param shape_out: Determines how to handle the dimensions of DZ_3D. If len(time) = 1 (one timestep), DZ_3D is returned as [1, lev, lat, lon] as opposed to [lev, lat, lon] - :type shape_out: float - - :raises: - + :type shape_out: float :return: ``DZ`` Layer thickness in altitude units (m) - :rtype: array [time, lev, lat, lon] + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs +.. py:function:: compute_DZ_full_pstd(pstd, T, ftype='average') -.. py:function:: compute_DZ_full_pstd(pstd, temp, ftype='average') + Calculate layer thickness. - Calculate the thickness of a layer from the midpoint of the - standard pressure levels (``pstd``). + Computes from the midpoint of the standard pressure levels (``pstd``). In this context, ``pfull=pstd`` with the layer interfaces defined somewhere in between successive layers:: @@ -185,54 +242,54 @@ Attributes / / / / :param pstd: Vertical coordinate (pstd; Pa) - :type pstd: array [lev] - - :param temp: Temperature (K) - :type temp: array [time, lev, lat, lon] - + :type pstd: array [lev] + :param T: Temperature (K) + :type T: array [time, lev, lat, lon] :param f_type: The FV3 file type: diurn, daily, or average - :type f_stype: str - - :raises: - + :type f_stype: str :return: DZ_full_pstd, Layer thicknesses (Pa) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs + :raises RuntimeError: If the layer thickness calculation fails + :raises ZeroDivisionError: If the temperature is zero .. py:function:: compute_Ek(ucomp, vcomp) - Calculate wave kinetic energ:: + Calculate wave kinetic energy + + Calculation:: Ek = 1/2 (u'**2+v'**2) :param ucomp: Zonal wind (m/s) - :type ucomp: array [time, lev, lat, lon] - + :type ucomp: array [time, lev, lat, lon] :param vcomp: Meridional wind (m/s) - :type vcomp: array [time, lev, lat, lon] - - :raises: - + :type vcomp: array [time, lev, lat, lon] :return: ``Ek`` Wave kinetic energy (J/kg) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs .. py:function:: compute_Ep(temp) - Calculate wave potential energy:: + Calculate wave potential energy. + + Calculation:: Ep = 1/2 (g/N)^2 (temp'/temp)^2 :param temp: Temperature (K) - :type temp: array [time, lev, lat, lon] - - :raises: - + :type temp: array [time, lev, lat, lon] :return: ``Ep`` Wave potential energy (J/kg) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs .. py:function:: compute_MF(UVcomp, w) @@ -240,16 +297,14 @@ Attributes Calculate zonal or meridional momentum fluxes. :param UVcomp: Zonal or meridional wind (ucomp or vcomp)(m/s) - :type UVcomp: array - + :type UVcomp: array :param w: Vertical wind (m/s) - :type w: array [time, lev, lat, lon] - - :raises: - + :type w: array [time, lev, lat, lon] :return: ``u'w'`` or ``v'w'``, Zonal/meridional momentum flux (J/kg) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs .. py:function:: compute_N(theta, zfull) @@ -257,57 +312,56 @@ Attributes Calculate the Brunt Vaisala freqency. :param theta: Potential temperature (K) - :type theta: array [time, lev, lat, lon] - + :type theta: array [time, lev, lat, lon] :param zfull: Altitude above ground level at the layer midpoint (m) - :type zfull: array [time, lev, lat, lon] - - :raises: - + :type zfull: array [time, lev, lat, lon] :return: ``N``, Brunt Vaisala freqency [rad/s] - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs .. py:function:: compute_Tco2(P_3D) Calculate the frost point of CO2. + Adapted from Fannale (1982) - Mars: The regolith-atmosphere cap system and climate change. Icarus. :param P_3D: The full 3D pressure array (Pa) - :type p_3D: array [time, lev, lat, lon] - - :raises: - + :type p_3D: array [time, lev, lat, lon] :return: CO2 frost point [K] - :rtype: array [time, lev, lat, lon] + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs - -.. py:function:: compute_Vg_sed(xTau, nTau, temp) +.. py:function:: compute_Vg_sed(xTau, nTau, T) Calculate the sedimentation rate of the dust. + [Courtney Batterson, 2023] :param xTau: Dust or ice MASS mixing ratio (ppm) - :type xTau: array [time, lev, lat, lon] - + :type xTau: array [time, lev, lat, lon] :param nTau: Dust or ice NUMBER mixing ratio (None) - :type nTau: array [time, lev, lat, lon] - - :param temp: Temperature (K) - :type temp: array [time, lev, lat, lon] - - :raises: - + :type nTau: array [time, lev, lat, lon] + :param T: Temperature (K) + :type T: array [time, lev, lat, lon] :return: ``Vg`` Dust sedimentation rate (m/s) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs + :raises RuntimeError: If the sedimentation rate calculation fails .. py:function:: compute_WMFF(MF, rho, lev, interp_type) - Calculate the zonal or meridional wave-mean flow forcing:: + Calculate the zonal or meridional wave-mean flow forcing. + + Calculation:: ax = -1/rho d(rho u'w')/dz ay = -1/rho d(rho v'w')/dz @@ -322,50 +376,90 @@ Attributes [du/dz = (du/dp).(-rho*g)] > [du/dz = -rho*g * (du/dp)] :param MF: Zonal/meridional momentum flux (J/kg) - :type MF: array [time, lev, lat, lon] - + :type MF: array [time, lev, lat, lon] :param rho: Atmospheric density (kg/m^3) - :type rho: array [time, lev, lat, lon] - + :type rho: array [time, lev, lat, lon] :param lev: Array for the vertical grid (zagl, zstd, pstd, or pfull) - :type lev: array [lev] - + :type lev: array [lev] :param interp_type: The vertical grid type (``zagl``, ``zstd``, ``pstd``, or ``pfull``) - :type interp_type: str + :type interp_type: str + :return: The zonal or meridional wave-mean flow forcing (m/s2) + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs + :raises RuntimeError: If the wave-mean flow forcing calculation fails + :raises ZeroDivisionError: If rho is zero - :raises: - :return: The zonal or meridional wave-mean flow forcing (m/s2) - :rtype: array [time, lev, lat, lon] +.. py:function:: compute_dustref_per_pa(dustref, delp) + + Computes visible dust opacity per Pascal from dustref and delp. + dustref is visible dust opacity per level (model layer). + + opacity per Pa = opacity per layer / layer thickness in Pa:: + + dustref_per_pa = dustref/delp + + [Courtney Batterson, 2025] + + :param dustref: Visible dust opacity [op/model layer] + :type dustref: array [time, lev, lat, lon] + :param delp: Layer thickness [Pa] + :type delp: array [time, lev, lat, lon] + :return: `dustref_per_pa` Visible dust opacity [op/Pa] + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs + + +.. py:function:: compute_dustref_per_z(dustref, delz) + + Computes visible dust opacity per kilometer from dustref and delz. + dustref is visible dust opacity per level (model layer). + + opacity per km = opacity per layer / layer thickness in m * 1000:: + dustref_per_z = dustref/delz*1000 + + [Courtney Batterson, 2025] + + :param dustref: Visible dust opacity [op/model layer] + :type dustref: array [time, lev, lat, lon] + :param delz: Layer thickness [m] + :type delz: array [time, lev, lat, lon] + :return: `dustref_per_z` Visible dust opacity [op/km] + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs .. py:function:: compute_mmr(xTau, temp, lev, const, f_type) Compute the dust or ice mixing ratio. + Adapted from Heavens et al. (2011) observations from MCS (JGR). + [Courtney Batterson, 2023] :param xTau: Dust or ice extinction rate (km-1) - :type xTau: array [time, lev, lat, lon] - + :type xTau: array [time, lev, lat, lon] :param temp: Temperature (K) - :type temp: array [time, lev, lat, lon] - + :type temp: array [time, lev, lat, lon] :param lev: Vertical coordinate (e.g., pstd) (e.g., Pa) - :type lev: array [lev] - + :type lev: array [lev] :param const: Dust or ice constant - :type const: array - + :type const: array :param f_type: The FV3 file type: diurn, daily, or average - :type f_stype: str - - :raises: - + :type f_stype: str :return: ``q``, Dust or ice mass mixing ratio (ppm) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs + :raises RuntimeError: If the mixing ratio calculation fails .. py:function:: compute_p_3D(ps, ak, bk, shape_out) @@ -373,24 +467,21 @@ Attributes Compute the 3D pressure at layer midpoints. :param ps: Surface pressure (Pa) - :type ps: array [time, lat, lon] - + :type ps: array [time, lat, lon] :param ak: Vertical coordinate pressure value (Pa) - :type ak: array [phalf] - + :type ak: array [phalf] :param bk: Vertical coordinate sigma value (None) - :type bk: array [phalf] - + :type bk: array [phalf] :param shape_out: Determines how to handle the dimensions of p_3D. If ``len(time) = 1`` (one timestep), ``p_3D`` is returned as [1, lev, lat, lon] as opposed to [lev, lat, lon] - :type shape_out: float - - :raises: - + :type shape_out: float :return: ``p_3D`` The full 3D pressure array (Pa) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs + :raises RuntimeError: If the pressure calculation fails .. py:function:: compute_rho(p_3D, temp) @@ -398,16 +489,14 @@ Attributes Compute density. :param p_3D: Pressure (Pa) - :type p_3D: array [time, lev, lat, lon] - + :type p_3D: array [time, lev, lat, lon] :param temp: Temperature (K) - :type temp: array [time, lev, lat, lon] - - :raises: - + :type temp: array [time, lev, lat, lon] :return: Density (kg/m^3) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs .. py:function:: compute_scorer(N, ucomp, zfull) @@ -415,42 +504,35 @@ Attributes Calculate the Scorer wavelength. :param N: Brunt Vaisala freqency (rad/s) - :type N: float [time, lev, lat, lon] - + :type N: float [time, lev, lat, lon] :param ucomp: Zonal wind (m/s) - :type ucomp: array [time, lev, lat, lon] - + :type ucomp: array [time, lev, lat, lon] :param zfull: Altitude above ground level at the layer midpoint (m) - :type zfull: array [time, lev, lat, lon] - - :raises: - + :type zfull: array [time, lev, lat, lon] :return: ``scorer_wl`` Scorer horizontal wavelength (m) - :rtype: array [time, lev, lat, lon] + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs - -.. py:function:: compute_theta(p_3D, ps, temp, f_type) +.. py:function:: compute_theta(p_3D, ps, T, f_type) Compute the potential temperature. :param p_3D: The full 3D pressure array (Pa) - :type p_3D: array [time, lev, lat, lon] - + :type p_3D: array [time, lev, lat, lon] :param ps: Surface pressure (Pa) - :type ps: array [time, lat, lon] - - :param temp: Temperature (K) - :type temp: array [time, lev, lat, lon] - + :type ps: array [time, lat, lon] + :param T: Temperature (K) + :type T: array [time, lev, lat, lon] :param f_type: The FV3 file type: diurn, daily, or average - :type f_type: str - - :raises: - + :type f_type: str :return: Potential temperature (K) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs .. py:function:: compute_w(rho, omega) @@ -471,112 +553,347 @@ Attributes omega = -rho * g * w :param rho: Atmospheric density (kg/m^3) - :type rho: array [time, lev, lat, lon] - + :type rho: array [time, lev, lat, lon] :param omega: Rate of change in pressure at layer midpoint (Pa/s) - :type omega: array [time, lev, lat, lon] - - :raises: - + :type omega: array [time, lev, lat, lon] :return: vertical wind (m/s) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs + :raises RuntimeError: If the vertical wind calculation fails + :raises ZeroDivisionError: If rho or omega is zero + :raises OverflowError: If the calculation results in an overflow + :raises Exception: If any other error occurs .. py:function:: compute_w_net(Vg, wvar) - Computes the net vertical wind, which is the vertical wind (w) - minus the sedimentation rate (``Vg_sed``):: + Computes the net vertical wind. + + w_net = vertical wind (w) - sedimentation rate (``Vg_sed``):: w_net = w - Vg_sed - :param Vg: Dust sedimentation rate (m/s) - :type Vg: array [time, lev, lat, lon] + [Courtney Batterson, 2023] + :param Vg: Dust sedimentation rate (m/s) + :type Vg: array [time, lev, lat, lon] :param wvar: Vertical wind (m/s) - :type wvar: array [time, lev, lat, lon] - - :raises: - + :type wvar: array [time, lev, lat, lon] :return: `w_net` Net vertical wind speed (m/s) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs .. py:function:: compute_xzTau(q, temp, lev, const, f_type) Compute the dust or ice extinction rate. + Adapted from Heavens et al. (2011) observations from MCS (JGR). + [Courtney Batterson, 2023] :param q: Dust or ice mass mixing ratio (ppm) - :type q: array [time, lev, lat, lon] - + :type q: array [time, lev, lat, lon] :param temp: Temperature (K) - :type temp: array [time, lev, lat, lon] - + :type temp: array [time, lev, lat, lon] :param lev: Vertical coordinate (e.g., pstd) (e.g., Pa) - :type lev: array [lev] - + :type lev: array [lev] :param const: Dust or ice constant - :type const: array - + :type const: array :param f_type: The FV3 file type: diurn, daily, or average - :type f_stype: str - - :raises: - + :type f_stype: str :return: ``xzTau`` Dust or ice extinction rate (km-1) - :rtype: array [time, lev, lat, lon] - + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs + :raises RuntimeError: If the extinction rate calculation fails -.. py:function:: compute_zfull(ps, ak, bk, temp) +.. py:function:: compute_zfull(ps, ak, bk, T) Calculate the altitude of the layer midpoints above ground level. :param ps: Surface pressure (Pa) - :type ps: array [time, lat, lon] - + :type ps: array [time, lat, lon] :param ak: Vertical coordinate pressure value (Pa) - :type ak: array [phalf] - + :type ak: array [phalf] :param bk: Vertical coordinate sigma value (None) - :type bk: array [phalf] - - :param temp: Temperature (K) - :type temp: array [time, lev, lat, lon] - - :raises: - + :type bk: array [phalf] + :param T: Temperature (K) + :type T: array [time, lev, lat, lon] :return: ``zfull`` (m) - :rtype: array [time, lev, lat, lon] + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs - -.. py:function:: compute_zhalf(ps, ak, bk, temp) +.. py:function:: compute_zhalf(ps, ak, bk, T) Calculate the altitude of the layer interfaces above ground level. :param ps: Surface pressure (Pa) - :type ps: array [time, lat, lon] - + :type ps: array [time, lat, lon] :param ak: Vertical coordinate pressure value (Pa) - :type ak: array [phalf] - + :type ak: array [phalf] :param bk: Vertical coordinate sigma value (None) - :type bk: array [phalf] + :type bk: array [phalf] + :param T: Temperature (K) + :type T: array [time, lev, lat, lon] + :return: ``zhalf`` (m) + :rtype: array [time, lev, lat, lon] + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs + + +.. py:function:: debug_wrapper(func) + + A decorator that wraps a function with error handling + based on the --debug flag. + If the --debug flag is set, it prints the full traceback + of any exception that occurs. Otherwise, it prints a + simplified error message. + + :param func: The function to wrap. + :type func: function + :return: The wrapped function. + :rtype: function + :raises Exception: If an error occurs during the function call. + :raises TypeError: If the function is not callable. + :raises ValueError: If the function is not found. + :raises NameError: If the function is not defined. + :raises AttributeError: If the function does not have the + specified attribute. + :raises ImportError: If the function cannot be imported. + :raises RuntimeError: If the function cannot be run. + :raises KeyError: If the function does not have the + specified key. + :raises IndexError: If the function does not have the + specified index. + :raises IOError: If the function cannot be opened. + :raises OSError: If the function cannot be accessed. + :raises EOFError: If the function cannot be read. + :raises MemoryError: If the function cannot be allocated. + :raises OverflowError: If the function cannot be overflowed. + :raises ZeroDivisionError: If the function cannot be divided by zero. + :raises StopIteration: If the function cannot be stopped. + :raises KeyboardInterrupt: If the function cannot be interrupted. + :raises SystemExit: If the function cannot be exited. + :raises AssertionError: If the function cannot be asserted. + + +.. py:function:: ensure_file_closed(filepath, delay=0.5) + + Try to ensure a file is not being accessed by the system. + + This is especially helpful for Windows environments. + + :param filepath: Path to the file + :param delay: Delay in seconds to wait for handles to release + :return: None + :rtype: None + :raises FileNotFoundError: If the file does not exist + :raises OSError: If the file is locked or cannot be accessed + :raises Exception: If any other error occurs + :raises TypeError: If the filepath is not a string + :raises ValueError: If the filepath is empty + :raises RuntimeError: If the file cannot be closed + + +.. py:function:: force_close_netcdf_files(file_or_dir, delay=1.0) + + Aggressively try to ensure netCDF files are closed on Windows systems. + + :param file_or_dir: Path to the file or directory to process + :param delay: Delay in seconds after forcing closure + :return: None + :rtype: None + :raises FileNotFoundError: If the file or directory does not exist + :raises OSError: If the file is locked or cannot be accessed + :raises Exception: If any other error occurs + :raises TypeError: If the file_or_dir is not a string + :raises ValueError: If the file_or_dir is empty + :raises RuntimeError: If the file or directory cannot be processed + :raises ImportError: If the netCDF4 module is not available + :raises AttributeError: If the netCDF4 module does not have the required attributes + :raises Exception: If any other error occurs + + +.. py:function:: get_existing_var_name(var_name, file_vars) + + Get the actual variable name that exists in the file. + + Considers alternative naming conventions. + + :param var_name: Variable name to check + :type var_name: str + :param file_vars: Set of variable names in the file + :type file_vars: set + :return: Actual variable name in the file + :rtype: str + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs - :param temp: Temperature (K) - :type temp: array [time, lev, lat, lon] - :raises: +.. py:function:: main() - :return: ``zhalf`` (m) - :rtype: array [time, lev, lat, lon] + Main function for variable manipulations in NetCDF files. + + This function performs a sequence of operations on one or more + NetCDF files, as specified by command-line arguments. The operations + include removing variables, extracting variables, adding variables, + vertical differentiation, zonal detrending, opacity conversions, + column integration, and editing variable metadata or values. + + Workflow: + - Iterates over all input NetCDF files. + - For each file, performs the following operations as requested + by arguments: + * Remove specified variables and update the file. + * Extract specified variables into a new file. + * Add new variables using provided methods. + * Compute vertical derivatives of variables with respect to + height or pressure. + * Remove zonal mean (detrend) from specified variables. + * Convert variables between dp/dz and dz/dp representations. + * Perform column integration of variables. + * Edit variable metadata (name, long_name, units) or scale + values. + + Arguments: + args: Namespace + Parsed command-line arguments specifying which operations + to perform and their parameters. + master_list: list + List of available variables and their properties (used for + adding variables). + debug: bool + If True, prints detailed error messages and stack traces. + Notes: + - Handles both Unix and Windows file operations for safe file + replacement. + - Uses helper functions for NetCDF file manipulation, variable + existence checks, and error handling. + - Assumes global constants and utility functions (e.g., Dataset, + Ncdf, check_file_tape, etc.) are defined elsewhere. + - Uses global variables lev_T and lev_T_out for axis + manipulation in vertical operations. + + Raises: + Exceptions are caught and logged for each operation; files are + cleaned up on error. + + +.. py:function:: patched_print_message(self, message, file=None) + + Patched version of _print_message that handles Unicode encoding errors. + + :param self: The ArgumentParser instance + :param message: The message to print + :param file: The file to print to (default is sys.stdout) + :type file: file-like object + :return: None + :rtype: None + :raises UnicodeEncodeError: If the message cannot be encoded + :raises TypeError: If the message is not a string + :raises ValueError: If the message is empty + :raises Exception: If any other error occurs + :raises RuntimeError: If the message cannot be printed + + +.. py:function:: process_add_variables(file_name, add_list, master_list, debug=False) + + Process a list of variables to add. + + Dependent variables are added in the correct order. + If a variable is already in the file, it is skipped. + If a variable cannot be added, an error message is printed. + + :param file_name: Input file path + :param add_list: List of variables to add + :param master_list: Dictionary of supported variables and their dependencies + :param debug: Whether to show debug information + :type debug: bool + :raises ValueError: If the input dimensions are not compatible + :raises TypeError: If the input types are not compatible + :raises Exception: If any other error occurs + :raises RuntimeError: If the variable cannot be added + + +.. py:function:: safe_copy_replace(src_file, dst_file, max_attempts=5, delay=1.0) + + Windows-specific approach to copy file contents and replace destination. + + This avoids move operations which are more likely to fail with locking + + :param src_file: Source file path + :param dst_file: Destination file path + :param max_attempts: Maximum number of retry attempts + :param delay: Base delay between attempts (increases with retries) + :return: True if successful, False otherwise + :rtype: bool + :raises FileNotFoundError: If the source file does not exist + :raises OSError: If the file is locked or cannot be accessed + :raises Exception: If any other error occurs + :raises TypeError: If the src_file or dst_file is not a string + :raises ValueError: If the src_file or dst_file is empty + :raises RuntimeError: If the file cannot be copied or replaced + + +.. py:function:: safe_move_file(src_file, dst_file, max_attempts=5, delay=1) + + Safely move a file with retries for Windows file locking issues. + + :param src_file: Source file path + :param dst_file: Destination file path + :param max_attempts: Number of attempts to make + :param delay: Delay between attempts in seconds + :return: True if successful, False otherwise + :rtype: bool + :raises FileNotFoundError: If the source file does not exist + :raises OSError: If the file is locked or cannot be accessed + :raises Exception: If any other error occurs + :raises TypeError: If the src_file or dst_file is not a string + :raises ValueError: If the src_file or dst_file is empty + :raises RuntimeError: If the file cannot be moved + + +.. py:function:: safe_print(text) + + Print text safely, handling encoding issues on Windows. + + :param text: Text to print + :type text: str + :return: None + :rtype: None + :raises UnicodeEncodeError: If the text cannot be encoded + :raises TypeError: If the text is not a string + :raises ValueError: If the text is empty + :raises Exception: If any other error occurs + :raises RuntimeError: If the text cannot be printed +.. py:function:: safe_remove_file(filepath, max_attempts=5, delay=1) -.. py:function:: main() + Safely remove a file with retries for Windows file locking issues + + :param filepath: Path to the file to remove + :param max_attempts: Number of attempts to make + :param delay: Delay between attempts in seconds + :return: True if successful, False otherwise + :rtype: bool + :raises FileNotFoundError: If the file does not exist + :raises OSError: If the file is locked or cannot be accessed + :raises Exception: If any other error occurs + :raises TypeError: If the filepath is not a string + :raises ValueError: If the filepath is empty + :raises RuntimeError: If the file cannot be removed .. py:data:: C_dst @@ -673,6 +990,14 @@ Attributes +.. py:data:: debug + + + +.. py:data:: exit_code + + + .. py:data:: filepath @@ -699,6 +1024,10 @@ Attributes +.. py:data:: original_print_message + + + .. py:data:: parser diff --git a/docs/source/autoapi/bin/index.rst b/docs/source/autoapi/bin/index.rst index 121d61e2..29a90bc7 100644 --- a/docs/source/autoapi/bin/index.rst +++ b/docs/source/autoapi/bin/index.rst @@ -14,6 +14,7 @@ Submodules MarsFiles/index.rst MarsFormat/index.rst MarsInterp/index.rst + MarsNest/index.rst MarsPlot/index.rst MarsPull/index.rst MarsVars/index.rst diff --git a/docs/source/installation.rst b/docs/source/installation.rst index 2bd98770..bbad6bbb 100644 --- a/docs/source/installation.rst +++ b/docs/source/installation.rst @@ -470,6 +470,8 @@ The pip installation method is less recommended as it requires manual installati # Install CAP with spectral analysis support pip install "amescap[spectral] @ git+https://github.com/NASA-Planetary-Science/AmesCAP.git" + # OR for a specific branch, e.g., the devel branch + pip install "amescap[spectral] @ git+https://github.com/NASA-Planetary-Science/AmesCAP.git@devel" # Don't forget to copy the profile file to your home directory cp amescap/mars_templates/amescap_profile ~/.amescap_profile diff --git a/tests/create_gcm_files.py b/tests/create_gcm_files.py index 1c5a74d0..1ac31d21 100644 --- a/tests/create_gcm_files.py +++ b/tests/create_gcm_files.py @@ -7,8 +7,10 @@ import numpy as np from netCDF4 import Dataset -import os -import datetime + +# ---------------------------------------------------------------------- +# EMARS Dummy File +# ---------------------------------------------------------------------- def create_emars_test(): """Create emars_test.nc with the exact variables and structure as real EMARS files.""" @@ -475,6 +477,10 @@ def create_var(name, dimensions, units, longname, min_val, max_val, data_type=np nc_file.close() print("Created emars_test.nc") +# ---------------------------------------------------------------------- +# OpenMARS Dummy File +# ---------------------------------------------------------------------- + def create_openmars_test(): """Create openmars_test.nc with the exact variables and structure as real OpenMARS files.""" nc_file = Dataset('openmars_test.nc', 'w', format='NETCDF4') @@ -532,6 +538,10 @@ def create_var(name, dimensions, units, min_val, max_val, data_type=np.float32): nc_file.close() print("Created openmars_test.nc") +# ---------------------------------------------------------------------- +# PCM Dummy File +# ---------------------------------------------------------------------- + def create_pcm_test(): """Create pcm_test.nc with the exact variables and structure as real PCM files.""" nc_file = Dataset('pcm_test.nc', 'w', format='NETCDF4') @@ -775,6 +785,10 @@ def create_var(name, dimensions, units, longname, min_val, max_val, data_type=np nc_file.close() print("Created pcm_test.nc") +# ---------------------------------------------------------------------- +# MarsWRF Dummy File +# ---------------------------------------------------------------------- + def create_marswrf_test(): """Create marswrf_test.nc with the exact variables and structure as real MarsWRF files.""" nc_file = Dataset('marswrf_test.nc', 'w', format='NETCDF4') diff --git a/tests/test_marsfiles.py b/tests/test_marsfiles.py index 1eebbecd..ba8b70fd 100644 --- a/tests/test_marsfiles.py +++ b/tests/test_marsfiles.py @@ -502,6 +502,14 @@ def test_temporal_filters(self): """Test all temporal filtering operations""" # High-pass filter result = self.run_mars_files(['01336.atmos_daily.nc', '-hpt', '10', '-incl', 'temp']) + + # TEMPORARY DEBUG OUTPUT + if result.returncode != 0: + print("\n=== HIGH-PASS TEMPORAL FILTER FAILED ===") + print("STDOUT:", result.stdout) + print("STDERR:", result.stderr) + print("========================================\n") + self.assertEqual(result.returncode, 0, "High-pass temporal filter command failed") high_pass_file = self.check_file_exists('01336.atmos_daily_hpt.nc') self.verify_netcdf_has_variable(high_pass_file, 'temp') @@ -529,6 +537,14 @@ def test_spatial_filters(self): # High-pass filter result = self.run_mars_files(['01336.atmos_daily.nc', '-hps', '10', '-incl', 'temp']) + + # TEMPORARY DEBUG OUTPUT + if result.returncode != 0: + print("\n=== HIGH-PASS SPATIAL FILTER FAILED ===") + print("STDOUT:", result.stdout) + print("STDERR:", result.stderr) + print("=======================================\n") + self.assertEqual(result.returncode, 0, "High-pass spatial filter command failed") high_pass_file = self.check_file_exists('01336.atmos_daily_hps.nc') self.verify_netcdf_has_variable(high_pass_file, 'temp') @@ -686,6 +702,13 @@ def test_regrid(self): '-regrid', '01336.atmos_average_pstd_c48.nc' ]) + + # TEMPORARY DEBUG OUTPUT + if result.returncode != 0: + print("\n=== REGRID COMMAND FAILED ===") + print("STDOUT:", result.stdout) + print("STDERR:", result.stderr) + print("============================\n") # Check for successful execution self.assertEqual(result.returncode, 0, "Regrid command failed")