From c6a452c440825ed52849d03d4f3df623cae04e4c Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Tue, 26 Mar 2024 13:28:51 -0600
Subject: [PATCH 01/20] Updated Violin Plot with bigger label fonts and ticks
 consistent with other plots in aircraftplots

---
 ...m_aircraft_Latestfor_develop_aircraft.yaml |  2 +-
 ...elop_aircraft_testmultiplemodelgroups.yaml |  2 +-
 melodies_monet/plots/aircraftplots.py         | 20 +++++++++----------
 3 files changed, 12 insertions(+), 12 deletions(-)

diff --git a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft.yaml b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft.yaml
index 0c22109d..5280f4c7 100644
--- a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft.yaml
+++ b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft.yaml
@@ -153,7 +153,7 @@ plots:
     fig_kwargs:
       figsize: [10, 8]
     text_kwargs:
-      fontsize: 25.
+      fontsize: 20.
     domain_type: ['all']
     domain_name: ['CONUS']
     data: ['firexaq_wrfchem_v4.2']
diff --git a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
index 01619ffc..8e4db45c 100644
--- a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
+++ b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
@@ -190,7 +190,7 @@ plots:
     fig_kwargs:
       figsize: [10, 8]
     text_kwargs:
-      fontsize: 25.
+      fontsize: 20.
     domain_type: ['all']
     domain_name: ['CONUS']
     data: ['firexaq_wrfchem_v4.2','firexaq_wrfchem_v4.2_test']
diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index 3ee5a03a..ffebf422 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -776,23 +776,23 @@ def make_violin_plot(comb_violin, label_violin, outname='plot',
     # Convert the DataFrame to long-form or "tidy" format suitable for sns.violinplot
     melted_comb_violin = pd.melt(comb_violin, var_name='group', value_name='value')
     
+    # Set labels and title using text_dict
+    # Set default text size, modify here for bigger text
+    def_text = dict(fontsize=14)  # can increase fontsize for default text (> 14) 
+    text_kwargs = {**def_text, **text_dict} if text_dict else def_text
+
+
     # Create the violin plot
-    
-    # Updated sns.violinplot call
     # Use 'hue' parameter and set 'orient' to 'v' for vertical orientation
     sns.violinplot(x='group', y='value', data=melted_comb_violin, hue='group', palette=palette, cut=0, orient='v', density_norm='width', inner='quartile')
 
 
-
-
-
-    # Set labels and title using text_dict
-    
-    # Set labels and title
-    def_text = dict(fontsize=14)
-    text_kwargs = {**def_text, **text_dict} if text_dict else def_text
+    # Set labels and title with increased size
     plt.xlabel('', weight='bold', fontsize=text_kwargs['fontsize'])
     plt.ylabel(ylabel if ylabel else 'Value', weight='bold', fontsize=text_kwargs['fontsize'])
+
+    # Increase tick label size
+    plt.tick_params(axis='both', labelsize=text_kwargs['fontsize'])
     
     # Set y-axis limits if provided
     if vmin is not None and vmax is not None:

From 5f32335c3b70303bcd2d3a3ab93122daf01ea8a2 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Tue, 9 Apr 2024 09:37:22 -0600
Subject: [PATCH 02/20] Violin plot title and ticklabel fontsizes fixed

---
 melodies_monet/driver.py              | 10 +++++++++-
 melodies_monet/plots/aircraftplots.py | 14 +++++++++++---
 2 files changed, 20 insertions(+), 4 deletions(-)

diff --git a/melodies_monet/driver.py b/melodies_monet/driver.py
index b39c10ef..3b7d4beb 100644
--- a/melodies_monet/driver.py
+++ b/melodies_monet/driver.py
@@ -1353,7 +1353,7 @@ def plotting(self):
                                 elif filter_op == 'isnotin':
                                     pairdf_all.query(f'{column} != {filter_vals}', inplace=True)
                                 else:
-                                    pairdf_all.query(f'{column} {filter_op} {filter_vals}', inplace=True)
+                                    pairdf_all.query(f'{column} {filter_op} {filter_vals}', inplace=True) 
                         elif 'filter_string' in grp_dict['data_proc']:
                             pairdf_all.query(grp_dict['data_proc']['filter_string'], inplace=True)
 
@@ -1554,6 +1554,8 @@ def plotting(self):
                                 savefig(outname + '.png', logo_height=150)
                                 del (ax, fig_dict, plot_dict, text_dict, obs_dict, obs_plot_dict)  # Clear axis for next plot.
                                 
+
+
                             # At the end save the plot.
                             ##if p_index == len(pair_labels) - 1:
                                 #Adding Altitude variable as secondary y-axis to timeseries (for, model vs aircraft) qzr++
@@ -1565,6 +1567,10 @@ def plotting(self):
                                   ##  ax = airplots.add_yax2_altitude(ax, pairdf, altitude_variable, altitude_ticks, text_kwargs)
                                 ##savefig(outname + '.png', logo_height=150)
                                 ##del (ax, fig_dict, plot_dict, text_dict, obs_dict, obs_plot_dict) #Clear axis for next plot.
+
+                                                                     
+
+                           
                                 
                         #qzr++ Added vertprofile plotype for aircraft vs model comparisons         
                         elif plot_type.lower() == 'vertprofile':
@@ -1629,6 +1635,8 @@ def plotting(self):
                                 del (ax, fig_dict, plot_dict, text_dict, obs_dict, obs_plot_dict) # Clear axis for next plot.
 
                         
+
+                        
                         elif plot_type.lower() == 'violin':
                             if set_yaxis:
                                 if all(k in obs_plot_dict for k in ('vmin_plot', 'vmax_plot')):
diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index ffebf422..2ce533dc 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -24,6 +24,11 @@
 from matplotlib.colors import TwoSlopeNorm, ListedColormap, LinearSegmentedColormap
 from matplotlib.patches import Rectangle
 from matplotlib.ticker import FuncFormatter
+
+import matplotlib.dates as mdates
+from matplotlib.colors import BoundaryNorm
+from matplotlib.colors import Normalize
+
 from monet.util.tools import get_epa_region_bounds as get_epa_bounds 
 import math
 from ..plots import savefig
@@ -222,7 +227,9 @@ def add_yax2_altitude(ax, pairdf, altitude_yax2, text_kwargs, vmin_y2, vmax_y2):
 
     return ax
 
-                              
+
+
+
 ####NEW vertprofile has option for both shading (for interquartile range) or box (interquartile range)-whisker (10th-90th percentile bounds) (qzr++)
 def make_vertprofile(df, column=None, label=None, ax=None, bins=None, altitude_variable=None, ylabel=None,
                      vmin=None, vmax=None, 
@@ -792,7 +799,8 @@ def make_violin_plot(comb_violin, label_violin, outname='plot',
     plt.ylabel(ylabel if ylabel else 'Value', weight='bold', fontsize=text_kwargs['fontsize'])
 
     # Increase tick label size
-    plt.tick_params(axis='both', labelsize=text_kwargs['fontsize'])
+    plt.tick_params(axis='both', labelsize=text_kwargs['fontsize']*0.8)
+    
     
     # Set y-axis limits if provided
     if vmin is not None and vmax is not None:
@@ -800,7 +808,7 @@ def make_violin_plot(comb_violin, label_violin, outname='plot',
     
     # Set title
     if domain_type and domain_name:
-        plt.title(f"Violin Plot for {domain_type} - {domain_name}", fontweight='bold')
+        plt.title(f"Violin Plot for {domain_type} - {domain_name}", fontweight='bold',fontsize=text_kwargs['fontsize'])
 
     # Finalize and save plot
     plt.tight_layout()

From 7b87cdaeecc9deec3d81fdf9e192249e7fb7d7a1 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Mon, 15 Apr 2024 23:38:35 -0600
Subject: [PATCH 03/20] Curtain pLot of Model with Observations scatter
 overlay: separate for each model-obs pair like scatter plots

---
 ...elop_aircraft_testmultiplemodelgroups.yaml |  36 ++++-
 melodies_monet/driver.py                      |  56 ++++++-
 melodies_monet/plots/aircraftplots.py         | 152 +++++++++++++++++-
 3 files changed, 237 insertions(+), 7 deletions(-)

diff --git a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
index 8e4db45c..c29bbe0c 100644
--- a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
+++ b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
@@ -161,6 +161,7 @@ plots:
           altitude:
             oper: "between"
             value: [2000,8000] # values are [vmin_y2, vmax_y2] 
+            
   plot_grp2:
     type: 'vertprofile'  # plot type
     fig_kwargs:  # Opt to define figure options
@@ -197,7 +198,6 @@ plots:
     data_proc:
       rem_obs_nan: True 
       set_axis: False
-
     
   plot_grp4:
     type: 'scatter_density'
@@ -235,7 +235,35 @@ plots:
     # Add other seaborn.kdeplot keyword arguments here as needed
     vcenter: #0  # Central value for TwoSlopeNorm
     extensions: ['min', 'max']  # Extensions for the colorbar  
+
   plot_grp5:
+    type: 'curtain'
+    fig_kwargs:
+      figsize: [12, 8]
+    text_kwargs:
+      fontsize: 18
+    domain_type: ['all']
+    domain_name: ['CONUS']
+    data: ['firexaq_wrfchem_v4.2','firexaq_wrfchem_v4.2_test']
+    data_proc:
+      rem_obs_nan: True
+      set_axis: True
+      ts_select_time: 'time'
+      #ts_avg_window: # 'H'  #Uncomment to use time averaging
+      altitude_variable: 'altitude'
+      vmin_plot: 1000
+      vmax_plot: 10000
+    # Optionally specify a custom colormap and number of levels
+    ##If want to use standard matplotlib colorscale comment the True option and colors and uncomment the following:
+    ##color_map_custom: True
+    ##colors: ['blue','royalblue', 'cyan', 'yellow', 'orange', 'red'] ##['#440154', '#3b528b', '#21908d', '#5dc963', '#fde725'] #['purple', 'blue', 'green', 'yellow'] # Example gradient
+    color_map_custom: False
+    color_map: 'jet'  # This tells the plotting function to use the built-in 'rainbow' colormap
+    color_levels: 15 #10  # Define the number of distinct colors in the color bar
+    # Optionally specify dynamic color bar ranges based on data
+    dynamic_color_range: True  # New key to use dynamic ranges based on model and obs data
+    
+  plot_grp6:
     type: 'taylor' # plot type
     fig_kwargs: #Opt to define figure options
       figsize: [8,8] # figure size if multiple plots
@@ -250,7 +278,8 @@ plots:
     data_proc:
       rem_obs_nan: True # True: Remove all points where model or obs variable is NaN. False: Remove only points where model variable is NaN.
       set_axis: True #If select True, add ty_scale for each variable in obs.
-  plot_grp6:
+      
+  plot_grp7:
     type: 'boxplot' # plot type
     fig_kwargs: #Opt to define figure options
       figsize: [8, 6] # figure size 
@@ -262,6 +291,8 @@ plots:
     data_proc:
       rem_obs_nan: True # True: Remove all points where model or obs variable is NaN. False: Remove only points where model variable is NaN.
       set_axis: False #If select True, add vmin_plot and vmax_plot for each variable in obs.
+
+      
 stats:
   #Stats require positive numbers, so if you want to calculate temperature use Kelvin!
   #Wind direction has special calculations for AirNow if obs name is 'WD'
@@ -279,4 +310,5 @@ stats:
   domain_name: ['CONUS'] #List of domain names. If domain_type = all domain_name is used in plot title.
   data: ['firexaq_wrfchem_v4.2','firexaq_wrfchem_v4.2_test'] # make this a list of pairs in obs_model where the obs is the obs label and model is the model_label
 
+  
 
diff --git a/melodies_monet/driver.py b/melodies_monet/driver.py
index 3b7d4beb..6a7a7f93 100644
--- a/melodies_monet/driver.py
+++ b/melodies_monet/driver.py
@@ -1568,9 +1568,61 @@ def plotting(self):
                                 ##savefig(outname + '.png', logo_height=150)
                                 ##del (ax, fig_dict, plot_dict, text_dict, obs_dict, obs_plot_dict) #Clear axis for next plot.
 
-                                                                     
+                        
+                        elif plot_type.lower() == 'curtain':
+                            curtain_config = grp_dict
+                            vmin = curtain_config['data_proc'].get('vmin_plot', None)
+                            vmax = curtain_config['data_proc'].get('vmax_plot', None)
+                            
+                            model_label = p.model
+                            obs_label = p.obs
+                        
+                            if 'altitude' in pairdf.columns and vmin is not None and vmax is not None:
+                                pairdf_filtered = pairdf[(pairdf['altitude'] >= vmin) & (pairdf['altitude'] <= vmax)]
+                            else:
+                                pairdf_filtered = pairdf
+                        
+                            try:
+                                model_data_pivot = pairdf_filtered.pivot(index='time', columns='altitude', values=modvar)
+                                obs_data_pivot = pairdf_filtered.pivot(index='time', columns='altitude', values=obsvar)
+                        
+                                # Convert datetime index to float seconds since the first timestamp
+                                time_values = model_data_pivot.index.values
+                                time_as_float = (time_values - time_values[0]).astype('timedelta64[s]').astype(float)
+                                
+                                altitude_values = model_data_pivot.columns.values
+                                model_data_2d = model_data_pivot.values
+                                obs_data_2d = obs_data_pivot.values
+                        
+                                model_data_2d_cleaned = airplots.interpolate_and_clean_model_data(time_as_float, altitude_values, model_data_2d)
+                                #Create Curtain Plots
+                                print(f"Processing curtain plot for model '{model_label}' and observation '{obs_label}'...")
+                        
+                                airplots.make_curtain_plot(
+                                    time=time_values, #time_as_float, #for datetime format on x-axis use time_values
+                                    altitude=altitude_values,
+                                    model_data_2d=model_data_2d_cleaned,
+                                    obs_data_2d=obs_data_2d,
+                                    model_var=modvar,
+                                    obs_var=obsvar,
+                                    grp_dict=curtain_config,
+                                    vmin=vmin,
+                                    vmax=vmax,
+                                    domain_type=domain_type,
+                                    domain_name=domain_name
+                                )
 
-                           
+                                # Save the scatter density plot for the current pair immediately
+                                outname_pair = f"{outname}_{obs_label}_vs_{model_label}.png"
+                                print(f"Saving Curtain plot (Model) with Scatter Overlay (Observation) to {outname_pair}...")  # Debugging print statement
+                                plt.savefig(outname_pair)
+                                plt.close()  # Close the current figure
+                                
+                        
+                            except Exception as e:
+                                print(f"Error generating curtain plot for {modvar} vs {obsvar}: {e}")
+
+                            
                                 
                         #qzr++ Added vertprofile plotype for aircraft vs model comparisons         
                         elif plot_type.lower() == 'vertprofile':
diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index 2ce533dc..428509e3 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -21,13 +21,14 @@
 from numpy import corrcoef
 sns.set_context('paper')
 from monet.plots.taylordiagram import TaylorDiagram as td
-from matplotlib.colors import TwoSlopeNorm, ListedColormap, LinearSegmentedColormap
+from matplotlib.colors import TwoSlopeNorm, ListedColormap, LinearSegmentedColormap, Normalize
 from matplotlib.patches import Rectangle
 from matplotlib.ticker import FuncFormatter
 
 import matplotlib.dates as mdates
-from matplotlib.colors import BoundaryNorm
-from matplotlib.colors import Normalize
+from scipy.interpolate import griddata
+
+
 
 from monet.util.tools import get_epa_region_bounds as get_epa_bounds 
 import math
@@ -227,6 +228,151 @@ def add_yax2_altitude(ax, pairdf, altitude_yax2, text_kwargs, vmin_y2, vmax_y2):
 
     return ax
 
+###NEW curtain plot qzr++  
+
+#### Function to interpolate and clean model data (FOR CURTAIN PLOT):
+
+def interpolate_and_clean_model_data(time, altitude, model_data):
+    """
+    Interpolates and cleans model data for creating curtain plots.
+
+    Parameters:
+    ----------
+    time : numpy.ndarray
+        Array of time points converted to floats (e.g., seconds since the start).
+    altitude : numpy.ndarray
+        Array of altitude points.
+    model_data : numpy.ndarray
+        2D array of model data values to be cleaned and interpolated.
+
+    Returns:
+    -------
+    numpy.ndarray
+        The cleaned and interpolated model data.
+    """
+    # Create meshgrid for interpolation based on the provided float time and altitude
+    time_grid, altitude_grid = np.meshgrid(np.unique(time), np.unique(altitude), indexing='ij')
+
+    # Flatten the grids and model data for interpolation
+    model_data_flat = model_data.flatten()
+    points = np.column_stack((time_grid.flatten(), altitude_grid.flatten()))
+
+    # Filter valid data points (non-NaN)
+    valid_mask = ~np.isnan(model_data_flat)
+    points_valid = points[valid_mask]
+    values_valid = model_data_flat[valid_mask]
+
+    # Perform interpolation using 'nearest' method to fill in NaNs
+    model_data_interpolated = griddata(points_valid, values_valid, points, method='nearest')
+
+    # Replace any remaining NaNs or Infs with appropriate limits
+    model_data_cleaned = np.nan_to_num(model_data_interpolated, nan=np.nanmean(model_data_interpolated), posinf=np.max(model_data_interpolated), neginf=np.min(model_data_interpolated))
+
+    # Reshape interpolated data back to the original grid shape
+    return model_data_cleaned.reshape(time_grid.shape)
+
+#### Main curtain plot function:
+def make_curtain_plot(time, altitude, model_data_2d, obs_data_2d, model_var, obs_var, grp_dict, vmin, vmax, domain_type, domain_name):
+    """
+    Generates a curtain plot comparing model data with observations across altitude over time,
+    with the ability to customize the appearance through a configuration dictionary.
+
+    Parameters:
+    ----------
+    time : numpy.ndarray
+        Array of time points, expected to be numerical values suitable for plotting.
+    altitude : numpy.ndarray
+        Array of altitude points.
+    model_data_2d : numpy.ndarray
+        2D array of cleaned and interpolated model data.
+    obs_data_2d : numpy.ndarray
+        2D array of observational data.
+    model_var : str
+        Model variable name for labeling.
+    obs_var : str
+        Observation variable name for labeling.
+    grp_dict : dict
+        Plot configuration options including aesthetics and normalization parameters.
+    vmin : float or None
+        Minimum value for color normalization, if applicable.
+    vmax : float or None
+        Maximum value for color normalization, if applicable.
+    domain_type : str
+        Type of domain being plotted (e.g., 'region', 'global').
+    domain_name : str
+        Name of the domain being plotted.
+
+    Returns:
+    -------
+    None
+        Saves the generated plot to a file.
+    """
+    # Create a figure with the specified dimensions
+    fig, ax = plt.subplots(figsize=grp_dict.get('fig_kwargs', {}).get('figsize', (12, 6)))
+
+    
+    # Calculate global min and max for color normalization
+    data_min = min(np.nanmin(model_data_2d), np.nanmin(obs_data_2d))
+    data_max = max(np.nanmax(model_data_2d), np.nanmax(obs_data_2d))
+
+    # Use custom colormap from YAML if provided, else default to 'viridis'
+    # Check if a custom colormap is specified or a named colormap should be used
+    color_levels = grp_dict.get('color_levels', 100)  # Default or from YAML
+    if grp_dict.get('color_map_custom', False):
+        # Use custom colors to make a colormap
+        colors = grp_dict['colors']
+        cmap = LinearSegmentedColormap.from_list("custom_cmap", colors, N=color_levels)
+    else:
+        # Use a named colormap
+        cmap = plt.get_cmap(grp_dict.get('color_map', 'viridis'))
+
+    norm = Normalize(vmin=data_min, vmax=data_max)
+
+    # Meshgrid for time and altitude to facilitate contour plotting
+    time_grid, altitude_grid = np.meshgrid(np.unique(time), np.unique(altitude), indexing='ij')    
+
+    # Example for contour plot:
+    contour = ax.contourf(time_grid, altitude_grid, model_data_2d, cmap=cmap, norm=norm)
+
+    # Observations overlayed on Model Curtain as Scatter
+    # Create a scatter plot using time and altitude
+    scatter = ax.scatter(time_grid.flatten(), altitude_grid.flatten(), c=obs_data_2d.flatten(), cmap=cmap, norm=norm, edgecolors='none')
+
+    cbar = plt.colorbar(contour, ax=ax)  #Colorbar for model and obs
+    # Get current font size from the 'text_kwargs' in grp_dict or use a default size
+    current_font_size = grp_dict.get('text_kwargs', {}).get('fontsize', 12) 
+    cbar.set_label(f'{model_var}', fontsize=current_font_size*0.6)  # You can adjust the fontsize as needed.
+    cbar.ax.tick_params(labelsize=current_font_size*0.6)  # Adjust the fontsize for the tick labels as needed.
+
+    # Set titles and labels with customized font settings
+    ax.set_title(f"Model Curtain Plot with Observations Scatter Overlaid: {model_var} vs {obs_var}", **grp_dict.get('text_kwargs', {'fontsize': 12}))
+    ax.set_xlabel('Time', **grp_dict.get('text_kwargs', {'fontsize': 12}))
+    ax.set_ylabel('Altitude', **grp_dict.get('text_kwargs', {'fontsize': 12}))
+    ##ax.tick_params(axis='both', which='major', labelsize=grp_dict.get('text_kwargs', {}).get('fontsize', 12))
+       
+    # Set tick parameters
+    ax.tick_params(axis='both', which='major', labelsize=current_font_size * 0.8)
+
+    # Format the x-axis to display dates and times
+    ax.xaxis_date()  # Interpret the x-axis values as dates
+    # Set format for date-time on the x-axis
+    ax.xaxis.set_major_formatter(mdates.DateFormatter('%Y-%m-%d %H:%M'))
+    fig.autofmt_xdate()  # Rotate date labels to prevent overlap
+
+    # Display legend and apply layout adjustments #Placeholder
+    # Plot legend
+    ##plt.legend(loc='upper left', bbox_to_anchor=(1.15, 1.15), borderaxespad=0.) #placeholder for any legend
+
+    # Adjust layout to fit the legend and prevent it from being cut off
+    plt.tight_layout(rect=[0, 0, 0.85, 1])
+
+    # Save the plot to a specified file path
+    plt.show()
+    ##plt.close() #close after showing to free up memory #if needed (placeholder)
+
+
+
+
 
 
 

From 684d8a1dc5972ad3618e4ca4c501516e0d4e4955 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Mon, 15 Apr 2024 23:50:30 -0600
Subject: [PATCH 04/20] Resolve Sphinx-Build CI error by adding hyphen(-) equal
 to the doctstring words such as, Parameters & Returns

---
 melodies_monet/plots/aircraftplots.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index 428509e3..5b76bebb 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -277,7 +277,7 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_data_2d, model_var, obs
     Generates a curtain plot comparing model data with observations across altitude over time,
     with the ability to customize the appearance through a configuration dictionary.
 
-    Parameters:
+    Parameters
     ----------
     time : numpy.ndarray
         Array of time points, expected to be numerical values suitable for plotting.
@@ -302,7 +302,7 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_data_2d, model_var, obs
     domain_name : str
         Name of the domain being plotted.
 
-    Returns:
+    Returns
     -------
     None
         Saves the generated plot to a file.

From 7edacd80cd2a60564ddaffe0cb95a9ef5810e645 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Mon, 15 Apr 2024 23:58:35 -0600
Subject: [PATCH 05/20] Fix Sphinx-Build related CI error for interpolate
 function docstring

---
 melodies_monet/plots/aircraftplots.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index 5b76bebb..afcc1fc1 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -236,7 +236,7 @@ def interpolate_and_clean_model_data(time, altitude, model_data):
     """
     Interpolates and cleans model data for creating curtain plots.
 
-    Parameters:
+    Parameters
     ----------
     time : numpy.ndarray
         Array of time points converted to floats (e.g., seconds since the start).
@@ -245,7 +245,7 @@ def interpolate_and_clean_model_data(time, altitude, model_data):
     model_data : numpy.ndarray
         2D array of model data values to be cleaned and interpolated.
 
-    Returns:
+    Returns
     -------
     numpy.ndarray
         The cleaned and interpolated model data.

From 0285d6a1f7057e042e90decde21249be43e456e6 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Tue, 11 Jun 2024 16:44:19 -0600
Subject: [PATCH 06/20] Subplots: 1) Model Curtain with Model scatter overlay
 and 2) Obs scatter

---
 ...elop_aircraft_testmultiplemodelgroups.yaml |  52 +++--
 melodies_monet/driver.py                      | 194 ++++++++++++++----
 melodies_monet/plots/aircraftplots.py         | 187 +++++++++--------
 3 files changed, 289 insertions(+), 144 deletions(-)

diff --git a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
index c29bbe0c..476a8438 100644
--- a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
+++ b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
@@ -6,12 +6,12 @@
 analysis:
   start_time: '2019-09-05-12:00:00' #UTC
   end_time: '2019-09-06-00:00:00' #UTC
-  output_dir: /wrk/qrasool/output_testaircraft_PR189_uptoviolin_multiplemodels #/wrk/charkins/melodies_monet/NOAA_CSL_main/test_pull_189_aircraft/output #/wrk/charkins/melodies_monet/aircraft/analysis #Opt if not specified plots stored in code directory.
+  output_dir: /wrk/qrasool/output_testaircraft_PR189_uptoviolin_multiplemodels_racmesrlvcpFTPnew  ##/wrk/qrasool/output_testaircraft_PR189_uptoviolin_multiplemodels #/wrk/charkins/melodies_monet/NOAA_CSL_main/test_pull_189_aircraft/output #/wrk/charkins/melodies_monet/aircraft/analysis #Opt if not specified plots stored in code directory.
   debug: True
 model:
   wrfchem_v4.2: # model label
     # files: /wrk/d2/rschwantes/wrf/firex_mech/qzhu/run_CONUS_fv19_BEIS_1.0xISO_RACM_v4.2.2_racm_berk_vcp_noI_phot_soa/0905/*
-    files: /wrk/qrasool/firex_mech_qzhu/run_CONUS_fv19_BEIS_1.0xISO_RACM_v4.2.2_racm_berk_vcp_noI_phot_soa/0905/* #/wrk/d2/rschwantes/wrf/firex_mech/qzhu/run_CONUS_fv19_BEIS_1.0xISO_RACM_v4.2.2_racm_berk_vcp_noI_phot_soa/0905/*
+    files: /wrk/qrasool/FIREX20190905wrfout_Newftp2024/racm_esrl_vcp/* #/wrk/d2/rschwantes/wrf/firex_mech/racm_esrl_vcp/*  (above path same as this) ##/wrk/qrasool/firex_mech_qzhu/run_CONUS_fv19_BEIS_1.0xISO_RACM_v4.2.2_racm_berk_vcp_noI_phot_soa/0905/* #/wrk/d2/rschwantes/wrf/firex_mech/qzhu/run_CONUS_fv19_BEIS_1.0xISO_RACM_v4.2.2_racm_berk_vcp_noI_phot_soa/0905/*
     mod_type: 'wrfchem'
     mod_kwargs: 
       mech: 'racm_esrl_vcp'
@@ -38,7 +38,7 @@ model:
       linestyle: ':' 
   wrfchem_v4.2_test: # model label
     # files: /wrk/d2/rschwantes/wrf/firex_mech/qzhu/run_CONUS_fv19_BEIS_1.0xISO_RACM_v4.2.2_racm_berk_vcp_noI_phot_soa/0905/*
-    files: /wrk/qrasool/firex_mech_qzhu/run_CONUS_fv19_BEIS_1.0xISO_RACM_v4.2.2_racm_berk_vcp_noI_phot_soa/0905/* #/wrk/d2/rschwantes/wrf/firex_mech/qzhu/run_CONUS_fv19_BEIS_1.0xISO_RACM_v4.2.2_racm_berk_vcp_noI_phot_soa/0905/*
+    files: /wrk/qrasool/FIREX20190905wrfout_Newftp2024/racm_esrl_vcp/* ##/wrk/qrasool/firex_mech_qzhu/run_CONUS_fv19_BEIS_1.0xISO_RACM_v4.2.2_racm_berk_vcp_noI_phot_soa/0905/* #/wrk/d2/rschwantes/wrf/firex_mech/qzhu/run_CONUS_fv19_BEIS_1.0xISO_RACM_v4.2.2_racm_berk_vcp_noI_phot_soa/0905/*
     mod_type: 'wrfchem'
     mod_kwargs: 
       mech: 'racm_esrl_vcp'
@@ -157,10 +157,10 @@ plots:
           color: g
         altitude_unit: m
         altitude_scaling_factor: 1
-        filter_dict: #Default is min and max if filter_dict doesn't define the options below (or if they are commented out)
-          altitude:
-            oper: "between"
-            value: [2000,8000] # values are [vmin_y2, vmax_y2] 
+        #filter_dict: #Default is min and max if filter_dict doesn't define the options below (or if they are commented out)
+         # altitude:
+          #  oper: "between"
+           # value: [2000,8000] # values are [vmin_y2, vmax_y2] 
             
   plot_grp2:
     type: 'vertprofile'  # plot type
@@ -239,9 +239,13 @@ plots:
   plot_grp5:
     type: 'curtain'
     fig_kwargs:
-      figsize: [12, 8]
+      figsize: [24,14] #[12, 8]
+    default_plot_kwargs:
+      #linewidth: 4.0
+      markersize: 30.
     text_kwargs:
       fontsize: 18
+      fontweight: 'bold'
     domain_type: ['all']
     domain_name: ['CONUS']
     data: ['firexaq_wrfchem_v4.2','firexaq_wrfchem_v4.2_test']
@@ -249,19 +253,35 @@ plots:
       rem_obs_nan: True
       set_axis: True
       ts_select_time: 'time'
-      #ts_avg_window: # 'H'  #Uncomment to use time averaging
+      ts_avg_window:  'H'  #Uncomment to use time averaging
       altitude_variable: 'altitude'
-      vmin_plot: 1000
-      vmax_plot: 10000
+      vmin_plot: #1000
+      vmax_plot: #10000
     # Optionally specify a custom colormap and number of levels
     ##If want to use standard matplotlib colorscale comment the True option and colors and uncomment the following:
-    ##color_map_custom: True
-    ##colors: ['blue','royalblue', 'cyan', 'yellow', 'orange', 'red'] ##['#440154', '#3b528b', '#21908d', '#5dc963', '#fde725'] #['purple', 'blue', 'green', 'yellow'] # Example gradient
-    color_map_custom: False
-    color_map: 'jet'  # This tells the plotting function to use the built-in 'rainbow' colormap
-    color_levels: 15 #10  # Define the number of distinct colors in the color bar
+    color_map_custom: False #True
+    colors: ["#ff8cff", "#dd6ff2", "#bb52e5", "#9935d8", "#7718cb",
+        "#0000bb", "#002ccc", "#0058dd", "#0084ee", "#00afff",
+        "#00ebff", "#27ffd7", "#63ff9b", "#a3ff5b", "#d3ff2b",
+        "#ffff00", "#ffcf00", "#ff9f00", "#ff6f00", "#ff3f00",
+        "#ff0000", "#d8000f", "#b2001f", "#8c002f", "#66003f",
+        "#343434", "#606060", "#8c8c8c", "#b8b8b8", "#e4e4e4"] #,
+        #"#ffffff"] #['blue','royalblue', 'cyan', 'yellow', 'orange', 'red'] ##['#440154', '#3b528b', '#21908d', '#5dc963', '#fde725'] #['purple', 'blue', 'green', 'yellow'] # Example gradient
+    ##use_interpolation: False #True #False  # Set to True for interpolated plot, False for block-like plot #Obsolete NOW
+    color_map: 'Spectral_r' #'jet'  # This tells the plotting function to use the built-in 'rainbow' colormap
+    color_levels: 30 #31 #15 #10  # Define the number of distinct colors in the color bar
     # Optionally specify dynamic color bar ranges based on data
     dynamic_color_range: True  # New key to use dynamic ranges based on model and obs data
+    variable_limits:
+      NO2_CL_RYERSON:
+        vmin: 0
+        vmax: 1 #10
+      NO_CL_RYERSON:
+        vmin: 0
+        vmax: 1 #15
+      O3_CL_RYERSON:
+        vmin: 0
+        vmax: 100
     
   plot_grp6:
     type: 'taylor' # plot type
diff --git a/melodies_monet/driver.py b/melodies_monet/driver.py
index 6a7a7f93..67655e5f 100644
--- a/melodies_monet/driver.py
+++ b/melodies_monet/driver.py
@@ -12,6 +12,7 @@
 import numpy as np
 import datetime
 
+
 from .util import write_util
 
 __all__ = (
@@ -1186,6 +1187,11 @@ def plotting(self):
         -------
         None
         """
+        import stratify
+        from .util.tools import resample_stratify
+        from .util.tools import vert_interp
+        import matplotlib.dates as mdates
+        import numpy as np
         import matplotlib.pyplot as plt
         pair_keys = list(self.paired.keys())
         if self.paired[pair_keys[0]].type.lower() in ['sat_grid_clm','sat_swath_clm']:
@@ -1245,6 +1251,7 @@ def plotting(self):
                     # Then loop through each of the pairs to add to the plot.
                     for p_index, p_label in enumerate(pair_labels):
                         p = self.paired[p_label]
+                        
                         # find the pair model label that matches the obs var
                         index = p.obs_vars.index(obsvar)
                         modvar = p.model_vars[index]
@@ -1567,60 +1574,173 @@ def plotting(self):
                                   ##  ax = airplots.add_yax2_altitude(ax, pairdf, altitude_variable, altitude_ticks, text_kwargs)
                                 ##savefig(outname + '.png', logo_height=150)
                                 ##del (ax, fig_dict, plot_dict, text_dict, obs_dict, obs_plot_dict) #Clear axis for next plot.
-
                         
+
+
+
                         elif plot_type.lower() == 'curtain':
                             curtain_config = grp_dict
-                            vmin = curtain_config['data_proc'].get('vmin_plot', None)
-                            vmax = curtain_config['data_proc'].get('vmax_plot', None)
+                        
+                            # Ensure we use the correct observation and model objects from pairing
+                            obs = self.obs[p.obs]
+                            mod = self.models[p.model]  # Get the model associated with this pair
+                            model_obj = mod.obj  # This assumes `mod.obj` contains the necessary model data
+                            model_data = model_obj[modvar].values
+                            obs_data = obs.obj[obsvar].values
+                        
+                            # Debugging: print dimensions
+                            ##print(f"Dimensions of model data: {model_data.shape}, ndims: {model_data.ndim}")
+                            ##print(f"Dimensions of observation data: {obs_data.shape}, ndims: {obs_data.ndim}")
+                        
+                            if not isinstance(obs.obj, pd.DataFrame):
+                                obs.obj = obs.obj.to_dataframe()
+                        
+                            # Drop any variables where coords NaN
+                            obs.obj = obs.obj.reset_index().dropna(subset=['pressure_obs', 'latitude', 'longitude']).set_index('time')
+                        
+                            # Convert to get something useful for MONET
+                            new_ds_obs = obs.obj.rename_axis('time_obs').reset_index().monet._df_to_da().set_coords(['time_obs', 'pressure_obs'])
+                        
+                            # Nearest neighbor approach to find closest grid cell to each point
+                            ds_model = m.util.combinetool.combine_da_to_da(model_obj, new_ds_obs, merge=False)
+                        
+                            # Interpolate based on time in the observations
+                            ds_model = ds_model.interp(time=ds_model.time_obs.squeeze())
+                        
+                            # Fetch the observation configuration for colorbar labels and plot limits
+                            obs_plot_dict = self.control_dict['obs'][p.obs]['variables'][obsvar]
+                        
+                            # Print ds_model and pressure_model values
+                            ##print(f"ds_model: {ds_model}")
+                            ##print(f"pressure_model values: {ds_model['pressure_model'].values}")
+                        
+                            # Define target pressures for interpolation based on the range of pressure_model
+                            min_pressure = ds_model['pressure_model'].min().compute()
+                            max_pressure = ds_model['pressure_model'].max().compute()
+                            interval = 100  # Interval in Pa
+                        
+                            print(f"Pressure MIN:{min_pressure}, max: {max_pressure}, interval: {interval}")
+                        
+                            target_pressures = np.linspace(max_pressure, min_pressure, interval)
                             
-                            model_label = p.model
-                            obs_label = p.obs
+                            # Debugging: print target pressures
+                            ##print(f"Generated target pressures: {target_pressures}, shape: {target_pressures.shape}")
                         
-                            if 'altitude' in pairdf.columns and vmin is not None and vmax is not None:
-                                pairdf_filtered = pairdf[(pairdf['altitude'] >= vmin) & (pairdf['altitude'] <= vmax)]
-                            else:
-                                pairdf_filtered = pairdf
+                            # Check for NaN values before interpolation
+                            ##print(f"NaNs in model_data before interpolation: {np.isnan(ds_model[modvar]).sum().compute()}")
+                            ##print(f"NaNs in pressure_model before interpolation: {np.isnan(ds_model['pressure_model']).sum().compute()}")
                         
-                            try:
-                                model_data_pivot = pairdf_filtered.pivot(index='time', columns='altitude', values=modvar)
-                                obs_data_pivot = pairdf_filtered.pivot(index='time', columns='altitude', values=obsvar)
+                            # Resample model data to target pressures using stratify
+                            da_wrf_const = resample_stratify(ds_model[modvar], target_pressures, ds_model['pressure_model'], axis=1, interpolation='linear', extrapolation='nan')
+                            da_wrf_const.name = modvar
                         
-                                # Convert datetime index to float seconds since the first timestamp
-                                time_values = model_data_pivot.index.values
-                                time_as_float = (time_values - time_values[0]).astype('timedelta64[s]').astype(float)
-                                
-                                altitude_values = model_data_pivot.columns.values
-                                model_data_2d = model_data_pivot.values
-                                obs_data_2d = obs_data_pivot.values
+                            # Create target_pressures DataArray
+                            da_target_pressures = xr.DataArray(target_pressures, dims=('z'))
+                            da_target_pressures.name = 'target_pressures'
+                        
+                            # Merge DataArrays into a single Dataset
+                            ds_wrf_const = xr.merge([da_wrf_const, da_target_pressures])
+                            ds_wrf_const = ds_wrf_const.set_coords('target_pressures')
+                        
+                            # Debugging: print merged dataset
+                            ##print(ds_wrf_const)
+                        
+                            # Ensure model_data_2d is properly reshaped for the contourf plot
+                            model_data_2d = ds_wrf_const[modvar].squeeze()
                         
-                                model_data_2d_cleaned = airplots.interpolate_and_clean_model_data(time_as_float, altitude_values, model_data_2d)
-                                #Create Curtain Plots
-                                print(f"Processing curtain plot for model '{model_label}' and observation '{obs_label}'...")
+                            # Debugging: print reshaped model data shape
+                            ##print(f"Reshaped model data shape: {model_data_2d.shape}")
+                        
+                            # Align the observation data with the time array length for scatter plot
+                            time = ds_wrf_const['time'].values
+                            obs_pressure = new_ds_obs['pressure_obs'].values.flatten()
+                            obs_data_2d = obs_data.flatten()
+                        
+                            ##print(f"time shape: {time.shape}, ndims: {time.ndim}")
+                            ##print(f"obs_pressure shape: {obs_pressure.shape}, ndims: {obs_pressure.ndim}")
+                            ##print(f"obs_data_2d shape: {obs_data_2d.shape}, ndims: {obs_data_2d.ndim}")
+                        
+                            time_dates = mdates.date2num(pd.to_datetime(time))
+                        
+                            if len(obs_pressure) != len(time_dates):
+                                obs_pressure = np.interp(time_dates, mdates.date2num(pd.to_datetime(new_ds_obs['time_obs'].values)), obs_pressure)
+                                obs_data_2d = np.interp(time_dates, mdates.date2num(pd.to_datetime(new_ds_obs['time_obs'].values)), obs_data_2d)
+                                ##print(f"Length of obs_pressure after interpolation: {len(obs_pressure)}")
+                                ##print(f"Length of obs_data_2d after interpolation: {len(obs_data_2d)}")
+                        
+                            # Fetch the observation configuration for colorbar labels
+                            obs_label_config = self.control_dict['obs'][p.obs]['variables']
+                            
+                            # Ensure the observation data is a DataFrame and reset index
+                            df_obs = new_ds_obs.to_dataframe().reset_index()
+                            df_obs['time'] = df_obs['time_obs']
+                        
+                            # Create var_name_list dynamically based on modvar and pressure_model
+                            var_name_list = [modvar, 'pressure_model'] if modvar != 'pressure_model' else [modvar]
+                            available_vars = [var for var in var_name_list if var in ds_model]
+                            
+                            # Print diagnostics for available variables
+                            ##print(f"Available variables for interpolation: {available_vars}")
+                            
+                            # Convert ds_model to DataFrame and reset index
+                            df_model = ds_model[available_vars].to_dataframe().reset_index()
+                            
+                            # Print diagnostics for the observation DataFrame
+                            ##print(f"Type of df_obs: {type(df_obs)}")
+                            ##print(f"Head of df_obs:\n{df_obs.head()}")
+                            
+                            # Print diagnostics for the model DataFrame
+                            ##print(f"Type of df_model: {type(df_model)}")
+                            ##print(f"Head of df_model:\n{df_model.head()}")
+                            
+                            # Ensure both DataFrames have the necessary columns for merging
+                            common_columns = set(df_obs.columns) & set(df_model.columns)
+                            ##print(f"Common columns between observation and model DataFrames: {common_columns}")
+                            if not {'latitude', 'longitude', 'pressure_obs', 'time'}.issubset(common_columns):
+                                raise ValueError("Missing necessary columns for merging. Columns needed: 'latitude', 'longitude', 'pressure_obs', 'time'")
+                            
+                            # Perform vertical interpolation only for available variables
+                            df_wrf = vert_interp(ds_model, df_obs, available_vars)
+                            
+                            # Print diagnostics for the resulting DataFrame
+                            ##print(f"Type of df_wrf: {type(df_wrf)}")
+                            ##print(f"Head of df_wrf:\n{df_wrf.head()}")
+                            
+                            # Ensure the modvar exists in the resulting DataFrame
+                            if modvar not in df_wrf.columns:
+                                raise ValueError(f"{modvar} not found in the resulting DataFrame columns")
+                            
+                            model_scatter_data = df_wrf[modvar].values
+                            
+                            # Print diagnostics for model scatter data
+                            ##print(f"model_scatter_data shape: {model_scatter_data.shape}, ndims: {model_scatter_data.ndim}")
+                        
+                            # Generate the curtain plot using airplots.make_curtain_plot
+                            try:
+                                ##print(f"Length of time_dates: {len(time_dates)}")
+                                ##print(f"Length of obs_pressure: {len(obs_pressure)}")
+                                ##print(f"Length of obs_data_2d: {len(obs_data_2d)}")
+                                ##print(f"model_scatter_data shape: {model_scatter_data.shape}")
                         
                                 airplots.make_curtain_plot(
-                                    time=time_values, #time_as_float, #for datetime format on x-axis use time_values
-                                    altitude=altitude_values,
-                                    model_data_2d=model_data_2d_cleaned,
-                                    obs_data_2d=obs_data_2d,
+                                    time=pd.to_datetime(time),
+                                    altitude=target_pressures,  # Use target_pressures for interpolation
+                                    model_data_2d=model_data_2d,  # Already reshaped to match the expected shape
+                                    obs_pressure=obs_pressure,  # Pressure_obs for obs scatter plot
+                                    obs_data_2d=obs_data_2d,  # Use original observation data for scatter plot
                                     model_var=modvar,
                                     obs_var=obsvar,
                                     grp_dict=curtain_config,
-                                    vmin=vmin,
-                                    vmax=vmax,
                                     domain_type=domain_type,
-                                    domain_name=domain_name
+                                    domain_name=domain_name,
+                                    obs_label_config=obs_label_config,
+                                    model_scatter_data=model_scatter_data  # Pass model scatter data for plotting
                                 )
-
-                                # Save the scatter density plot for the current pair immediately
-                                outname_pair = f"{outname}_{obs_label}_vs_{model_label}.png"
-                                print(f"Saving Curtain plot (Model) with Scatter Overlay (Observation) to {outname_pair}...")  # Debugging print statement
-                                plt.savefig(outname_pair)
-                                plt.close()  # Close the current figure
-                                
-                        
                             except Exception as e:
                                 print(f"Error generating curtain plot for {modvar} vs {obsvar}: {e}")
+                            finally:
+                                plt.close('all')  # Clean up matplotlib resources
+
 
                             
                                 
diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index afcc1fc1..4d3a3707 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -21,12 +21,13 @@
 from numpy import corrcoef
 sns.set_context('paper')
 from monet.plots.taylordiagram import TaylorDiagram as td
-from matplotlib.colors import TwoSlopeNorm, ListedColormap, LinearSegmentedColormap, Normalize
+from matplotlib.colors import TwoSlopeNorm, ListedColormap, LinearSegmentedColormap, Normalize, LogNorm
 from matplotlib.patches import Rectangle
 from matplotlib.ticker import FuncFormatter
 
+from matplotlib.dates import DateFormatter
 import matplotlib.dates as mdates
-from scipy.interpolate import griddata
+#from scipy.interpolate import griddata
 
 
 
@@ -228,54 +229,19 @@ def add_yax2_altitude(ax, pairdf, altitude_yax2, text_kwargs, vmin_y2, vmax_y2):
 
     return ax
 
-###NEW curtain plot qzr++  
 
-#### Function to interpolate and clean model data (FOR CURTAIN PLOT):
 
-def interpolate_and_clean_model_data(time, altitude, model_data):
-    """
-    Interpolates and cleans model data for creating curtain plots.
-
-    Parameters
-    ----------
-    time : numpy.ndarray
-        Array of time points converted to floats (e.g., seconds since the start).
-    altitude : numpy.ndarray
-        Array of altitude points.
-    model_data : numpy.ndarray
-        2D array of model data values to be cleaned and interpolated.
 
-    Returns
-    -------
-    numpy.ndarray
-        The cleaned and interpolated model data.
+###NEW curtain plot qzr++  (NEW CURTAIN model plot with model overlay, shared x-axis 
+def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d, model_var, obs_var, grp_dict, vmin=None, vmax=None, domain_type=None, domain_name=None, obs_label_config=None, model_scatter_data=None):
     """
-    # Create meshgrid for interpolation based on the provided float time and altitude
-    time_grid, altitude_grid = np.meshgrid(np.unique(time), np.unique(altitude), indexing='ij')
-
-    # Flatten the grids and model data for interpolation
-    model_data_flat = model_data.flatten()
-    points = np.column_stack((time_grid.flatten(), altitude_grid.flatten()))
-
-    # Filter valid data points (non-NaN)
-    valid_mask = ~np.isnan(model_data_flat)
-    points_valid = points[valid_mask]
-    values_valid = model_data_flat[valid_mask]
-
-    # Perform interpolation using 'nearest' method to fill in NaNs
-    model_data_interpolated = griddata(points_valid, values_valid, points, method='nearest')
-
-    # Replace any remaining NaNs or Infs with appropriate limits
-    model_data_cleaned = np.nan_to_num(model_data_interpolated, nan=np.nanmean(model_data_interpolated), posinf=np.max(model_data_interpolated), neginf=np.min(model_data_interpolated))
-
-    # Reshape interpolated data back to the original grid shape
-    return model_data_cleaned.reshape(time_grid.shape)
-
-#### Main curtain plot function:
-def make_curtain_plot(time, altitude, model_data_2d, obs_data_2d, model_var, obs_var, grp_dict, vmin, vmax, domain_type, domain_name):
-    """
-    Generates a curtain plot comparing model data with observations across altitude over time,
+    Generates a curtain plot comparing model data with obs across altitude (Pressure, right now) over time,
     with the ability to customize the appearance through a configuration dictionary.
+    ##Two Subplots: 1) model data contourf plot with model scaatter overlay, 
+    2) another for observation data using scatter plot. 
+    This layout ensures that both datasets can be analyzed without visual interference from each other.
+    Shared X-Axis: The time axis is shared between the two plots for better comparison.
+    Titles and Labels: Each subplot has a title specific to the data it displays.
 
     Parameters
     ----------
@@ -285,6 +251,8 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_data_2d, model_var, obs
         Array of altitude points.
     model_data_2d : numpy.ndarray
         2D array of cleaned and interpolated model data.
+    obs_pressure : numpy.ndarray
+        Array of pressure points corresponding to observations.
     obs_data_2d : numpy.ndarray
         2D array of observational data.
     model_var : str
@@ -301,77 +269,114 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_data_2d, model_var, obs
         Type of domain being plotted (e.g., 'region', 'global').
     domain_name : str
         Name of the domain being plotted.
+    obs_label_config : dict
+        Configuration dictionary for observation labels.
+    model_scatter_data : numpy.ndarray or None
+        Data for model scatter plot overlay.
 
     Returns
     -------
     None
         Saves the generated plot to a file.
     """
-    # Create a figure with the specified dimensions
-    fig, ax = plt.subplots(figsize=grp_dict.get('fig_kwargs', {}).get('figsize', (12, 6)))
+    if model_data_2d.size == 0 or obs_data_2d.size == 0:
+        print("Warning: Model or observation data is empty. Skipping plot.")
+        return
 
+    # Determine vmin and vmax from YAML config if provided
+    if 'variable_limits' in grp_dict:
+        var_limits = grp_dict['variable_limits'].get(obs_var, {})
+        vmin = var_limits.get('vmin', vmin)
+        vmax = var_limits.get('vmax', vmax)
+
+    # Debugging: print shapes and ndims
+    print(f"time shape: {time.shape}, ndims: {time.ndim}")
+    print(f"altitude shape: {altitude.shape}, ndims: {altitude.ndim}")
+    print(f"model_data_2d shape: {model_data_2d.shape}, ndims: {model_data_2d.ndim}")
+    print(f"obs_data_2d shape: {obs_data_2d.shape}, ndims: {obs_data_2d.ndim}")
+
+    # Determine vmin and vmax dynamically if not provided
+    if vmin is None or vmax is None:
+        vmin = min(np.nanmin(model_data_2d), np.nanmin(obs_data_2d))
+        vmax = max(np.nanmax(model_data_2d), np.nanmax(obs_data_2d))
     
-    # Calculate global min and max for color normalization
-    data_min = min(np.nanmin(model_data_2d), np.nanmin(obs_data_2d))
-    data_max = max(np.nanmax(model_data_2d), np.nanmax(obs_data_2d))
+    # Set colorbar min to zero if the minimum value is less than zero
+    if vmin < 0:
+        vmin = 0
+
+    print(f"vmin and vmax set dynamically based on model and observation data to {vmin} and {vmax}")
 
-    # Use custom colormap from YAML if provided, else default to 'viridis'
-    # Check if a custom colormap is specified or a named colormap should be used
-    color_levels = grp_dict.get('color_levels', 100)  # Default or from YAML
-    if grp_dict.get('color_map_custom', False):
-        # Use custom colors to make a colormap
+    fig, axs = plt.subplots(nrows=2, figsize=grp_dict.get('fig_kwargs', {}).get('figsize', (20, 8)), sharex=True, gridspec_kw={'height_ratios': [1, 1]})
+
+    # Handling custom color maps
+    if grp_dict.get('color_map_custom', False) and 'colors' in grp_dict:
         colors = grp_dict['colors']
-        cmap = LinearSegmentedColormap.from_list("custom_cmap", colors, N=color_levels)
+        cmap = LinearSegmentedColormap.from_list("custom_cmap", colors, N=grp_dict.get('color_levels', 100))
     else:
-        # Use a named colormap
         cmap = plt.get_cmap(grp_dict.get('color_map', 'viridis'))
 
-    norm = Normalize(vmin=data_min, vmax=data_max)
+    norm = Normalize(vmin=vmin, vmax=vmax)
 
-    # Meshgrid for time and altitude to facilitate contour plotting
-    time_grid, altitude_grid = np.meshgrid(np.unique(time), np.unique(altitude), indexing='ij')    
+    time_dates = mdates.date2num(time)
+    time_mesh, altitude_mesh = np.meshgrid(time_dates, altitude, indexing='ij')
 
-    # Example for contour plot:
-    contour = ax.contourf(time_grid, altitude_grid, model_data_2d, cmap=cmap, norm=norm)
+    levels = np.linspace(vmin, vmax, 200)
+    contourf_plot = axs[0].contourf(time_mesh, altitude_mesh, model_data_2d, levels=levels, cmap=cmap, norm=norm, extend='both')
 
-    # Observations overlayed on Model Curtain as Scatter
-    # Create a scatter plot using time and altitude
-    scatter = ax.scatter(time_grid.flatten(), altitude_grid.flatten(), c=obs_data_2d.flatten(), cmap=cmap, norm=norm, edgecolors='none')
+    # Adjust the layout to place the colorbar at the bottom
+    cbar = fig.colorbar(contourf_plot, ax=axs.ravel().tolist(), orientation='horizontal', pad=0.1, aspect=50)
+    cbar.set_label(obs_label_config.get(obs_var, {}).get('ylabel_plot', 'Variable (units)'), fontsize=12, fontweight='bold')
+    cbar.ax.tick_params(labelsize=10)
 
-    cbar = plt.colorbar(contour, ax=ax)  #Colorbar for model and obs
-    # Get current font size from the 'text_kwargs' in grp_dict or use a default size
-    current_font_size = grp_dict.get('text_kwargs', {}).get('fontsize', 12) 
-    cbar.set_label(f'{model_var}', fontsize=current_font_size*0.6)  # You can adjust the fontsize as needed.
-    cbar.ax.tick_params(labelsize=current_font_size*0.6)  # Adjust the fontsize for the tick labels as needed.
+    axs[0].xaxis.set_major_formatter(mdates.DateFormatter('%m-%d %H'))
+    axs[0].xaxis.set_major_locator(mdates.AutoDateLocator())
+    fig.autofmt_xdate(rotation=45, ha='right')
 
-    # Set titles and labels with customized font settings
-    ax.set_title(f"Model Curtain Plot with Observations Scatter Overlaid: {model_var} vs {obs_var}", **grp_dict.get('text_kwargs', {'fontsize': 12}))
-    ax.set_xlabel('Time', **grp_dict.get('text_kwargs', {'fontsize': 12}))
-    ax.set_ylabel('Altitude', **grp_dict.get('text_kwargs', {'fontsize': 12}))
-    ##ax.tick_params(axis='both', which='major', labelsize=grp_dict.get('text_kwargs', {}).get('fontsize', 12))
-       
-    # Set tick parameters
-    ax.tick_params(axis='both', which='major', labelsize=current_font_size * 0.8)
+    #axs[0].set_title(f"Model vs Observation Curtain Plot: {model_var} vs {obs_var}", fontsize=12, fontweight='bold')
+    # Set the main title and subplot titles
+    fig.suptitle(f"Model vs Observation Curtain Plot: {model_var} vs {obs_var}", fontsize=16, fontweight='bold')
+    axs[0].set_title("Model Curtain with Model Scatter Overlay", fontsize=12, fontweight='bold')
+    
+    axs[0].set_ylabel('Pressure (Pa)', fontsize=12, fontweight='bold')
+    axs[0].invert_yaxis()  # Invert y-axis to have max pressure at the bottom
+    axs[0].tick_params(axis='both', labelsize=10)
 
-    # Format the x-axis to display dates and times
-    ax.xaxis_date()  # Interpret the x-axis values as dates
-    # Set format for date-time on the x-axis
-    ax.xaxis.set_major_formatter(mdates.DateFormatter('%Y-%m-%d %H:%M'))
-    fig.autofmt_xdate()  # Rotate date labels to prevent overlap
+    # Separate subplot for the observation scatter plot
+    axs[1].set_title("Observation Scatter", fontsize=12, fontweight='bold')
+    scatter = axs[1].scatter(time_dates, obs_pressure, c=obs_data_2d, cmap=cmap, norm=norm, edgecolor=(0, 0, 0, 0), linewidth=0.3, alpha=0.5)
+    axs[1].xaxis.set_major_formatter(mdates.DateFormatter('%m-%d %H'))
+    axs[1].xaxis.set_major_locator(mdates.AutoDateLocator())
+    fig.autofmt_xdate(rotation=45, ha='right')
 
-    # Display legend and apply layout adjustments #Placeholder
-    # Plot legend
-    ##plt.legend(loc='upper left', bbox_to_anchor=(1.15, 1.15), borderaxespad=0.) #placeholder for any legend
+    axs[1].set_ylabel('Pressure (Pa)', fontsize=12, fontweight='bold')
+    axs[1].invert_yaxis()  # Invert y-axis to have max pressure at the bottom
+    axs[1].tick_params(axis='both', labelsize=10)
+    axs[1].set_xlabel('Time', fontsize=12, fontweight='bold')
 
-    # Adjust layout to fit the legend and prevent it from being cut off
-    plt.tight_layout(rect=[0, 0, 0.85, 1])
+    # Set the same y-tick labels for both subplots and exclude the minimum value
+    y_ticks = axs[0].get_yticks()
+    axs[1].set_yticks(y_ticks)
+    axs[1].set_yticklabels([str(int(tick)) if tick != y_ticks.min() else '' for tick in y_ticks])
 
-    # Save the plot to a specified file path
-    plt.show()
-    ##plt.close() #close after showing to free up memory #if needed (placeholder)
+    # Plot model scatter data on top of the model curtain plot
+    if model_scatter_data is not None:
+        axs[0].scatter(time_dates, obs_pressure, c=model_scatter_data, cmap=cmap, norm=norm, edgecolor='black', linewidth=0.5, alpha=0.7)
+
+        # Print diagnostics for model scatter data
+        print(f"model_scatter_data shape: {model_scatter_data.shape}, ndims: {model_scatter_data.ndim}")
 
+    plt.show()
 
 
+    #Diagnostic Histogram
+    #plt.figure(figsize=(10, 4))
+    #plt.hist(model_data_2d.values.flatten(), bins=100, alpha=0.15, color='blue', label='Model', range=(vmin, vmax))
+    #plt.hist(obs_data_2d.flatten(), bins=100, alpha=0.5, color='green', label='Observation', range=(vmin, vmax))
+    #plt.legend()
+    #plt.title('Distribution of Model and Observation Data')
+    #plt.xlabel('Value')
+    #plt.ylabel('Frequency')
+    #plt.show()
 
 
 

From de7aff9d00b373ab12959385bff80e32c4a952e1 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Mon, 22 Jul 2024 18:06:36 -0600
Subject: [PATCH 07/20] Some Fixes for curtain plot YAML and other minor fixes
 in driver and aircraftplots

---
 ...elop_aircraft_testmultiplemodelgroups.yaml | 34 ++++----
 melodies_monet/driver.py                      | 36 ++++++++-
 melodies_monet/plots/aircraftplots.py         | 77 ++++++++++++-------
 3 files changed, 96 insertions(+), 51 deletions(-)

diff --git a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
index 476a8438..47906b47 100644
--- a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
+++ b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
@@ -239,12 +239,12 @@ plots:
   plot_grp5:
     type: 'curtain'
     fig_kwargs:
-      figsize: [24,14] #[12, 8]
+      figsize: [20,14] #[12, 8]
     default_plot_kwargs:
       #linewidth: 4.0
-      markersize: 30.
+      markersize: 40.
     text_kwargs:
-      fontsize: 18
+      fontsize: 25 #18
       fontweight: 'bold'
     domain_type: ['all']
     domain_name: ['CONUS']
@@ -252,14 +252,9 @@ plots:
     data_proc:
       rem_obs_nan: True
       set_axis: True
-      ts_select_time: 'time'
-      ts_avg_window:  'H'  #Uncomment to use time averaging
       altitude_variable: 'altitude'
-      vmin_plot: #1000
-      vmax_plot: #10000
-    # Optionally specify a custom colormap and number of levels
-    ##If want to use standard matplotlib colorscale comment the True option and colors and uncomment the following:
-    color_map_custom: False #True
+    ##To use a custom Matplotlib colormap, set “color_map_custom”=True and specify “colors” and "color_levels" like the example below. To use a standard Matplotlib colormap, set “color_map_custom” = False and specify a “color_map”:
+    color_map_custom: True #False #True
     colors: ["#ff8cff", "#dd6ff2", "#bb52e5", "#9935d8", "#7718cb",
         "#0000bb", "#002ccc", "#0058dd", "#0084ee", "#00afff",
         "#00ebff", "#27ffd7", "#63ff9b", "#a3ff5b", "#d3ff2b",
@@ -267,21 +262,20 @@ plots:
         "#ff0000", "#d8000f", "#b2001f", "#8c002f", "#66003f",
         "#343434", "#606060", "#8c8c8c", "#b8b8b8", "#e4e4e4"] #,
         #"#ffffff"] #['blue','royalblue', 'cyan', 'yellow', 'orange', 'red'] ##['#440154', '#3b528b', '#21908d', '#5dc963', '#fde725'] #['purple', 'blue', 'green', 'yellow'] # Example gradient
-    ##use_interpolation: False #True #False  # Set to True for interpolated plot, False for block-like plot #Obsolete NOW
-    color_map: 'Spectral_r' #'jet'  # This tells the plotting function to use the built-in 'rainbow' colormap
     color_levels: 30 #31 #15 #10  # Define the number of distinct colors in the color bar
-    # Optionally specify dynamic color bar ranges based on data
-    dynamic_color_range: True  # New key to use dynamic ranges based on model and obs data
+    color_map: 'Spectral_r' #'jet'  # This tells the plotting function to use the built-in 'rainbow' colormap   
     variable_limits:
       NO2_CL_RYERSON:
-        vmin: 0
-        vmax: 1 #10
+        cmin: 0
+        cmax: 1 #10
       NO_CL_RYERSON:
-        vmin: 0
-        vmax: 1 #15
+        cmin: 0
+        cmax: 1 #15
       O3_CL_RYERSON:
-        vmin: 0
-        vmax: 100
+        cmin: 0
+        cmax: 100
+    vmin: 120000 #hPa  #0  #Optional (y-axis limits)
+    vmax: 5000 # Optional #need to be edited as per min and max altitude (i.e., vmin and vmax)
     
   plot_grp6:
     type: 'taylor' # plot type
diff --git a/melodies_monet/driver.py b/melodies_monet/driver.py
index 67655e5f..1a867b32 100644
--- a/melodies_monet/driver.py
+++ b/melodies_monet/driver.py
@@ -1187,7 +1187,7 @@ def plotting(self):
         -------
         None
         """
-        import stratify
+        
         from .util.tools import resample_stratify
         from .util.tools import vert_interp
         import matplotlib.dates as mdates
@@ -1579,7 +1579,28 @@ def plotting(self):
 
 
                         elif plot_type.lower() == 'curtain':
+                            if set_yaxis:
+                                if all(k in obs_plot_dict for k in ('vmin_plot', 'vmax_plot')):
+                                    vmin = obs_plot_dict['vmin_plot']
+                                    vmax = obs_plot_dict['vmax_plot']
+                                else:
+                                    print('Warning: vmin_plot and vmax_plot not specified for ' + obsvar + ', so default used.')
+                                    vmin = None
+                                    vmax = None
+                            else:
+                                vmin = None
+                                vmax = None
                             curtain_config = grp_dict
+
+                            # Determine cmin and cmax (colorbar min max) from YAML config if provided
+                            if 'variable_limits' in curtain_config:
+                                var_limits = curtain_config['variable_limits'].get(obsvar, {})
+                                cmin = var_limits.get('cmin')
+                                cmax = var_limits.get('cmax')
+
+                            # Inside your loop for processing each pair
+                            obs_label = p.obs
+                            model_label = p.model
                         
                             # Ensure we use the correct observation and model objects from pairing
                             obs = self.obs[p.obs]
@@ -1730,12 +1751,23 @@ def plotting(self):
                                     obs_data_2d=obs_data_2d,  # Use original observation data for scatter plot
                                     model_var=modvar,
                                     obs_var=obsvar,
+                                    cmin=cmin,
+                                    cmax=cmax,
+                                    vmin=vmin,
+                                    vmax=vmax,
+                                    plot_dict=plot_dict,
                                     grp_dict=curtain_config,
                                     domain_type=domain_type,
                                     domain_name=domain_name,
                                     obs_label_config=obs_label_config,
-                                    model_scatter_data=model_scatter_data  # Pass model scatter data for plotting
+                                    model_scatter_data=model_scatter_data,  # Pass model scatter data for plotting
+                                    debug=self.debug  # Pass debug flag
                                 )
+                            
+                                # Save the curtain plot for the current pair immediately
+                                outname_pair = f"{outname}_{obs_label}_vs_{model_label}.png"
+                                print(f"Saving curtain plot to {outname_pair}...")
+                                plt.savefig(outname_pair)
                             except Exception as e:
                                 print(f"Error generating curtain plot for {modvar} vs {obsvar}: {e}")
                             finally:
diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index 4d3a3707..2fa49da8 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -21,13 +21,12 @@
 from numpy import corrcoef
 sns.set_context('paper')
 from monet.plots.taylordiagram import TaylorDiagram as td
-from matplotlib.colors import TwoSlopeNorm, ListedColormap, LinearSegmentedColormap, Normalize, LogNorm
+from matplotlib.colors import TwoSlopeNorm, ListedColormap, LinearSegmentedColormap, Normalize
 from matplotlib.patches import Rectangle
 from matplotlib.ticker import FuncFormatter
 
 from matplotlib.dates import DateFormatter
 import matplotlib.dates as mdates
-#from scipy.interpolate import griddata
 
 
 
@@ -233,11 +232,11 @@ def add_yax2_altitude(ax, pairdf, altitude_yax2, text_kwargs, vmin_y2, vmax_y2):
 
 
 ###NEW curtain plot qzr++  (NEW CURTAIN model plot with model overlay, shared x-axis 
-def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d, model_var, obs_var, grp_dict, vmin=None, vmax=None, domain_type=None, domain_name=None, obs_label_config=None, model_scatter_data=None):
+def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d, model_var, obs_var, grp_dict, vmin=None, vmax=None, cmin=None, cmax=None, plot_dict=None, domain_type=None, domain_name=None, obs_label_config=None, model_scatter_data=None, debug=False):
     """
     Generates a curtain plot comparing model data with obs across altitude (Pressure, right now) over time,
     with the ability to customize the appearance through a configuration dictionary.
-    ##Two Subplots: 1) model data contourf plot with model scaatter overlay, 
+    ##Two Subplots: 1) model data contourf plot with model scatter overlay, 
     2) another for observation data using scatter plot. 
     This layout ensures that both datasets can be analyzed without visual interference from each other.
     Shared X-Axis: The time axis is shared between the two plots for better comparison.
@@ -261,10 +260,17 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
         Observation variable name for labeling.
     grp_dict : dict
         Plot configuration options including aesthetics and normalization parameters.
-    vmin : float or None
+    vmin : float
+        Min value to use on y-axis.
+    vmax : float
+        Max value to use on y-axis.
+    cmin : float or None
         Minimum value for color normalization, if applicable.
-    vmax : float or None
+    cmax : float or None
         Maximum value for color normalization, if applicable.
+    plot_dict : dictionary
+        Dictionary containing information about plotting for each pair
+        (e.g., color, linestyle, markerstyle).
     domain_type : str
         Type of domain being plotted (e.g., 'region', 'global').
     domain_name : str
@@ -273,21 +279,27 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
         Configuration dictionary for observation labels.
     model_scatter_data : numpy.ndarray or None
         Data for model scatter plot overlay.
+    debug : bool
+        Whether to plot interactively (True) or not (False). Flag for submitting jobs to supercomputer turn off interactive mode.
 
     Returns
     -------
     None
         Saves the generated plot to a file.
     """
+
+    if vmin is not None and vmax is not None:
+        plot_dict['ylim'] = [vmin, vmax]
+    
     if model_data_2d.size == 0 or obs_data_2d.size == 0:
         print("Warning: Model or observation data is empty. Skipping plot.")
         return
 
-    # Determine vmin and vmax from YAML config if provided
+    # Determine cmin and cmax (colorbar min max) from YAML config if provided
     if 'variable_limits' in grp_dict:
         var_limits = grp_dict['variable_limits'].get(obs_var, {})
-        vmin = var_limits.get('vmin', vmin)
-        vmax = var_limits.get('vmax', vmax)
+        cmin = var_limits.get('cmin', cmin)
+        cmax = var_limits.get('cmax', cmax)
 
     # Debugging: print shapes and ndims
     print(f"time shape: {time.shape}, ndims: {time.ndim}")
@@ -295,16 +307,16 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
     print(f"model_data_2d shape: {model_data_2d.shape}, ndims: {model_data_2d.ndim}")
     print(f"obs_data_2d shape: {obs_data_2d.shape}, ndims: {obs_data_2d.ndim}")
 
-    # Determine vmin and vmax dynamically if not provided
-    if vmin is None or vmax is None:
-        vmin = min(np.nanmin(model_data_2d), np.nanmin(obs_data_2d))
-        vmax = max(np.nanmax(model_data_2d), np.nanmax(obs_data_2d))
+    # Determine cmin and cmax dynamically if not provided
+    if cmin is None or cmax is None:
+        cmin = min(np.nanmin(model_data_2d), np.nanmin(obs_data_2d))
+        cmax = max(np.nanmax(model_data_2d), np.nanmax(obs_data_2d))
     
     # Set colorbar min to zero if the minimum value is less than zero
-    if vmin < 0:
-        vmin = 0
+    if cmin < 0:
+        cmin = 0
 
-    print(f"vmin and vmax set dynamically based on model and observation data to {vmin} and {vmax}")
+    print(f"cmin and cmax set dynamically based on model and observation data to {cmin} and {cmax}")
 
     fig, axs = plt.subplots(nrows=2, figsize=grp_dict.get('fig_kwargs', {}).get('figsize', (20, 8)), sharex=True, gridspec_kw={'height_ratios': [1, 1]})
 
@@ -315,18 +327,18 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
     else:
         cmap = plt.get_cmap(grp_dict.get('color_map', 'viridis'))
 
-    norm = Normalize(vmin=vmin, vmax=vmax)
+    norm = Normalize(vmin=cmin, vmax=cmax)
 
     time_dates = mdates.date2num(time)
     time_mesh, altitude_mesh = np.meshgrid(time_dates, altitude, indexing='ij')
 
-    levels = np.linspace(vmin, vmax, 200)
+    levels = np.linspace(cmin, cmax, 200)
     contourf_plot = axs[0].contourf(time_mesh, altitude_mesh, model_data_2d, levels=levels, cmap=cmap, norm=norm, extend='both')
 
     # Adjust the layout to place the colorbar at the bottom
     cbar = fig.colorbar(contourf_plot, ax=axs.ravel().tolist(), orientation='horizontal', pad=0.1, aspect=50)
-    cbar.set_label(obs_label_config.get(obs_var, {}).get('ylabel_plot', 'Variable (units)'), fontsize=12, fontweight='bold')
-    cbar.ax.tick_params(labelsize=10)
+    cbar.set_label(obs_label_config.get(obs_var, {}).get('ylabel_plot', 'Variable (units)'), fontsize=18, fontweight='bold')
+    cbar.ax.tick_params(labelsize=14)
 
     axs[0].xaxis.set_major_formatter(mdates.DateFormatter('%m-%d %H'))
     axs[0].xaxis.set_major_locator(mdates.AutoDateLocator())
@@ -335,27 +347,31 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
     #axs[0].set_title(f"Model vs Observation Curtain Plot: {model_var} vs {obs_var}", fontsize=12, fontweight='bold')
     # Set the main title and subplot titles
     fig.suptitle(f"Model vs Observation Curtain Plot: {model_var} vs {obs_var}", fontsize=16, fontweight='bold')
-    axs[0].set_title("Model Curtain with Model Scatter Overlay", fontsize=12, fontweight='bold')
+    axs[0].set_title("Model Curtain with Model Scatter Overlay", fontsize=18, fontweight='bold')
     
-    axs[0].set_ylabel('Pressure (Pa)', fontsize=12, fontweight='bold')
+    axs[0].set_ylabel('Pressure (Pa)', fontsize=18, fontweight='bold')
     axs[0].invert_yaxis()  # Invert y-axis to have max pressure at the bottom
-    axs[0].tick_params(axis='both', labelsize=10)
+    axs[0].tick_params(axis='both', labelsize=14)
 
     # Separate subplot for the observation scatter plot
-    axs[1].set_title("Observation Scatter", fontsize=12, fontweight='bold')
+    axs[1].set_title("Observation Scatter", fontsize=18, fontweight='bold')
     scatter = axs[1].scatter(time_dates, obs_pressure, c=obs_data_2d, cmap=cmap, norm=norm, edgecolor=(0, 0, 0, 0), linewidth=0.3, alpha=0.5)
     axs[1].xaxis.set_major_formatter(mdates.DateFormatter('%m-%d %H'))
     axs[1].xaxis.set_major_locator(mdates.AutoDateLocator())
     fig.autofmt_xdate(rotation=45, ha='right')
 
-    axs[1].set_ylabel('Pressure (Pa)', fontsize=12, fontweight='bold')
+    axs[1].set_ylabel('Pressure (Pa)', fontsize=18, fontweight='bold')
     axs[1].invert_yaxis()  # Invert y-axis to have max pressure at the bottom
-    axs[1].tick_params(axis='both', labelsize=10)
-    axs[1].set_xlabel('Time', fontsize=12, fontweight='bold')
+    axs[1].tick_params(axis='both', labelsize=14)
+    axs[1].set_xlabel('Time', fontsize=18, fontweight='bold')
+
+    
 
     # Set the same y-tick labels for both subplots and exclude the minimum value
     y_ticks = axs[0].get_yticks()
+    y_lim = axs[0].get_ylim()
     axs[1].set_yticks(y_ticks)
+    axs[1].set_ylim(y_lim)
     axs[1].set_yticklabels([str(int(tick)) if tick != y_ticks.min() else '' for tick in y_ticks])
 
     # Plot model scatter data on top of the model curtain plot
@@ -367,11 +383,14 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
 
     plt.show()
 
+    # Only close the plot if not in debug mode
+    if not debug:
+        plt.close()
 
     #Diagnostic Histogram
     #plt.figure(figsize=(10, 4))
-    #plt.hist(model_data_2d.values.flatten(), bins=100, alpha=0.15, color='blue', label='Model', range=(vmin, vmax))
-    #plt.hist(obs_data_2d.flatten(), bins=100, alpha=0.5, color='green', label='Observation', range=(vmin, vmax))
+    #plt.hist(model_data_2d.values.flatten(), bins=100, alpha=0.15, color='blue', label='Model', range=(cmin, cmax))
+    #plt.hist(obs_data_2d.flatten(), bins=100, alpha=0.5, color='green', label='Observation', range=(cmin, cmax))
     #plt.legend()
     #plt.title('Distribution of Model and Observation Data')
     #plt.xlabel('Value')

From 21b893789b1d8cd89a96725bfc6f5de556bbd499 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Wed, 24 Jul 2024 11:35:25 -0600
Subject: [PATCH 08/20] Saving plots for scatter and curtain plots for
 different model-obs pair separately as png now

---
 melodies_monet/driver.py              | 20 +++++++++++---------
 melodies_monet/plots/aircraftplots.py | 18 +++++++++++++-----
 2 files changed, 24 insertions(+), 14 deletions(-)

diff --git a/melodies_monet/driver.py b/melodies_monet/driver.py
index 1a867b32..5fdf7e33 100644
--- a/melodies_monet/driver.py
+++ b/melodies_monet/driver.py
@@ -1742,8 +1742,12 @@ def plotting(self):
                                 ##print(f"Length of obs_pressure: {len(obs_pressure)}")
                                 ##print(f"Length of obs_data_2d: {len(obs_data_2d)}")
                                 ##print(f"model_scatter_data shape: {model_scatter_data.shape}")
+
+                                outname_pair = f"{outname}_{obs_label}_vs_{model_label}.png"
+
+                                print(f"Saving curtain plot to {outname_pair}...")
                         
-                                airplots.make_curtain_plot(
+                                ax = airplots.make_curtain_plot(
                                     time=pd.to_datetime(time),
                                     altitude=target_pressures,  # Use target_pressures for interpolation
                                     model_data_2d=model_data_2d,  # Already reshaped to match the expected shape
@@ -1756,6 +1760,7 @@ def plotting(self):
                                     vmin=vmin,
                                     vmax=vmax,
                                     plot_dict=plot_dict,
+                                    outname=outname_pair,
                                     grp_dict=curtain_config,
                                     domain_type=domain_type,
                                     domain_name=domain_name,
@@ -1764,10 +1769,7 @@ def plotting(self):
                                     debug=self.debug  # Pass debug flag
                                 )
                             
-                                # Save the curtain plot for the current pair immediately
-                                outname_pair = f"{outname}_{obs_label}_vs_{model_label}.png"
-                                print(f"Saving curtain plot to {outname_pair}...")
-                                plt.savefig(outname_pair)
+                                
                             except Exception as e:
                                 print(f"Error generating curtain plot for {modvar} vs {obsvar}: {e}")
                             finally:
@@ -1966,7 +1968,9 @@ def plotting(self):
                                 del kwargs['shade_lowest']
 
 
+                            outname_pair = f"{outname}_{obs_label}_vs_{model_label}.png"
 
+                            print(f"Saving scatter density plot to {outname_pair}...")
                             
                             # Create the scatter density plot
                             print(f"Processing scatter density plot for model '{model_label}' and observation '{obs_label}'...")
@@ -1983,13 +1987,11 @@ def plotting(self):
                                 vmax_x=vmax_x,
                                 vmin_y=vmin_y,
                                 vmax_y=vmax_y,
+                                outname=outname_pair,
                                 **kwargs                            
                             )
 
-                            # Save the scatter density plot for the current pair immediately
-                            outname_pair = f"{outname}_{obs_label}_vs_{model_label}.png"
-                            print(f"Saving scatter density plot to {outname_pair}...")  # Debugging print statement
-                            plt.savefig(outname_pair)
+                            
                             plt.close()  # Close the current figure
                             
                         elif plot_type.lower() == 'boxplot':
diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index 2fa49da8..68c7b1a2 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -232,7 +232,7 @@ def add_yax2_altitude(ax, pairdf, altitude_yax2, text_kwargs, vmin_y2, vmax_y2):
 
 
 ###NEW curtain plot qzr++  (NEW CURTAIN model plot with model overlay, shared x-axis 
-def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d, model_var, obs_var, grp_dict, vmin=None, vmax=None, cmin=None, cmax=None, plot_dict=None, domain_type=None, domain_name=None, obs_label_config=None, model_scatter_data=None, debug=False):
+def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d, model_var, obs_var, grp_dict, vmin=None, vmax=None, cmin=None, cmax=None, plot_dict=None, outname='plot', domain_type=None, domain_name=None, obs_label_config=None, model_scatter_data=None, debug=False):
     """
     Generates a curtain plot comparing model data with obs across altitude (Pressure, right now) over time,
     with the ability to customize the appearance through a configuration dictionary.
@@ -271,6 +271,8 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
     plot_dict : dictionary
         Dictionary containing information about plotting for each pair
         (e.g., color, linestyle, markerstyle).
+    outname : str
+        File location and name of plot.
     domain_type : str
         Type of domain being plotted (e.g., 'region', 'global').
     domain_name : str
@@ -380,7 +382,10 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
 
         # Print diagnostics for model scatter data
         print(f"model_scatter_data shape: {model_scatter_data.shape}, ndims: {model_scatter_data.ndim}")
-
+        
+    # Save the curtain plot for the current pair immediately
+    print(f"Saving curtain plot to {outname}...")
+    savefig(f"{outname}", loc=4, logo_height=100, dpi=300)
     plt.show()
 
     # Only close the plot if not in debug mode
@@ -398,8 +403,6 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
     #plt.show()
 
 
-
-
 ####NEW vertprofile has option for both shading (for interquartile range) or box (interquartile range)-whisker (10th-90th percentile bounds) (qzr++)
 def make_vertprofile(df, column=None, label=None, ax=None, bins=None, altitude_variable=None, ylabel=None,
                      vmin=None, vmax=None, 
@@ -716,7 +719,7 @@ def make_vertprofile(df, column=None, label=None, ax=None, bins=None, altitude_v
 
 
 ##NEW Scatter Density Plot for model obs pairs (matplotlib scatter plot if fill=False or seaborn kde sactter density plot if fill= True)
-def make_scatter_density_plot(df, mod_var=None, obs_var=None, ax=None, color_map='viridis', xlabel=None, ylabel=None, title=None, fill=False, vmin_x=None, vmax_x=None, vmin_y=None, vmax_y=None, **kwargs):
+def make_scatter_density_plot(df, mod_var=None, obs_var=None, ax=None, color_map='viridis', xlabel=None, ylabel=None, title=None, fill=False, vmin_x=None, vmax_x=None, vmin_y=None, vmax_y=None, outname='plot', **kwargs):
     
     """  
     Creates a scatter density plot for the specified column (variable) in the paired DataFrame (df).
@@ -741,6 +744,8 @@ def make_scatter_density_plot(df, mod_var=None, obs_var=None, ax=None, color_map
         Title for the plot (optional)
     fill: bool
         Fill set to True for seaborn kde plot
+    outname : str
+        File location and name of plot.
     **kwargs: dict 
         Additional keyword arguments for customization
 
@@ -844,6 +849,9 @@ def make_scatter_density_plot(df, mod_var=None, obs_var=None, ax=None, color_map
         cbar = plt.colorbar(mappable, ax=ax)
     cbar.set_label(colorbar_label)
 
+    # Save the scatter density plot for the current pair immediately
+    print(f"Saving scatter density plot to {outname}...")
+    savefig(f"{outname}", loc=4, logo_height=100, dpi=300)
     plt.show()
 
     return ax

From 7f1814c45971a118a28f99483c5f499e266cdb89 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Mon, 29 Jul 2024 10:22:52 -0600
Subject: [PATCH 09/20] Comment color_map in yaml out because it's not used
 when color_map_custom is used.

---
 ...raft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
index 47906b47..6f50fb01 100644
--- a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
+++ b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
@@ -263,7 +263,7 @@ plots:
         "#343434", "#606060", "#8c8c8c", "#b8b8b8", "#e4e4e4"] #,
         #"#ffffff"] #['blue','royalblue', 'cyan', 'yellow', 'orange', 'red'] ##['#440154', '#3b528b', '#21908d', '#5dc963', '#fde725'] #['purple', 'blue', 'green', 'yellow'] # Example gradient
     color_levels: 30 #31 #15 #10  # Define the number of distinct colors in the color bar
-    color_map: 'Spectral_r' #'jet'  # This tells the plotting function to use the built-in 'rainbow' colormap   
+    #color_map: 'Spectral_r' #'jet'  # This can be set instead of color_map_custom, color_levels and color_levels to set a colormap defined in matplotlob.  
     variable_limits:
       NO2_CL_RYERSON:
         cmin: 0

From 160c1bb8b4b146c3f71397834c4ee2286389308f Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Tue, 30 Jul 2024 14:07:47 -0600
Subject: [PATCH 10/20] Some redundancies related to interpolation removed in
 driver, custom_color_map True statement in aircratplots more explicit,
 colobar minima not forced to be zero now

---
 melodies_monet/driver.py              | 5 -----
 melodies_monet/plots/aircraftplots.py | 8 ++++----
 2 files changed, 4 insertions(+), 9 deletions(-)

diff --git a/melodies_monet/driver.py b/melodies_monet/driver.py
index 5fdf7e33..20b187db 100644
--- a/melodies_monet/driver.py
+++ b/melodies_monet/driver.py
@@ -1683,11 +1683,6 @@ def plotting(self):
                         
                             time_dates = mdates.date2num(pd.to_datetime(time))
                         
-                            if len(obs_pressure) != len(time_dates):
-                                obs_pressure = np.interp(time_dates, mdates.date2num(pd.to_datetime(new_ds_obs['time_obs'].values)), obs_pressure)
-                                obs_data_2d = np.interp(time_dates, mdates.date2num(pd.to_datetime(new_ds_obs['time_obs'].values)), obs_data_2d)
-                                ##print(f"Length of obs_pressure after interpolation: {len(obs_pressure)}")
-                                ##print(f"Length of obs_data_2d after interpolation: {len(obs_data_2d)}")
                         
                             # Fetch the observation configuration for colorbar labels
                             obs_label_config = self.control_dict['obs'][p.obs]['variables']
diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index 68c7b1a2..7b24e42f 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -314,19 +314,19 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
         cmin = min(np.nanmin(model_data_2d), np.nanmin(obs_data_2d))
         cmax = max(np.nanmax(model_data_2d), np.nanmax(obs_data_2d))
     
-    # Set colorbar min to zero if the minimum value is less than zero
-    if cmin < 0:
-        cmin = 0
+    
 
     print(f"cmin and cmax set dynamically based on model and observation data to {cmin} and {cmax}")
 
     fig, axs = plt.subplots(nrows=2, figsize=grp_dict.get('fig_kwargs', {}).get('figsize', (20, 8)), sharex=True, gridspec_kw={'height_ratios': [1, 1]})
 
     # Handling custom color maps
-    if grp_dict.get('color_map_custom', False) and 'colors' in grp_dict:
+    if 'color_map_custom' in grp_dict and grp_dict['color_map_custom'] is True and 'colors' in grp_dict:
+        # Custom color map logic
         colors = grp_dict['colors']
         cmap = LinearSegmentedColormap.from_list("custom_cmap", colors, N=grp_dict.get('color_levels', 100))
     else:
+    # Default color map logic
         cmap = plt.get_cmap(grp_dict.get('color_map', 'viridis'))
 
     norm = Normalize(vmin=cmin, vmax=cmax)

From 9681061e46c565bd78142a6c6a00e960a1015584 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Tue, 30 Jul 2024 16:38:40 -0600
Subject: [PATCH 11/20] colobar min and max now taken from vmin_plot ad
 vmax_plot from obs key list in YAML for different variables plotted

---
 ...velop_aircraft_testmultiplemodelgroups.yaml | 18 ++++--------------
 melodies_monet/driver.py                       | 10 +++++-----
 melodies_monet/plots/aircraftplots.py          |  6 ------
 3 files changed, 9 insertions(+), 25 deletions(-)

diff --git a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
index 6f50fb01..a50ae7fe 100644
--- a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
+++ b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
@@ -86,8 +86,8 @@ obs:
         LLOD_value: -888888 # Opt Set this value to LLOD_setvalue
         LLOD_setvalue: 0.0 # Opt Set LLOD_value=LLOD_setvalue, applied AFTER unit_scale and obs_unit
         ylabel_plot: 'Ozone (ppbv)'
-        vmin_plot: 15.0 #Opt Min for y-axis during plotting. To apply to a plot, change restrict_yaxis = True.
-        vmax_plot: 55.0 #Opt Max for y-axis during plotting. To apply to a plot, change restrict_yaxis = True.
+        vmin_plot: 0 #15.0 #Opt Min for y-axis during plotting. To apply to a plot, change restrict_yaxis = True.
+        vmax_plot: 100 #55.0 #Opt Max for y-axis during plotting. To apply to a plot, change restrict_yaxis = True.
         vdiff_plot: 20.0 #Opt +/- range to use in bias plots. To apply to a plot, change restrict_yaxis = True.
         # nlevels_plot: 21 #Opt number of levels used in colorbar for contourf plot.
       'NO_CL_RYERSON':
@@ -96,7 +96,7 @@ obs:
         LLOD_setvalue: 0.0 # Opt Set LLOD_value=LLOD_setvalue, applied AFTER unit_scale and obs_unit
         ylabel_plot: 'NO (ppbv)' #Optional to set ylabel so can include units and/or instr etc.
         vmin_plot: 0.0 #Opt Min for y-axis during plotting. To apply to a plot, change restrict_yaxis = True.
-        vmax_plot: 20.0 #Opt Max for y-axis during plotting. To apply to a plot, change restrict_yaxis = True.
+        vmax_plot: 1 #20.0 #Opt Max for y-axis during plotting. To apply to a plot, change restrict_yaxis = True.
         vdiff_plot: 15.0 #Opt +/- range to use in bias plots. To apply to a plot, change restrict_yaxis = True.
         nlevels_plot: 21 #Opt number of levels used in colorbar for contourf plot.
       'NO2_CL_RYERSON':
@@ -105,7 +105,7 @@ obs:
         LLOD_setvalue: 0.0 # Opt Set LLOD_value=LLOD_setvalue, applied AFTER unit_scale and obs_unit
         ylabel_plot: 'NO2 (ppbv)' #Optional to set ylabel so can include units and/or instr etc.
         vmin_plot: 0.0 #Opt Min for y-axis during plotting. To apply to a plot, change restrict_yaxis = True.
-        vmax_plot: 20.0 #Opt Max for y-axis during plotting. To apply to a plot, change restrict_yaxis = True.
+        vmax_plot: 1 #20.0 #Opt Max for y-axis during plotting. To apply to a plot, change restrict_yaxis = True.
         vdiff_plot: 15.0 #Opt +/- range to use in bias plots. To apply to a plot, change restrict_yaxis = True.
         nlevels_plot: 21 #Opt number of levels used in colorbar for contourf plot.
       'Latitude_YANG':
@@ -264,16 +264,6 @@ plots:
         #"#ffffff"] #['blue','royalblue', 'cyan', 'yellow', 'orange', 'red'] ##['#440154', '#3b528b', '#21908d', '#5dc963', '#fde725'] #['purple', 'blue', 'green', 'yellow'] # Example gradient
     color_levels: 30 #31 #15 #10  # Define the number of distinct colors in the color bar
     #color_map: 'Spectral_r' #'jet'  # This can be set instead of color_map_custom, color_levels and color_levels to set a colormap defined in matplotlob.  
-    variable_limits:
-      NO2_CL_RYERSON:
-        cmin: 0
-        cmax: 1 #10
-      NO_CL_RYERSON:
-        cmin: 0
-        cmax: 1 #15
-      O3_CL_RYERSON:
-        cmin: 0
-        cmax: 100
     vmin: 120000 #hPa  #0  #Optional (y-axis limits)
     vmax: 5000 # Optional #need to be edited as per min and max altitude (i.e., vmin and vmax)
     
diff --git a/melodies_monet/driver.py b/melodies_monet/driver.py
index 20b187db..8392b704 100644
--- a/melodies_monet/driver.py
+++ b/melodies_monet/driver.py
@@ -1592,11 +1592,6 @@ def plotting(self):
                                 vmax = None
                             curtain_config = grp_dict
 
-                            # Determine cmin and cmax (colorbar min max) from YAML config if provided
-                            if 'variable_limits' in curtain_config:
-                                var_limits = curtain_config['variable_limits'].get(obsvar, {})
-                                cmin = var_limits.get('cmin')
-                                cmax = var_limits.get('cmax')
 
                             # Inside your loop for processing each pair
                             obs_label = p.obs
@@ -1608,6 +1603,11 @@ def plotting(self):
                             model_obj = mod.obj  # This assumes `mod.obj` contains the necessary model data
                             model_data = model_obj[modvar].values
                             obs_data = obs.obj[obsvar].values
+
+                            # Determine cmin and cmax (colorbar min max) from observation config if provided
+                            cmin = obs_plot_dict.get('vmin_plot', None)
+                            cmax = obs_plot_dict.get('vmax_plot', None)
+
                         
                             # Debugging: print dimensions
                             ##print(f"Dimensions of model data: {model_data.shape}, ndims: {model_data.ndim}")
diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index 7b24e42f..ed252665 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -297,12 +297,6 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
         print("Warning: Model or observation data is empty. Skipping plot.")
         return
 
-    # Determine cmin and cmax (colorbar min max) from YAML config if provided
-    if 'variable_limits' in grp_dict:
-        var_limits = grp_dict['variable_limits'].get(obs_var, {})
-        cmin = var_limits.get('cmin', cmin)
-        cmax = var_limits.get('cmax', cmax)
-
     # Debugging: print shapes and ndims
     print(f"time shape: {time.shape}, ndims: {time.ndim}")
     print(f"altitude shape: {altitude.shape}, ndims: {altitude.ndim}")

From 59cbbfeb710e76bf08c9a42fdc42fe97e8548e1c Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Wed, 31 Jul 2024 21:15:50 -0600
Subject: [PATCH 12/20] Model and Obs Scatter overlay now plotted using pairdf
 as in other plot types, Also in aircraftplots make_curtain_plot() now
 fontsize, fontweight and labelsize all linked via YAML

---
 ...elop_aircraft_testmultiplemodelgroups.yaml |   4 +-
 melodies_monet/driver.py                      | 131 ++++++------------
 melodies_monet/plots/aircraftplots.py         |  56 ++++----
 3 files changed, 71 insertions(+), 120 deletions(-)

diff --git a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
index a50ae7fe..df7b44ee 100644
--- a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
+++ b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
@@ -246,6 +246,7 @@ plots:
     text_kwargs:
       fontsize: 25 #18
       fontweight: 'bold'
+      labelsize: 16
     domain_type: ['all']
     domain_name: ['CONUS']
     data: ['firexaq_wrfchem_v4.2','firexaq_wrfchem_v4.2_test']
@@ -264,8 +265,9 @@ plots:
         #"#ffffff"] #['blue','royalblue', 'cyan', 'yellow', 'orange', 'red'] ##['#440154', '#3b528b', '#21908d', '#5dc963', '#fde725'] #['purple', 'blue', 'green', 'yellow'] # Example gradient
     color_levels: 30 #31 #15 #10  # Define the number of distinct colors in the color bar
     #color_map: 'Spectral_r' #'jet'  # This can be set instead of color_map_custom, color_levels and color_levels to set a colormap defined in matplotlob.  
-    vmin: 120000 #hPa  #0  #Optional (y-axis limits)
+    vmin: 120000 #Pressure in Pa  #0  #Optional (y-axis limits)
     vmax: 5000 # Optional #need to be edited as per min and max altitude (i.e., vmin and vmax)
+    interval : 500 #Pressure in Pa; target pressure interpolation steps for model curtain 2d plot
     
   plot_grp6:
     type: 'taylor' # plot type
diff --git a/melodies_monet/driver.py b/melodies_monet/driver.py
index 8392b704..47fffb68 100644
--- a/melodies_monet/driver.py
+++ b/melodies_monet/driver.py
@@ -1057,9 +1057,23 @@ def pair_data(self, time_interval=None):
                     ds_model = m.util.combinetool.combine_da_to_da(model_obj,new_ds_obs,merge=False)
                     #Interpolate based on time in the observations
                     ds_model = ds_model.interp(time=ds_model.time_obs.squeeze())
+
+                    # Debugging: Print the variables in ds_model to verify 'pressure_model' is included  ##qzr++
+                    #print("Variables in ds_model after combine_da_to_da and interp:", ds_model.variables)
+                    
+                    # Ensure 'pressure_model' is included in ds_model (checked it exists)
+                    #if 'pressure_model' not in ds_model:
+                     #   raise KeyError("'pressure_model' is missing in the model dataset")   #qzr++
                     
                     paired_data = vert_interp(ds_model,obs.obj,keys+mod_vars)
                     print('After pairing: ', paired_data)
+
+                    # Ensure 'pressure_model' is included in the DataFrame (pairdf) #qzr++
+                    #if 'pressure_model' not in paired_data.columns:
+                       # raise KeyError("'pressure_model' is missing in the paired_data")   #qzr++
+
+                      
+                    
                     # this outputs as a pandas dataframe.  Convert this to xarray obj
                     p = pair()
                     p.type = 'aircraft'
@@ -1590,35 +1604,29 @@ def plotting(self):
                             else:
                                 vmin = None
                                 vmax = None
-                            curtain_config = grp_dict
-
-
+                                
+                            curtain_config = grp_dict # Curtain plot grp YAML dict
                             # Inside your loop for processing each pair
                             obs_label = p.obs
                             model_label = p.model
                         
-                            # Ensure we use the correct observation and model objects from pairing
+                            
+                            #Ensure we use the correct observation and model objects from pairing
                             obs = self.obs[p.obs]
-                            mod = self.models[p.model]  # Get the model associated with this pair
-                            model_obj = mod.obj  # This assumes `mod.obj` contains the necessary model data
-                            model_data = model_obj[modvar].values
-                            obs_data = obs.obj[obsvar].values
+                            mod = self.models[p.model]
+                            model_obj = mod.obj
+                        
+                            # Fetch the observation configuration for colorbar labels
+                            obs_label_config = self.control_dict['obs'][obs_label]['variables']
+                        
+                            # Fetch the model and observation data from pairdf
+                            pairdf = pairdf_all.reset_index()
 
                             # Determine cmin and cmax (colorbar min max) from observation config if provided
                             cmin = obs_plot_dict.get('vmin_plot', None)
                             cmax = obs_plot_dict.get('vmax_plot', None)
-
-                        
-                            # Debugging: print dimensions
-                            ##print(f"Dimensions of model data: {model_data.shape}, ndims: {model_data.ndim}")
-                            ##print(f"Dimensions of observation data: {obs_data.shape}, ndims: {obs_data.ndim}")
-                        
-                            if not isinstance(obs.obj, pd.DataFrame):
-                                obs.obj = obs.obj.to_dataframe()
-                        
-                            # Drop any variables where coords NaN
-                            obs.obj = obs.obj.reset_index().dropna(subset=['pressure_obs', 'latitude', 'longitude']).set_index('time')
                         
+                            #### For model_data_2d for curtain/contourfill plot #####                       
                             # Convert to get something useful for MONET
                             new_ds_obs = obs.obj.rename_axis('time_obs').reset_index().monet._df_to_da().set_coords(['time_obs', 'pressure_obs'])
                         
@@ -1626,19 +1634,17 @@ def plotting(self):
                             ds_model = m.util.combinetool.combine_da_to_da(model_obj, new_ds_obs, merge=False)
                         
                             # Interpolate based on time in the observations
-                            ds_model = ds_model.interp(time=ds_model.time_obs.squeeze())
-                        
-                            # Fetch the observation configuration for colorbar labels and plot limits
-                            obs_plot_dict = self.control_dict['obs'][p.obs]['variables'][obsvar]
+                            ds_model = ds_model.interp(time=ds_model.time_obs.squeeze())                 
                         
-                            # Print ds_model and pressure_model values
+                            # Print ds_model and pressure_model values #Debugging
                             ##print(f"ds_model: {ds_model}")
                             ##print(f"pressure_model values: {ds_model['pressure_model'].values}")
                         
                             # Define target pressures for interpolation based on the range of pressure_model
                             min_pressure = ds_model['pressure_model'].min().compute()
                             max_pressure = ds_model['pressure_model'].max().compute()
-                            interval = 100  # Interval in Pa
+                            # Fetch the interval from curtain_config
+                            interval = curtain_config.get('interval', 100)  # Default to 100 (in Pa) if not provided 
                         
                             print(f"Pressure MIN:{min_pressure}, max: {max_pressure}, interval: {interval}")
                         
@@ -1650,6 +1656,7 @@ def plotting(self):
                             # Check for NaN values before interpolation
                             ##print(f"NaNs in model_data before interpolation: {np.isnan(ds_model[modvar]).sum().compute()}")
                             ##print(f"NaNs in pressure_model before interpolation: {np.isnan(ds_model['pressure_model']).sum().compute()}")
+
                         
                             # Resample model data to target pressures using stratify
                             da_wrf_const = resample_stratify(ds_model[modvar], target_pressures, ds_model['pressure_model'], axis=1, interpolation='linear', extrapolation='nan')
@@ -1666,70 +1673,20 @@ def plotting(self):
                             # Debugging: print merged dataset
                             ##print(ds_wrf_const)
                         
-                            # Ensure model_data_2d is properly reshaped for the contourf plot
+                            # Ensure model_data_2d is properly reshaped for the contourfill plot
                             model_data_2d = ds_wrf_const[modvar].squeeze()
                         
                             # Debugging: print reshaped model data shape
                             ##print(f"Reshaped model data shape: {model_data_2d.shape}")
-                        
-                            # Align the observation data with the time array length for scatter plot
-                            time = ds_wrf_const['time'].values
-                            obs_pressure = new_ds_obs['pressure_obs'].values.flatten()
-                            obs_data_2d = obs_data.flatten()
-                        
-                            ##print(f"time shape: {time.shape}, ndims: {time.ndim}")
-                            ##print(f"obs_pressure shape: {obs_pressure.shape}, ndims: {obs_pressure.ndim}")
-                            ##print(f"obs_data_2d shape: {obs_data_2d.shape}, ndims: {obs_data_2d.ndim}")
-                        
-                            time_dates = mdates.date2num(pd.to_datetime(time))
-                        
-                        
-                            # Fetch the observation configuration for colorbar labels
-                            obs_label_config = self.control_dict['obs'][p.obs]['variables']
-                            
-                            # Ensure the observation data is a DataFrame and reset index
-                            df_obs = new_ds_obs.to_dataframe().reset_index()
-                            df_obs['time'] = df_obs['time_obs']
-                        
-                            # Create var_name_list dynamically based on modvar and pressure_model
-                            var_name_list = [modvar, 'pressure_model'] if modvar != 'pressure_model' else [modvar]
-                            available_vars = [var for var in var_name_list if var in ds_model]
-                            
-                            # Print diagnostics for available variables
-                            ##print(f"Available variables for interpolation: {available_vars}")
-                            
-                            # Convert ds_model to DataFrame and reset index
-                            df_model = ds_model[available_vars].to_dataframe().reset_index()
-                            
-                            # Print diagnostics for the observation DataFrame
-                            ##print(f"Type of df_obs: {type(df_obs)}")
-                            ##print(f"Head of df_obs:\n{df_obs.head()}")
                             
-                            # Print diagnostics for the model DataFrame
-                            ##print(f"Type of df_model: {type(df_model)}")
-                            ##print(f"Head of df_model:\n{df_model.head()}")
-                            
-                            # Ensure both DataFrames have the necessary columns for merging
-                            common_columns = set(df_obs.columns) & set(df_model.columns)
-                            ##print(f"Common columns between observation and model DataFrames: {common_columns}")
-                            if not {'latitude', 'longitude', 'pressure_obs', 'time'}.issubset(common_columns):
-                                raise ValueError("Missing necessary columns for merging. Columns needed: 'latitude', 'longitude', 'pressure_obs', 'time'")
-                            
-                            # Perform vertical interpolation only for available variables
-                            df_wrf = vert_interp(ds_model, df_obs, available_vars)
-                            
-                            # Print diagnostics for the resulting DataFrame
-                            ##print(f"Type of df_wrf: {type(df_wrf)}")
-                            ##print(f"Head of df_wrf:\n{df_wrf.head()}")
-                            
-                            # Ensure the modvar exists in the resulting DataFrame
-                            if modvar not in df_wrf.columns:
-                                raise ValueError(f"{modvar} not found in the resulting DataFrame columns")
-                            
-                            model_scatter_data = df_wrf[modvar].values
-                            
-                            # Print diagnostics for model scatter data
-                            ##print(f"model_scatter_data shape: {model_scatter_data.shape}, ndims: {model_scatter_data.ndim}")
+                            #### model_data_2d for curtain plot ready ####
+
+
+                            # Fetch model pressure and other model and observation data from "pairdf" (for scatter plot overlay)
+                            time = pairdf['time']
+                            obs_pressure = pairdf['pressure_obs']
+                            obs_data = pairdf[obsvar]
+                            model_data = pairdf[modvar]        
                         
                             # Generate the curtain plot using airplots.make_curtain_plot
                             try:
@@ -1746,9 +1703,9 @@ def plotting(self):
                                     time=pd.to_datetime(time),
                                     altitude=target_pressures,  # Use target_pressures for interpolation
                                     model_data_2d=model_data_2d,  # Already reshaped to match the expected shape
+                                    pairdf=pairdf,  #use pairdf for scatter overlay (model and obs)
                                     obs_pressure=obs_pressure,  # Pressure_obs for obs scatter plot
-                                    obs_data_2d=obs_data_2d,  # Use original observation data for scatter plot
-                                    model_var=modvar,
+                                    mod_var=modvar,
                                     obs_var=obsvar,
                                     cmin=cmin,
                                     cmax=cmax,
@@ -1760,7 +1717,7 @@ def plotting(self):
                                     domain_type=domain_type,
                                     domain_name=domain_name,
                                     obs_label_config=obs_label_config,
-                                    model_scatter_data=model_scatter_data,  # Pass model scatter data for plotting
+                                    text_dict=text_dict,
                                     debug=self.debug  # Pass debug flag
                                 )
                             
diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index ed252665..05a3c3d3 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -232,7 +232,7 @@ def add_yax2_altitude(ax, pairdf, altitude_yax2, text_kwargs, vmin_y2, vmax_y2):
 
 
 ###NEW curtain plot qzr++  (NEW CURTAIN model plot with model overlay, shared x-axis 
-def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d, model_var, obs_var, grp_dict, vmin=None, vmax=None, cmin=None, cmax=None, plot_dict=None, outname='plot', domain_type=None, domain_name=None, obs_label_config=None, model_scatter_data=None, debug=False):
+def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, pairdf, mod_var, obs_var, grp_dict, vmin=None, vmax=None, cmin=None, cmax=None, plot_dict=None, outname='plot', domain_type=None, domain_name=None, obs_label_config=None, text_dict=None, debug=False):
     """
     Generates a curtain plot comparing model data with obs across altitude (Pressure, right now) over time,
     with the ability to customize the appearance through a configuration dictionary.
@@ -252,9 +252,9 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
         2D array of cleaned and interpolated model data.
     obs_pressure : numpy.ndarray
         Array of pressure points corresponding to observations.
-    obs_data_2d : numpy.ndarray
-        2D array of observational data.
-    model_var : str
+    pairdf : pandas.DataFrame
+        Paired DataFrame containing model and observation data for scatter plots.
+    mod_var : str
         Model variable name for labeling.
     obs_var : str
         Observation variable name for labeling.
@@ -279,8 +279,8 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
         Name of the domain being plotted.
     obs_label_config : dict
         Configuration dictionary for observation labels.
-    model_scatter_data : numpy.ndarray or None
-        Data for model scatter plot overlay.
+    text_dict : dict
+        Dictionary containing text properties (fontsize, fontweight, etc.).
     debug : bool
         Whether to plot interactively (True) or not (False). Flag for submitting jobs to supercomputer turn off interactive mode.
 
@@ -293,7 +293,7 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
     if vmin is not None and vmax is not None:
         plot_dict['ylim'] = [vmin, vmax]
     
-    if model_data_2d.size == 0 or obs_data_2d.size == 0:
+    if model_data_2d.size == 0 or pairdf.empty:
         print("Warning: Model or observation data is empty. Skipping plot.")
         return
 
@@ -301,15 +301,13 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
     print(f"time shape: {time.shape}, ndims: {time.ndim}")
     print(f"altitude shape: {altitude.shape}, ndims: {altitude.ndim}")
     print(f"model_data_2d shape: {model_data_2d.shape}, ndims: {model_data_2d.ndim}")
-    print(f"obs_data_2d shape: {obs_data_2d.shape}, ndims: {obs_data_2d.ndim}")
+    print(f"pairdf shape: {pairdf.shape}, ndims: {pairdf.ndim}")
 
     # Determine cmin and cmax dynamically if not provided
     if cmin is None or cmax is None:
-        cmin = min(np.nanmin(model_data_2d), np.nanmin(obs_data_2d))
-        cmax = max(np.nanmax(model_data_2d), np.nanmax(obs_data_2d))
+        cmin = min(np.nanmin(model_data_2d), np.nanmin(pairdf[obs_var].values))
+        cmax = max(np.nanmax(model_data_2d), np.nanmax(pairdf[obs_var].values))
     
-    
-
     print(f"cmin and cmax set dynamically based on model and observation data to {cmin} and {cmax}")
 
     fig, axs = plt.subplots(nrows=2, figsize=grp_dict.get('fig_kwargs', {}).get('figsize', (20, 8)), sharex=True, gridspec_kw={'height_ratios': [1, 1]})
@@ -320,7 +318,7 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
         colors = grp_dict['colors']
         cmap = LinearSegmentedColormap.from_list("custom_cmap", colors, N=grp_dict.get('color_levels', 100))
     else:
-    # Default color map logic
+        # Default color map logic
         cmap = plt.get_cmap(grp_dict.get('color_map', 'viridis'))
 
     norm = Normalize(vmin=cmin, vmax=cmax)
@@ -333,35 +331,32 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
 
     # Adjust the layout to place the colorbar at the bottom
     cbar = fig.colorbar(contourf_plot, ax=axs.ravel().tolist(), orientation='horizontal', pad=0.1, aspect=50)
-    cbar.set_label(obs_label_config.get(obs_var, {}).get('ylabel_plot', 'Variable (units)'), fontsize=18, fontweight='bold')
-    cbar.ax.tick_params(labelsize=14)
+    cbar.set_label(obs_label_config.get(obs_var, {}).get('ylabel_plot', 'Variable (units)'), fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))
+    cbar.ax.tick_params(labelsize=text_dict.get('labelsize', 14))
 
     axs[0].xaxis.set_major_formatter(mdates.DateFormatter('%m-%d %H'))
     axs[0].xaxis.set_major_locator(mdates.AutoDateLocator())
     fig.autofmt_xdate(rotation=45, ha='right')
 
-    #axs[0].set_title(f"Model vs Observation Curtain Plot: {model_var} vs {obs_var}", fontsize=12, fontweight='bold')
     # Set the main title and subplot titles
-    fig.suptitle(f"Model vs Observation Curtain Plot: {model_var} vs {obs_var}", fontsize=16, fontweight='bold')
-    axs[0].set_title("Model Curtain with Model Scatter Overlay", fontsize=18, fontweight='bold')
+    fig.suptitle(f"Model vs Observation Curtain Plot: {mod_var} vs {obs_var}", fontsize=text_dict.get('fontsize', 16), fontweight=text_dict.get('fontweight', 'bold'))
+    axs[0].set_title("Model Curtain with Model Scatter Overlay", fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))
     
-    axs[0].set_ylabel('Pressure (Pa)', fontsize=18, fontweight='bold')
+    axs[0].set_ylabel('Pressure (Pa)', fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))
     axs[0].invert_yaxis()  # Invert y-axis to have max pressure at the bottom
-    axs[0].tick_params(axis='both', labelsize=14)
+    axs[0].tick_params(axis='both', labelsize=text_dict.get('labelsize', 14))
 
     # Separate subplot for the observation scatter plot
-    axs[1].set_title("Observation Scatter", fontsize=18, fontweight='bold')
-    scatter = axs[1].scatter(time_dates, obs_pressure, c=obs_data_2d, cmap=cmap, norm=norm, edgecolor=(0, 0, 0, 0), linewidth=0.3, alpha=0.5)
+    axs[1].set_title("Observation Scatter", fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))
+    scatter = axs[1].scatter(time_dates, obs_pressure, c=pairdf[obs_var].values, cmap=cmap, norm=norm, edgecolor=(0, 0, 0, 0), linewidth=0.3, alpha=0.5)
     axs[1].xaxis.set_major_formatter(mdates.DateFormatter('%m-%d %H'))
     axs[1].xaxis.set_major_locator(mdates.AutoDateLocator())
     fig.autofmt_xdate(rotation=45, ha='right')
 
-    axs[1].set_ylabel('Pressure (Pa)', fontsize=18, fontweight='bold')
+    axs[1].set_ylabel('Pressure (Pa)', fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))
     axs[1].invert_yaxis()  # Invert y-axis to have max pressure at the bottom
-    axs[1].tick_params(axis='both', labelsize=14)
-    axs[1].set_xlabel('Time', fontsize=18, fontweight='bold')
-
-    
+    axs[1].tick_params(axis='both', labelsize=text_dict.get('labelsize', 14))
+    axs[1].set_xlabel('Time', fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))
 
     # Set the same y-tick labels for both subplots and exclude the minimum value
     y_ticks = axs[0].get_yticks()
@@ -371,12 +366,8 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
     axs[1].set_yticklabels([str(int(tick)) if tick != y_ticks.min() else '' for tick in y_ticks])
 
     # Plot model scatter data on top of the model curtain plot
-    if model_scatter_data is not None:
-        axs[0].scatter(time_dates, obs_pressure, c=model_scatter_data, cmap=cmap, norm=norm, edgecolor='black', linewidth=0.5, alpha=0.7)
+    axs[0].scatter(time_dates, obs_pressure, c=pairdf[mod_var].values, cmap=cmap, norm=norm, edgecolor='black', linewidth=0.5, alpha=0.7)
 
-        # Print diagnostics for model scatter data
-        print(f"model_scatter_data shape: {model_scatter_data.shape}, ndims: {model_scatter_data.ndim}")
-        
     # Save the curtain plot for the current pair immediately
     print(f"Saving curtain plot to {outname}...")
     savefig(f"{outname}", loc=4, logo_height=100, dpi=300)
@@ -386,6 +377,7 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, obs_data_2d,
     if not debug:
         plt.close()
 
+
     #Diagnostic Histogram
     #plt.figure(figsize=(10, 4))
     #plt.hist(model_data_2d.values.flatten(), bins=100, alpha=0.15, color='blue', label='Model', range=(cmin, cmax))

From a1bd14d3107433cfd156172df0fff542564da588 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Wed, 31 Jul 2024 21:56:54 -0600
Subject: [PATCH 13/20] removal of parts not really used such as, model_data
 and obs_data

---
 melodies_monet/driver.py | 15 +++++----------
 1 file changed, 5 insertions(+), 10 deletions(-)

diff --git a/melodies_monet/driver.py b/melodies_monet/driver.py
index 47fffb68..9a93bba5 100644
--- a/melodies_monet/driver.py
+++ b/melodies_monet/driver.py
@@ -1654,7 +1654,7 @@ def plotting(self):
                             ##print(f"Generated target pressures: {target_pressures}, shape: {target_pressures.shape}")
                         
                             # Check for NaN values before interpolation
-                            ##print(f"NaNs in model_data before interpolation: {np.isnan(ds_model[modvar]).sum().compute()}")
+                            ##print(f"NaNs in model data before interpolation: {np.isnan(ds_model[modvar]).sum().compute()}")
                             ##print(f"NaNs in pressure_model before interpolation: {np.isnan(ds_model['pressure_model']).sum().compute()}")
 
                         
@@ -1670,7 +1670,7 @@ def plotting(self):
                             ds_wrf_const = xr.merge([da_wrf_const, da_target_pressures])
                             ds_wrf_const = ds_wrf_const.set_coords('target_pressures')
                         
-                            # Debugging: print merged dataset
+                            # Debugging: print merged dataset for model curtain
                             ##print(ds_wrf_const)
                         
                             # Ensure model_data_2d is properly reshaped for the contourfill plot
@@ -1684,17 +1684,12 @@ def plotting(self):
 
                             # Fetch model pressure and other model and observation data from "pairdf" (for scatter plot overlay)
                             time = pairdf['time']
-                            obs_pressure = pairdf['pressure_obs']
-                            obs_data = pairdf[obsvar]
-                            model_data = pairdf[modvar]        
+                            obs_pressure = pairdf['pressure_obs']  
+                            ##print(f"Length of time: {len(time)}") #Debugging
+                            ##print(f"Length of obs_pressure: {len(obs_pressure)}") #Debugging
                         
                             # Generate the curtain plot using airplots.make_curtain_plot
                             try:
-                                ##print(f"Length of time_dates: {len(time_dates)}")
-                                ##print(f"Length of obs_pressure: {len(obs_pressure)}")
-                                ##print(f"Length of obs_data_2d: {len(obs_data_2d)}")
-                                ##print(f"model_scatter_data shape: {model_scatter_data.shape}")
-
                                 outname_pair = f"{outname}_{obs_label}_vs_{model_label}.png"
 
                                 print(f"Saving curtain plot to {outname_pair}...")

From 0ba4156a43ab7083d9e1acde2b53b92200ef8356 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Tue, 17 Sep 2024 16:49:57 -0600
Subject: [PATCH 14/20] Resolved issues remaining with yaxis (altitude variable
 ie. pressure) limits, Tested with unit change Pa to hPa, cleaned up codes and
 yaml

---
 ...elop_aircraft_testmultiplemodelgroups.yaml | 25 +++----
 melodies_monet/driver.py                      | 46 +++++++------
 melodies_monet/plots/aircraftplots.py         | 66 +++++++++++++++----
 3 files changed, 92 insertions(+), 45 deletions(-)

diff --git a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
index df7b44ee..8e774d99 100644
--- a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
+++ b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft_testmultiplemodelgroups.yaml
@@ -25,8 +25,8 @@ model:
     variables:
         'pres_pa_mid':
             rename: pressure_model # name to convert this variable to, use 'pressure_model' for aircraft obs
-            unit_scale: 1 #Opt Scaling factor 
-            unit_scale_method: '*' #Opt Multiply = '*' , Add = '+', subtract = '-', divide = '/'
+            unit_scale: 100 #1 #Opt Scaling factor 
+            unit_scale_method: '/' #convert to hPa (1 hPa=100 Pa) #'*' #Opt Multiply = '*' , Add = '+', subtract = '-', divide = '/'
         'temperature_k':
             rename: temp_model # name to convert this variable to
             unit_scale: 1 #Opt Scaling factor 
@@ -52,8 +52,8 @@ model:
     variables:
         'pres_pa_mid':
             rename: pressure_model # name to convert this variable to, use 'pressure_model' for aircraft obs
-            unit_scale: 1 #Opt Scaling factor 
-            unit_scale_method: '*' #Opt Multiply = '*' , Add = '+', subtract = '-', divide = '/'
+            unit_scale: 100 #1 #Opt Scaling factor 
+            unit_scale_method: '/' #convert to hPa (1 hPa=100 Pa) #'*' #Opt Multiply = '*' , Add = '+', subtract = '-', divide = '/'
         'temperature_k':
             rename: temp_model # name to convert this variable to
             unit_scale: 1 #Opt Scaling factor 
@@ -118,7 +118,7 @@ obs:
             unit_scale_method: '*' #Opt Multiply = '*' , Add = '+', subtract = '-', divide = '/'
       'P_BUI':
             rename: pressure_obs # name to convert this variable to
-            unit_scale: 100 #Opt Scaling factor 
+            unit_scale: 1 #remain in hPa #100 # hpa to Pa  #Opt Scaling factor 
             unit_scale_method: '*' #Opt Multiply = '*' , Add = '+', subtract = '-', divide = '/'
       'MSL_GPS_Altitude_YANG':
           rename: altitude # name to convert this variable to
@@ -255,19 +255,20 @@ plots:
       set_axis: True
       altitude_variable: 'altitude'
     ##To use a custom Matplotlib colormap, set “color_map_custom”=True and specify “colors” and "color_levels" like the example below. To use a standard Matplotlib colormap, set “color_map_custom” = False and specify a “color_map”:
-    color_map_custom: True #False #True
+    color_map_custom: True #False 
     colors: ["#ff8cff", "#dd6ff2", "#bb52e5", "#9935d8", "#7718cb",
         "#0000bb", "#002ccc", "#0058dd", "#0084ee", "#00afff",
         "#00ebff", "#27ffd7", "#63ff9b", "#a3ff5b", "#d3ff2b",
         "#ffff00", "#ffcf00", "#ff9f00", "#ff6f00", "#ff3f00",
         "#ff0000", "#d8000f", "#b2001f", "#8c002f", "#66003f",
-        "#343434", "#606060", "#8c8c8c", "#b8b8b8", "#e4e4e4"] #,
-        #"#ffffff"] #['blue','royalblue', 'cyan', 'yellow', 'orange', 'red'] ##['#440154', '#3b528b', '#21908d', '#5dc963', '#fde725'] #['purple', 'blue', 'green', 'yellow'] # Example gradient
-    color_levels: 30 #31 #15 #10  # Define the number of distinct colors in the color bar
+        "#343434", "#606060", "#8c8c8c", "#b8b8b8", "#e4e4e4"] # Example gradient
+    color_levels: 30   # Define the number of distinct colors in the color bar (see te numbers of colors in 'colors' key above is matching)
     #color_map: 'Spectral_r' #'jet'  # This can be set instead of color_map_custom, color_levels and color_levels to set a colormap defined in matplotlob.  
-    vmin: 120000 #Pressure in Pa  #0  #Optional (y-axis limits)
-    vmax: 5000 # Optional #need to be edited as per min and max altitude (i.e., vmin and vmax)
-    interval : 500 #Pressure in Pa; target pressure interpolation steps for model curtain 2d plot
+    vmin: 50 #in hPa #5000 #120000 #Pressure in Pa  #0  #Optional (y-axis limits)
+    vmax: 1200 # in hPa #120000 #5000 # Optional #need to be edited as per min and max altitude (i.e., vmin and vmax)
+    num_levels: 100   # Number of vertical levels for interpolation
+    interval: 100 #hPa #10000 #Pa   # Y-axis tick interval in Pa (e.g., ticks every 10,000 Pa)
+    pressure_units: 'hPa' #'Pa' (units for y-axis)
     
   plot_grp6:
     type: 'taylor' # plot type
diff --git a/melodies_monet/driver.py b/melodies_monet/driver.py
index 9a93bba5..599f7558 100644
--- a/melodies_monet/driver.py
+++ b/melodies_monet/driver.py
@@ -1203,9 +1203,6 @@ def plotting(self):
         """
         
         from .util.tools import resample_stratify
-        from .util.tools import vert_interp
-        import matplotlib.dates as mdates
-        import numpy as np
         import matplotlib.pyplot as plt
         pair_keys = list(self.paired.keys())
         if self.paired[pair_keys[0]].type.lower() in ['sat_grid_clm','sat_swath_clm']:
@@ -1593,17 +1590,27 @@ def plotting(self):
 
 
                         elif plot_type.lower() == 'curtain':
+                            # Set cmin and cmax from obs_plot_dict for colorbar limits
                             if set_yaxis:
                                 if all(k in obs_plot_dict for k in ('vmin_plot', 'vmax_plot')):
-                                    vmin = obs_plot_dict['vmin_plot']
-                                    vmax = obs_plot_dict['vmax_plot']
+                                    cmin = obs_plot_dict['vmin_plot']
+                                    cmax = obs_plot_dict['vmax_plot']
                                 else:
                                     print('Warning: vmin_plot and vmax_plot not specified for ' + obsvar + ', so default used.')
-                                    vmin = None
-                                    vmax = None
+                                    cmin = None
+                                    cmax = None
+                            else:
+                                cmin = None
+                                cmax = None
+                            
+                            # Set vmin and vmax from grp_dict for altitude limits
+                            if set_yaxis:
+                                vmin = grp_dict.get('vmin', None)
+                                vmax = grp_dict.get('vmax', None)
                             else:
                                 vmin = None
                                 vmax = None
+
                                 
                             curtain_config = grp_dict # Curtain plot grp YAML dict
                             # Inside your loop for processing each pair
@@ -1621,10 +1628,6 @@ def plotting(self):
                         
                             # Fetch the model and observation data from pairdf
                             pairdf = pairdf_all.reset_index()
-
-                            # Determine cmin and cmax (colorbar min max) from observation config if provided
-                            cmin = obs_plot_dict.get('vmin_plot', None)
-                            cmax = obs_plot_dict.get('vmax_plot', None)
                         
                             #### For model_data_2d for curtain/contourfill plot #####                       
                             # Convert to get something useful for MONET
@@ -1643,12 +1646,15 @@ def plotting(self):
                             # Define target pressures for interpolation based on the range of pressure_model
                             min_pressure = ds_model['pressure_model'].min().compute()
                             max_pressure = ds_model['pressure_model'].max().compute()
-                            # Fetch the interval from curtain_config
-                            interval = curtain_config.get('interval', 100)  # Default to 100 (in Pa) if not provided 
-                        
-                            print(f"Pressure MIN:{min_pressure}, max: {max_pressure}, interval: {interval}")
+                            
+                            # Fetch the interval and num_levels from curtain_config
+                            interval = curtain_config.get('interval', 10000)  # Default to 10,000 Pa if not provided      # Y-axis tick interval
+                            num_levels = curtain_config.get('num_levels', 100)   # Default to 100 levels if not provided
                         
-                            target_pressures = np.linspace(max_pressure, min_pressure, interval)
+                            print(f"Pressure MIN:{min_pressure}, max: {max_pressure}, ytick_interval: {interval}, interpolation_levels: {num_levels}  ")
+                            
+                            # Use num_levels to define target_pressures interpolation levels 
+                            target_pressures = np.linspace(max_pressure, min_pressure, num_levels)
                             
                             # Debugging: print target pressures
                             ##print(f"Generated target pressures: {target_pressures}, shape: {target_pressures.shape}")
@@ -1698,17 +1704,17 @@ def plotting(self):
                                     time=pd.to_datetime(time),
                                     altitude=target_pressures,  # Use target_pressures for interpolation
                                     model_data_2d=model_data_2d,  # Already reshaped to match the expected shape
-                                    pairdf=pairdf,  #use pairdf for scatter overlay (model and obs)
                                     obs_pressure=obs_pressure,  # Pressure_obs for obs scatter plot
+                                    pairdf=pairdf,  #use pairdf for scatter overlay (model and obs)
                                     mod_var=modvar,
                                     obs_var=obsvar,
-                                    cmin=cmin,
-                                    cmax=cmax,
+                                    grp_dict=curtain_config,
                                     vmin=vmin,
                                     vmax=vmax,
+                                    cmin=cmin,
+                                    cmax=cmax,
                                     plot_dict=plot_dict,
                                     outname=outname_pair,
-                                    grp_dict=curtain_config,
                                     domain_type=domain_type,
                                     domain_name=domain_name,
                                     obs_label_config=obs_label_config,
diff --git a/melodies_monet/plots/aircraftplots.py b/melodies_monet/plots/aircraftplots.py
index 05a3c3d3..3f14bd58 100755
--- a/melodies_monet/plots/aircraftplots.py
+++ b/melodies_monet/plots/aircraftplots.py
@@ -289,10 +289,7 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, pairdf, mod_v
     None
         Saves the generated plot to a file.
     """
-
-    if vmin is not None and vmax is not None:
-        plot_dict['ylim'] = [vmin, vmax]
-    
+  
     if model_data_2d.size == 0 or pairdf.empty:
         print("Warning: Model or observation data is empty. Skipping plot.")
         return
@@ -342,9 +339,59 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, pairdf, mod_v
     fig.suptitle(f"Model vs Observation Curtain Plot: {mod_var} vs {obs_var}", fontsize=text_dict.get('fontsize', 16), fontweight=text_dict.get('fontweight', 'bold'))
     axs[0].set_title("Model Curtain with Model Scatter Overlay", fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))
     
-    axs[0].set_ylabel('Pressure (Pa)', fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))
+    ##axs[0].set_ylabel('Pressure (Pa)', fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold')) #removed explicit y-axis label (made it flexible: see pressure_units via yaml)
     axs[0].invert_yaxis()  # Invert y-axis to have max pressure at the bottom
     axs[0].tick_params(axis='both', labelsize=text_dict.get('labelsize', 14))
+    axs[1].invert_yaxis()  # Invert y-axis to have max pressure at the bottom (FOR SECOND SUBPLOT)
+
+    # Set y-axis limits if vmin and vmax are provided
+    if vmin is not None and vmax is not None:
+        # Since y-axis is inverted already (see invert_yaxis() above), set limits accordingly
+        print(f"Setting y-axis limits: vmin={vmin}, vmax={vmax}")
+        axs[0].set_ylim(vmax, vmin)  # Note the order to match inversion
+        axs[1].set_ylim(vmax, vmin)  # Note the order to match inversion
+    else:
+        vmin, vmax = axs[0].get_ylim()
+
+
+
+    # Retrieve pressure units from grp_dict or default to 'Pa'
+    pressure_units = grp_dict.get('pressure_units', 'Pa')
+    
+    # Set y-tick labels and y-axis label based on pressure units
+    if pressure_units == 'hPa':
+        y_axis_label = 'Pressure (hPa)'
+    else:
+        y_axis_label = 'Pressure (Pa)'
+    
+    # Apply y-tick labels and y-axis labels (both subplots)
+    axs[0].set_ylabel(y_axis_label, fontsize=text_dict.get('fontsize', 18),
+                      fontweight=text_dict.get('fontweight', 'bold'))
+    axs[1].set_ylabel(y_axis_label, fontsize=text_dict.get('fontsize', 18),
+                      fontweight=text_dict.get('fontweight', 'bold'))
+   
+
+    # Set y-axis ticks at specified intervals for axs[0] (first subplot)
+    if 'interval' in grp_dict:
+        interval = grp_dict['interval']
+        print(f"Interval value: {interval}")  # This will print the interval value
+        y_ticks = np.arange(vmin, vmax + interval, interval)
+        y_ticks = y_ticks[y_ticks <= vmax]  # Ensure ticks do not exceed vmax
+        print(f"Calculated y_ticks: {y_ticks}")
+        axs[0].set_yticks(y_ticks)
+        # Format y-tick labels
+        axs[0].set_yticklabels([str(int(tick)) for tick in y_ticks])
+    else:
+        y_ticks = axs[0].get_yticks()
+    
+    # Synchronize y-ticks and y-limits for axs[1] (Second subplot) based on axs[0] (First subplot)
+    y_ticks = axs[0].get_yticks()
+    y_lim = axs[0].get_ylim()
+    axs[1].set_yticks(y_ticks)
+    axs[1].set_ylim(y_lim)
+    # Set y-tick labels for axs[1]
+    axs[1].set_yticklabels([str(int(tick)) for tick in y_ticks])  ##axs[1].set_yticklabels([str(int(tick)) if tick != y_ticks.min() else '' for tick in y_ticks]) #if needed to hide minimum 
+
 
     # Separate subplot for the observation scatter plot
     axs[1].set_title("Observation Scatter", fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))
@@ -353,17 +400,10 @@ def make_curtain_plot(time, altitude, model_data_2d, obs_pressure, pairdf, mod_v
     axs[1].xaxis.set_major_locator(mdates.AutoDateLocator())
     fig.autofmt_xdate(rotation=45, ha='right')
 
-    axs[1].set_ylabel('Pressure (Pa)', fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))
-    axs[1].invert_yaxis()  # Invert y-axis to have max pressure at the bottom
+    ##axs[1].set_ylabel('Pressure (Pa)', fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))  #made it flexible via YAML specified 'pressure_units'
     axs[1].tick_params(axis='both', labelsize=text_dict.get('labelsize', 14))
     axs[1].set_xlabel('Time', fontsize=text_dict.get('fontsize', 18), fontweight=text_dict.get('fontweight', 'bold'))
 
-    # Set the same y-tick labels for both subplots and exclude the minimum value
-    y_ticks = axs[0].get_yticks()
-    y_lim = axs[0].get_ylim()
-    axs[1].set_yticks(y_ticks)
-    axs[1].set_ylim(y_lim)
-    axs[1].set_yticklabels([str(int(tick)) if tick != y_ticks.min() else '' for tick in y_ticks])
 
     # Plot model scatter data on top of the model curtain plot
     axs[0].scatter(time_dates, obs_pressure, c=pairdf[mod_var].values, cmap=cmap, norm=norm, edgecolor='black', linewidth=0.5, alpha=0.7)

From b9ba0ef99e383d0c0ef6065ed3026f545b049040 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Tue, 17 Sep 2024 17:00:12 -0600
Subject: [PATCH 15/20] Added curtain plot options to the single model vs obs
 yaml as well

---
 ...m_aircraft_Latestfor_develop_aircraft.yaml | 40 ++++++++++++++++++-
 1 file changed, 39 insertions(+), 1 deletion(-)

diff --git a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft.yaml b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft.yaml
index 5280f4c7..317dab7f 100644
--- a/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft.yaml
+++ b/examples/yaml/control_wrfchem_aircraft_Latestfor_develop_aircraft.yaml
@@ -124,6 +124,7 @@ plots:
           altitude:
             oper: "between"
             value: [2000,8000] # values are [vmin_y2, vmax_y2] 
+            
   plot_grp2:
     type: 'vertprofile'  # plot type
     fig_kwargs:  # Opt to define figure options
@@ -198,7 +199,42 @@ plots:
     # Add other seaborn.kdeplot keyword arguments here as needed
     vcenter: #0  # Central value for TwoSlopeNorm
     extensions: ['min', 'max']  # Extensions for the colorbar  
+
   plot_grp5:
+    type: 'curtain'
+    fig_kwargs:
+      figsize: [20,14] #[12, 8]
+    default_plot_kwargs:
+      #linewidth: 4.0
+      markersize: 40.
+    text_kwargs:
+      fontsize: 25 #18
+      fontweight: 'bold'
+      labelsize: 16
+    domain_type: ['all']
+    domain_name: ['CONUS']
+    data: ['firexaq_wrfchem_v4.2','firexaq_wrfchem_v4.2_test']
+    data_proc:
+      rem_obs_nan: True
+      set_axis: True
+      altitude_variable: 'altitude'
+    ##To use a custom Matplotlib colormap, set “color_map_custom”=True and specify “colors” and "color_levels" like the example below. To use a standard Matplotlib colormap, set “color_map_custom” = False and specify a “color_map”:
+    color_map_custom: True #False 
+    colors: ["#ff8cff", "#dd6ff2", "#bb52e5", "#9935d8", "#7718cb",
+        "#0000bb", "#002ccc", "#0058dd", "#0084ee", "#00afff",
+        "#00ebff", "#27ffd7", "#63ff9b", "#a3ff5b", "#d3ff2b",
+        "#ffff00", "#ffcf00", "#ff9f00", "#ff6f00", "#ff3f00",
+        "#ff0000", "#d8000f", "#b2001f", "#8c002f", "#66003f",
+        "#343434", "#606060", "#8c8c8c", "#b8b8b8", "#e4e4e4"] # Example gradient
+    color_levels: 30   # Define the number of distinct colors in the color bar (see te numbers of colors in 'colors' key above is matching)
+    #color_map: 'Spectral_r' #'jet'  # This can be set instead of color_map_custom, color_levels and color_levels to set a colormap defined in matplotlob.  
+    vmin: 50 #in hPa #5000 #120000 #Pressure in Pa  #0  #Optional (y-axis limits)
+    vmax: 1200 # in hPa #120000 #5000 # Optional #need to be edited as per min and max altitude (i.e., vmin and vmax)
+    num_levels: 100   # Number of vertical levels for interpolation
+    interval: 100 #hPa #10000 #Pa   # Y-axis tick interval in Pa (e.g., ticks every 10,000 Pa)
+    pressure_units: 'hPa' #'Pa' (units for y-axis)
+    
+  plot_grp6:
     type: 'taylor' # plot type
     fig_kwargs: #Opt to define figure options
       figsize: [8,8] # figure size if multiple plots
@@ -213,7 +249,8 @@ plots:
     data_proc:
       rem_obs_nan: True # True: Remove all points where model or obs variable is NaN. False: Remove only points where model variable is NaN.
       set_axis: True #If select True, add ty_scale for each variable in obs.
-  plot_grp6:
+      
+  plot_grp7:
     type: 'boxplot' # plot type
     fig_kwargs: #Opt to define figure options
       figsize: [8, 6] # figure size 
@@ -225,6 +262,7 @@ plots:
     data_proc:
       rem_obs_nan: True # True: Remove all points where model or obs variable is NaN. False: Remove only points where model variable is NaN.
       set_axis: False #If select True, add vmin_plot and vmax_plot for each variable in obs.
+      
 stats:
   #Stats require positive numbers, so if you want to calculate temperature use Kelvin!
   #Wind direction has special calculations for AirNow if obs name is 'WD'

From d53037998bedff4b7f492d7629912f39f65d198c Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Wed, 18 Sep 2024 12:54:16 -0600
Subject: [PATCH 16/20] docs/develop/developers_guide.rst Edits for fixing
 MELODIES-MONET project board link with text edit works

---
 docs/develop/developers_guide.rst | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/docs/develop/developers_guide.rst b/docs/develop/developers_guide.rst
index f75a0637..0a005103 100644
--- a/docs/develop/developers_guide.rst
+++ b/docs/develop/developers_guide.rst
@@ -172,7 +172,7 @@ The generated HTML will be created in ``docs/_build/html``,
 with ``docs/_build/html/index.html`` the main page that can be
 viewed in any browser.
 
-Please see the `Documentation <https://github.com/NOAA-CSL/MELODIES-MONET/projects/2>`_ 
+Please see the `MELODIES-MONET Documentation <https://github.com/orgs/NOAA-CSL/projects/6>`_ 
 project on GitHub to learn about current and future development.
 
 

From 79ca49e306893a6ba1a731e4f02e57092dd41840 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Wed, 18 Sep 2024 13:10:34 -0600
Subject: [PATCH 17/20] TEST docs/background/supported_plots.rst Edits for
 fixing MELODIES-MONET project board link with text edit works

---
 docs/background/supported_plots.rst | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/docs/background/supported_plots.rst b/docs/background/supported_plots.rst
index 5cb0b7d5..5ad1f095 100644
--- a/docs/background/supported_plots.rst
+++ b/docs/background/supported_plots.rst
@@ -5,7 +5,7 @@ Model to Model Comparisons
 --------------------------
 Under development. 
 
-See the `Spatial Verification <https://github.com/NOAA-CSL/MELODIES-MONET/projects/3>`_ 
+.. See the `Spatial Verification <https://github.com/NOAA-CSL/MELODIES-MONET/projects/3>`_ 
 project on GitHub to learn about current and future development.
 
 Model to Observation Comparisons

From efbfe90ae27b5cd6195731e4dcaab3abd88c7d6e Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Wed, 18 Sep 2024 13:14:15 -0600
Subject: [PATCH 18/20] Revert "TEST docs/background/supported_plots.rst Edits
 for fixing MELODIES-MONET project board link with text edit works"

This reverts commit 79ca49e306893a6ba1a731e4f02e57092dd41840.
---
 docs/background/supported_plots.rst | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/docs/background/supported_plots.rst b/docs/background/supported_plots.rst
index 5ad1f095..5cb0b7d5 100644
--- a/docs/background/supported_plots.rst
+++ b/docs/background/supported_plots.rst
@@ -5,7 +5,7 @@ Model to Model Comparisons
 --------------------------
 Under development. 
 
-.. See the `Spatial Verification <https://github.com/NOAA-CSL/MELODIES-MONET/projects/3>`_ 
+See the `Spatial Verification <https://github.com/NOAA-CSL/MELODIES-MONET/projects/3>`_ 
 project on GitHub to learn about current and future development.
 
 Model to Observation Comparisons

From 0c9778159d6c60efeb12609c1425cbfe5b1d1003 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Wed, 18 Sep 2024 13:37:11 -0600
Subject: [PATCH 19/20] Renamed docs/background to docs/users_guide

---
 docs/{background => users_guide}/description.rst        | 0
 docs/{background => users_guide}/gridded_datasets.rst   | 0
 docs/{background => users_guide}/introduction.rst       | 0
 docs/{background => users_guide}/supported_analyses.rst | 0
 docs/{background => users_guide}/supported_datasets.rst | 0
 docs/{background => users_guide}/supported_plots.rst    | 0
 docs/{background => users_guide}/supported_stats.rst    | 0
 docs/{background => users_guide}/time_chunking.rst      | 0
 8 files changed, 0 insertions(+), 0 deletions(-)
 rename docs/{background => users_guide}/description.rst (100%)
 rename docs/{background => users_guide}/gridded_datasets.rst (100%)
 rename docs/{background => users_guide}/introduction.rst (100%)
 rename docs/{background => users_guide}/supported_analyses.rst (100%)
 rename docs/{background => users_guide}/supported_datasets.rst (100%)
 rename docs/{background => users_guide}/supported_plots.rst (100%)
 rename docs/{background => users_guide}/supported_stats.rst (100%)
 rename docs/{background => users_guide}/time_chunking.rst (100%)

diff --git a/docs/background/description.rst b/docs/users_guide/description.rst
similarity index 100%
rename from docs/background/description.rst
rename to docs/users_guide/description.rst
diff --git a/docs/background/gridded_datasets.rst b/docs/users_guide/gridded_datasets.rst
similarity index 100%
rename from docs/background/gridded_datasets.rst
rename to docs/users_guide/gridded_datasets.rst
diff --git a/docs/background/introduction.rst b/docs/users_guide/introduction.rst
similarity index 100%
rename from docs/background/introduction.rst
rename to docs/users_guide/introduction.rst
diff --git a/docs/background/supported_analyses.rst b/docs/users_guide/supported_analyses.rst
similarity index 100%
rename from docs/background/supported_analyses.rst
rename to docs/users_guide/supported_analyses.rst
diff --git a/docs/background/supported_datasets.rst b/docs/users_guide/supported_datasets.rst
similarity index 100%
rename from docs/background/supported_datasets.rst
rename to docs/users_guide/supported_datasets.rst
diff --git a/docs/background/supported_plots.rst b/docs/users_guide/supported_plots.rst
similarity index 100%
rename from docs/background/supported_plots.rst
rename to docs/users_guide/supported_plots.rst
diff --git a/docs/background/supported_stats.rst b/docs/users_guide/supported_stats.rst
similarity index 100%
rename from docs/background/supported_stats.rst
rename to docs/users_guide/supported_stats.rst
diff --git a/docs/background/time_chunking.rst b/docs/users_guide/time_chunking.rst
similarity index 100%
rename from docs/background/time_chunking.rst
rename to docs/users_guide/time_chunking.rst

From 51cbf21099e6ad8a3af9eafbca8478879170bca8 Mon Sep 17 00:00:00 2001
From: quaz115 <quazi.rasool@noaa.gov>
Date: Wed, 18 Sep 2024 14:09:22 -0600
Subject: [PATCH 20/20] Renamed supported_analyses.rst to
 supported_diagnostics.rst and updated content

---
 .../{supported_analyses.rst => supported_diagnostics.rst}       | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)
 rename docs/users_guide/{supported_analyses.rst => supported_diagnostics.rst} (98%)

diff --git a/docs/users_guide/supported_analyses.rst b/docs/users_guide/supported_diagnostics.rst
similarity index 98%
rename from docs/users_guide/supported_analyses.rst
rename to docs/users_guide/supported_diagnostics.rst
index 618fec9e..f774c4ca 100644
--- a/docs/users_guide/supported_analyses.rst
+++ b/docs/users_guide/supported_diagnostics.rst
@@ -1,4 +1,4 @@
-Supported Analyses
+Supported Diagnostics
 ==================
 
 Supported data analysis options in MELODIES MONET are explained below.