From d90cb674d4d77bbce7fbdc4d063f6975ae6bad45 Mon Sep 17 00:00:00 2001
From: Axel Lauer <axel.lauer@dlr.de>
Date: Tue, 17 Jan 2023 13:12:28 +0100
Subject: [PATCH 1/5] added support for ICON to lauer22jclim cloud diagnostics

---
 doc/sphinx/source/recipes/recipe_clouds.rst   |  6 +-
 esmvaltool/diag_scripts/clouds/clouds.ncl     | 41 +++++++++-----
 .../diag_scripts/clouds/clouds_dyn_matrix.ncl | 31 ++++++++---
 .../clouds/clouds_interannual.ncl             | 22 +++++---
 esmvaltool/diag_scripts/clouds/clouds_pdf.ncl | 55 ++++++++++---------
 .../diag_scripts/clouds/clouds_scatter.ncl    | 22 ++++++--
 .../clouds/clouds_seasonal_cycle.ncl          | 22 +++++---
 .../diag_scripts/clouds/clouds_taylor.ncl     | 19 ++++---
 .../diag_scripts/clouds/clouds_zonal.ncl      | 20 ++++++-
 9 files changed, 155 insertions(+), 83 deletions(-)

diff --git a/doc/sphinx/source/recipes/recipe_clouds.rst b/doc/sphinx/source/recipes/recipe_clouds.rst
index 6b606fa83a..9328c1af3c 100644
--- a/doc/sphinx/source/recipes/recipe_clouds.rst
+++ b/doc/sphinx/source/recipes/recipe_clouds.rst
@@ -109,6 +109,8 @@ User settings in recipe
      "CylindricalEquidistant")
    * showdiff: calculate and plot differences model - reference
      (default = false)
+   * showyears: add start and end years to the plot titles
+     (default = false)
    * rel_diff: if showdiff = true, then plot relative differences (%)
      (default = False)
    * ref_diff_min: lower cutoff value in case of calculating relative
@@ -297,7 +299,7 @@ User settings in recipe
 
    none
 
-8. Script clouds_lifrac_pdf.ncl
+8. Script clouds_pdf.ncl
 
    *Required settings (scripts)*
 
@@ -333,6 +335,8 @@ User settings in recipe
    * explicit_cn_levels: use these contour levels for plotting
    * filename_add: optionally add this string to plot filesnames
    * projection: map projection, e.g., Mollweide, Mercator
+   * showyears: add start and end years to the plot titles
+     (default = false)
    * var: short_name of variable to process (default = "" i.e. use
      first variable in variable list)
 
diff --git a/esmvaltool/diag_scripts/clouds/clouds.ncl b/esmvaltool/diag_scripts/clouds/clouds.ncl
index 9f933d3e89..928d7cc96d 100644
--- a/esmvaltool/diag_scripts/clouds/clouds.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds.ncl
@@ -30,7 +30,9 @@
 ;                       plot(s)
 ;   projection:         map projection for plotting (default =
 ;                       "CylindricalEquidistant")
-;   showdiff            calculate and plot differences (default = False)
+;   showdiff:           calculate and plot differences (default = False)
+;   showyears:          add start and end years to the plot titles (default =
+;                       false)
 ;   rel_diff:           if showdiff = True, then plot relative differences (%)
 ;                       (default = False)
 ;   rel_diff_min:       lower cutoff value in case of calculating relative
@@ -58,6 +60,7 @@
 ;   none
 ;
 ; Modification history
+;   20230117-lauer_axel: added support for ICON (code from Manuel)
 ;   20211021-lauer_axel: added output of basic statistics as ascii files
 ;   20211006-lauer_axel: removed write_plots
 ;   20210325-lauer_axel: added option to estimate observational uncertainty
@@ -98,6 +101,7 @@ begin
   set_default_att(diag_script_info, "rel_diff", False)
   set_default_att(diag_script_info, "rel_diff_min", -1.0e19)
   set_default_att(diag_script_info, "showdiff", False)
+  set_default_att(diag_script_info, "showyears", False)
   set_default_att(diag_script_info, "timemean", "annualclim")
   set_default_att(diag_script_info, "treat_var_as_error", False)
   set_default_att(diag_script_info, "var", "")
@@ -109,17 +113,18 @@ begin
   end if
 
   variables = metadata_att_as_array(variable_info, "short_name")
-  varidx = ind(variables .eq. var0)
-  if (ismissing(varidx)) then
-    errstr = "diagnostic " + diag + " requires the following variable: var0"
+  if (.not. any(variables .eq. var0)) then
+    errstr = "diagnostic " + diag + " requires the following variable: " + var0
     error_msg("f", DIAG_SCRIPT, "", errstr)
   end if
 
+  var0_info = select_metadata_by_name(variable_info, var0)
+  var0_info := var0_info[0]
   info0 = select_metadata_by_name(input_file_info, var0)
   dim_MOD = ListCount(info0)
 
-  if (isatt(variable_info[varidx], "reference_dataset")) then
-    refname = variable_info[varidx]@reference_dataset
+  if (isatt(var0_info, "reference_dataset")) then
+    refname = var0_info@reference_dataset
   end if
   names = metadata_att_as_array(info0, "dataset")
   projects = metadata_att_as_array(info0, "project")
@@ -453,8 +458,18 @@ begin
 
     ; annotation
 
+    ; if desired, add years to plot title
+    years_str = ""
+    if (diag_script_info@showyears) then
+      years_str = " (" + var0_info@start_year
+      if (var0_info@start_year .ne. var0_info@end_year) then
+        years_str = years_str + "-" + var0_info@end_year
+      end if
+      years_str = years_str + ")"
+    end if
+
 ;    res@tiMainOn             = False
-    res@tiMainString = names(imod)
+    res@tiMainString = names(imod) + years_str
     res@tiMainFontHeightF = 0.025
 
     res@cnLevelSelectionMode = "ExplicitLevels"
@@ -465,8 +480,8 @@ begin
 
     ; variable specific plotting settings
 
-    if (isatt(variable_info[varidx], "units")) then
-      data1@units = variable_info[varidx]@units
+    if (isatt(var0_info, "units")) then
+      data1@units = var0_info@units
     else
       data1@units = ""
     end if
@@ -576,8 +591,8 @@ begin
       data1@diag_script = (/DIAG_SCRIPT/)
     end if
 
-    if (isatt(variable_info[varidx], "long_name")) then
-      data1@long_name = variable_info[varidx]@long_name
+    if (isatt(var0_info, "long_name")) then
+      data1@long_name = var0_info@long_name
     end if
 
     data1@var = var0
@@ -784,7 +799,7 @@ begin
         dres@mpPerimOn      = perim    ; draw line around map
         dres@gsnStringFontHeightF = 0.02
 
-        dres@tiMainString =  names(imod) + " - " + refname
+        dres@tiMainString =  names(imod) + " - " + refname + years_str
         dres@tiMainFontHeightF = 0.025
 
         rmsd(imod, :) = rmsd@_FillValue
@@ -818,7 +833,7 @@ begin
         dres@gsnLeftString  = ""
         dres@gsnRightString = ""
         dres@gsnCenterString = ""
-        dres@tiMainString =  refname + " uncertainty"
+        dres@tiMainString =  refname + " uncertainty" + years_str
         rmsd(imod, :) = rmsd@_FillValue
         bias(imod, :) = bias@_FillValue
       else
diff --git a/esmvaltool/diag_scripts/clouds/clouds_dyn_matrix.ncl b/esmvaltool/diag_scripts/clouds/clouds_dyn_matrix.ncl
index 972847c427..9fcfadebbb 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_dyn_matrix.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_dyn_matrix.ncl
@@ -37,6 +37,7 @@
 ;   none
 ;
 ; Modification history
+;   20230117-lauer_axel: added support for ICON (code from Manuel)
 ;   20220126-lauer_axel: added optional variable labels for x- and y-axes
 ;   20211118-lauer_axel: added output of frequency distributions
 ;   20210607-lauer_axel: added multi-model-average (= average over all models)
@@ -57,7 +58,8 @@ begin
   enter_msg(DIAG_SCRIPT, "")
 
   diag = "clouds_dyn_matrix.ncl"
-  variables = metadata_att_as_array(variable_info, "short_name")
+  variables = get_unique_values(metadata_att_as_array(variable_info, \
+                                                      "short_name"))
 
   ; Check required diag_script_info attributes
   exit_if_missing_atts(diag_script_info, (/"var_x", "var_y", "var_z", \
@@ -162,21 +164,26 @@ begin
   ; is equal for each variable
 
   do i = 0, nVAR - 1
-    if (isatt(variable_info[idx(i)], "reference_dataset")) then
-      refname(i) = variable_info[idx(i)]@reference_dataset
+    var = variables(idx(i))
+    var_info = select_metadata_by_name(variable_info, var)
+    var_info := var_info[0]
+    if (isatt(var_info, "reference_dataset")) then
+      refname(i) = var_info@reference_dataset
     end if
-    info = select_metadata_by_name(input_file_info, variables(idx(i)))
+    info = select_metadata_by_name(input_file_info, var)
     if (i .eq. 0) then
       dim_MOD = ListCount(info)
     else
       dim_test = ListCount(info)
       if (dim_test .ne. dim_MOD) then
         error_msg("f", DIAG_SCRIPT, "", "number of datasets for variable " \
-                  + variables(i) + " does not match number of datasets for " \
-                  + variables(0))
+                  + var + " does not match number of datasets for " \
+                  + variables(idx(0)))
       end if
     end if
     delete(info)
+    delete(var)
+    delete(var_info)
   end do
 
   ; Set default values for non-required diag_script_info attributes
@@ -251,9 +258,11 @@ begin
 
   idxmod = get_mod(names_x, projects_x)
 
-  if (idxmod(0) .eq. -1) then
+  if (idxmod(0) .eq. -1) then  ; no model found
     flag_multimod = False
-  else
+  elseif (dimsizes(idxmod) .eq. 1) then  ; one model found
+    flag_multimod = False
+  else  ; more than one model found
     flag_multimod = True
   end if
 
@@ -814,7 +823,11 @@ begin
 
   ref = where(count(ref_ind, :, :) .gt. 1.0e-3, count(ref_ind, :, :), \
               count@_FillValue)
-  modavg = dim_avg_n_Wrap(count(idxmod, :, :), 0)
+  if (dimsizes(idxmod) .gt. 1) then
+    modavg = dim_avg_n_Wrap(count(idxmod, :, :), 0)
+  else
+    modavg = count(idxmod, :, :)
+  end if
   modavg = where(modavg .gt. 1.0e-3, modavg, modavg@_FillValue)
 
   diff = modavg
diff --git a/esmvaltool/diag_scripts/clouds/clouds_interannual.ncl b/esmvaltool/diag_scripts/clouds/clouds_interannual.ncl
index fa81697564..e27984a822 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_interannual.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_interannual.ncl
@@ -31,6 +31,7 @@
 ;   none
 ;
 ; Modification history
+;   20230117-lauer_axel: added support for ICON (code from Manuel)
 ;   20211006-lauer_axel: removed write_plots
 ;   20210413-lauer_axel: added multi-obs mean and average over all models
 ;                        as individual plots
@@ -79,17 +80,18 @@ begin
   end if
 
   variables = metadata_att_as_array(variable_info, "short_name")
-  varidx = ind(variables .eq. var0)
-  if (ismissing(varidx)) then
-    errstr = "diagnostic " + diag + " requires the following variable: var0"
+  if (.not. any(variables .eq. var0)) then
+    errstr = "diagnostic " + diag + " requires the following variable: " + var0
     error_msg("f", DIAG_SCRIPT, "", errstr)
   end if
 
+  var0_info = select_metadata_by_name(variable_info, var0)
+  var0_info := var0_info[0]
   info0 = select_metadata_by_name(input_file_info, var0)
   dim_MOD = ListCount(info0)
 
-  if (isatt(variable_info[varidx], "reference_dataset")) then
-    refname = variable_info[varidx]@reference_dataset
+  if (isatt(var0_info, "reference_dataset")) then
+    refname = var0_info@reference_dataset
   end if
 
   names = metadata_att_as_array(info0, "dataset")
@@ -141,9 +143,11 @@ begin
 
   idxmod = get_mod(names, projects)
 
-  if (idxmod(0) .eq. -1) then
+  if (idxmod(0) .eq. -1) then  ; no model found
     flag_multimod = False
-  else
+  elseif (dimsizes(idxmod) .eq. 1) then  ; one model found
+    flag_multimod = False
+  else  ; more than one model found
     flag_multimod = True
   end if
 
@@ -416,8 +420,8 @@ begin
       data1@diag_script = (/DIAG_SCRIPT/)
     end if
     data1@var = var0  ; Overwrite existing entry
-    if (isatt(variable_info[varidx], "long_name")) then
-      data1@long_name = variable_info[varidx]@long_name
+    if (isatt(var0_info, "long_name")) then
+      data1@long_name = var0_info@long_name
     else
       data1@long_name = var0
     end if
diff --git a/esmvaltool/diag_scripts/clouds/clouds_pdf.ncl b/esmvaltool/diag_scripts/clouds/clouds_pdf.ncl
index d8e24fc995..36f9f6fbef 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_pdf.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_pdf.ncl
@@ -9,7 +9,7 @@
 ;
 ; Required diag_script_info attributes (diagnostic specific)
 ;   xmin                min value for bins (x axis)
-;   xmax                max value for bins (y axis)
+;   xmax                max value for bins (x axis)
 ;
 ; Optional diag_script_info attributes (diagnostic specific)
 ;   filename_add:       optionally add this string to output filenames
@@ -513,36 +513,39 @@ begin
   end if
   idxmod = get_mod(names, projects)
   ; calculate mmm
-  mmmpdf = dim_avg_n_Wrap(data_all(idxmod, :), 0)
-
-  mmmpdf@diag_script = (/DIAG_SCRIPT/)
-  mmmpdf@var = var0 + "_pdf_mmm"
-  mmmpdf@datasets = str_join(names(idxmod), ", ")
-
-  if (var0 .eq. "clt") then
-    x = mmmpdf&$var0$
-    idx60 = ind(x .ge. 60.)
-    idx80 = ind(x .ge. 80.)
-    mmmpdf@cummulative_freq_ge_60 = sum(mmmpdf(idx60))
-    mmmpdf@cummulative_freq_ge_80 = sum(mmmpdf(idx80))
-  end if
+  if (dimsizes(idxmod) .gt. 1) then
+    mmmpdf = dim_avg_n_Wrap(data_all(idxmod, :), 0)
 
-  nc_filename2 = work_dir + "clouds_pdf_" + var0 + filename_add + "_mmm_ref.nc"
+    mmmpdf@diag_script = (/DIAG_SCRIPT/)
+    mmmpdf@var = var0 + "_pdf_mmm"
+    mmmpdf@datasets = str_join(names(idxmod), ", ")
 
-  if (ref_ind .ge. 0) then
-    refpdf = data_all(ref_ind, :)
-    refpdf@diag_script = (/DIAG_SCRIPT/)
-    refpdf@var = var0 + "_pdf_ref"
-    refpdf@datasets = names(ref_ind)
     if (var0 .eq. "clt") then
-      refpdf@cummulative_freq_ge_60 = sum(refpdf(idx60))
-      refpdf@cummulative_freq_ge_80 = sum(refpdf(idx80))
+      x = mmmpdf&$var0$
+      idx60 = ind(x .ge. 60.)
+      idx80 = ind(x .ge. 80.)
+      mmmpdf@cummulative_freq_ge_60 = sum(mmmpdf(idx60))
+      mmmpdf@cummulative_freq_ge_80 = sum(mmmpdf(idx80))
     end if
-    nc_outfile2 = ncdf_write(refpdf, nc_filename2)
-    nc_filename2@existing = "append"
-  end if
 
-  nc_outfile2 = ncdf_write(mmmpdf, nc_filename2)
+    nc_filename2 = work_dir + "clouds_pdf_" + var0 + filename_add + \
+      "_mmm_ref.nc"
+
+      if (ref_ind .ge. 0) then
+      refpdf = data_all(ref_ind, :)
+      refpdf@diag_script = (/DIAG_SCRIPT/)
+      refpdf@var = var0 + "_pdf_ref"
+      refpdf@datasets = names(ref_ind)
+      if (var0 .eq. "clt") then
+        refpdf@cummulative_freq_ge_60 = sum(refpdf(idx60))
+        refpdf@cummulative_freq_ge_80 = sum(refpdf(idx80))
+      end if
+      nc_outfile2 = ncdf_write(refpdf, nc_filename2)
+      nc_filename2@existing = "append"
+    end if
+
+    nc_outfile2 = ncdf_write(mmmpdf, nc_filename2)
+  end if
 
   ; ========================================================================
 
diff --git a/esmvaltool/diag_scripts/clouds/clouds_scatter.ncl b/esmvaltool/diag_scripts/clouds/clouds_scatter.ncl
index b019718dfb..35e217613b 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_scatter.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_scatter.ncl
@@ -28,6 +28,7 @@
 ;   none
 ;
 ; Modification history
+;   20230117-lauer_axel: added support for ICON (code from Manuel)
 ;   20210210-lauer_axel: written.
 ;
 ; ############################################################################
@@ -44,7 +45,8 @@ begin
   enter_msg(DIAG_SCRIPT, "")
 
   diag = "clouds_scatter.ncl"
-  variables = metadata_att_as_array(variable_info, "short_name")
+  variables = get_unique_values(metadata_att_as_array(variable_info, \
+                                                      "short_name"))
 
   ; Check required diag_script_info attributes
   exit_if_missing_atts(diag_script_info, (/"var_x", "var_y", "xmin", "xmax"/))
@@ -98,21 +100,26 @@ begin
   ; is equal for each variable
 
   do i = 0, nVAR - 1
-    if (isatt(variable_info[idx(i)], "reference_dataset")) then
-      refname(i) = variable_info[idx(i)]@reference_dataset
+    var = variables(idx(i))
+    var_info = select_metadata_by_name(variable_info, var)
+    var_info := var_info[0]
+    if (isatt(var_info, "reference_dataset")) then
+      refname(i) = var_info@reference_dataset
     end if
-    info = select_metadata_by_name(input_file_info, variables(idx(i)))
+    info = select_metadata_by_name(input_file_info, var)
     if (i .eq. 0) then
       dim_MOD = ListCount(info)
     else
       dim_test = ListCount(info)
       if (dim_test .ne. dim_MOD) then
         error_msg("f", DIAG_SCRIPT, "", "number of datasets for variable " \
-                  + variables(i) + " does not match number of datasets for " \
-                  + variables(0))
+                  + var + " does not match number of datasets for " \
+                  + variables(idx(0)))
       end if
     end if
     delete(info)
+    delete(var)
+    delete(var_info)
   end do
 
   ; Set default values for non-required diag_script_info attributes
@@ -347,6 +354,9 @@ begin
   if (idxmod(0) .eq. -1) then
     flag_multimod = False
     mm_ind = -1
+  elseif (dimsizes(idxmod) .eq. 1) then
+    flag_multimod = False
+    mm_ind = -1
   else
     flag_multimod = True
     mmavg = new((/1, nbins/), float)
diff --git a/esmvaltool/diag_scripts/clouds/clouds_seasonal_cycle.ncl b/esmvaltool/diag_scripts/clouds/clouds_seasonal_cycle.ncl
index f7f91e5122..d918e9eccb 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_seasonal_cycle.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_seasonal_cycle.ncl
@@ -29,6 +29,7 @@
 ;   none
 ;
 ; Modification history
+;   20230117-lauer_axel: added support for ICON (code from Manuel)
 ;   20210415-lauer_axel: written.
 ;
 ; ############################################################################
@@ -64,17 +65,18 @@ begin
   end if
 
   variables = metadata_att_as_array(variable_info, "short_name")
-  varidx = ind(variables .eq. var0)
-  if (ismissing(varidx)) then
-    errstr = "diagnostic " + diag + " requires the following variable: var0"
+  if (.not. any(variables .eq. var0)) then
+    errstr = "diagnostic " + diag + " requires the following variable: " + var0
     error_msg("f", DIAG_SCRIPT, "", errstr)
   end if
 
+  var0_info = select_metadata_by_name(variable_info, var0)
+  var0_info := var0_info[0]
   info0 = select_metadata_by_name(input_file_info, var0)
   dim_MOD = ListCount(info0)
 
-  if (isatt(variable_info[varidx], "reference_dataset")) then
-    refname = variable_info[varidx]@reference_dataset
+  if (isatt(var0_info, "reference_dataset")) then
+    refname = var0_info@reference_dataset
   end if
 
   names = metadata_att_as_array(info0, "dataset")
@@ -119,9 +121,11 @@ begin
 
   idxmod = get_mod(names, projects)
 
-  if (idxmod(0) .eq. -1) then
+  if (idxmod(0) .eq. -1) then  ; no model found
     flag_multimod = False
-  else
+  elseif (dimsizes(idxmod) .eq. 1) then  ; one model found
+    flag_multimod = False
+  else  ; more than one model found
     flag_multimod = True
   end if
 
@@ -375,8 +379,8 @@ begin
       data1@diag_script = (/DIAG_SCRIPT/)
     end if
     data1@var = var0  ; Overwrite existing entry
-    if (isatt(variable_info[varidx], "long_name")) then
-      data1@long_name = variable_info[varidx]@long_name
+    if (isatt(var0_info, "long_name")) then
+      data1@long_name = var0_info@long_name
     else
       data1@long_name = var0
     end if
diff --git a/esmvaltool/diag_scripts/clouds/clouds_taylor.ncl b/esmvaltool/diag_scripts/clouds/clouds_taylor.ncl
index 7a4763e852..e2de9117cd 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_taylor.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_taylor.ncl
@@ -79,6 +79,7 @@
 ;        variable has been calculated from) the *second* (third, ...)
 ;
 ; Modification history
+;   20230117-lauer_axel: added support for ICON (code from Manuel)
 ;   20211028-lauer_axel: added option to calculate tcwp as lwp + iwp
 ;   20211018-lauer_axel: added option to remove individual models from legend
 ;   20211006-lauer_axel: removed write_plots
@@ -128,7 +129,8 @@ begin
 
   estimate_obs_uncertainty = diag_script_info@estimate_obs_uncertainty
 
-  variables = metadata_att_as_array(variable_info, "short_name")
+  variables = get_unique_values(metadata_att_as_array(variable_info, \
+                                                      "short_name"))
   numvars = dimsizes(variables)
 
   ; find "main" variable and if present uncertainty estimates or auxiliary
@@ -159,17 +161,20 @@ begin
         end if
       end do
     end do
-    var0 = variable_info[mainvarind]@short_name
+    var0 = variables(mainvarind)
   else
     mainvarind = ind(variables .eq. diag_script_info@var)
     var0 = diag_script_info@var
   end if
 
   if (ismissing(mainvarind)) then
-    errstr = "diagnostic " + diag + " requires the following variable: var0"
+    errstr = "diagnostic " + diag + " requires the following variable: " + var0
     error_msg("f", DIAG_SCRIPT, "", errstr)
   end if
 
+  var0_info = select_metadata_by_name(variable_info, var0)
+  var0_info := var0_info[0]
+
   flag_multiobs_unc = diag_script_info@multiobs_uncertainty
   multiobs_exclude = diag_script_info@multiobs_exclude
   if (estimate_obs_uncertainty .and. flag_multiobs_unc) then
@@ -197,11 +202,11 @@ begin
 
   info0 = select_metadata_by_name(input_file_info, var0)
   dim_MOD = ListCount(info0)
-  if (isatt(variable_info[mainvarind], "reference_dataset")) then
-    refname = variable_info[mainvarind]@reference_dataset
+  if (isatt(var0_info, "reference_dataset")) then
+    refname = var0_info@reference_dataset
   end if
-  if (isatt(variable_info[mainvarind], "alternative_dataset")) then
-    refname2 = variable_info[mainvarind]@alternative_dataset
+  if (isatt(var0_info, "alternative_dataset")) then
+    refname2 = var0_info@alternative_dataset
   end if
   names = metadata_att_as_array(info0, "dataset")
   projects = metadata_att_as_array(info0, "project")
diff --git a/esmvaltool/diag_scripts/clouds/clouds_zonal.ncl b/esmvaltool/diag_scripts/clouds/clouds_zonal.ncl
index a77677771b..32b7befb7f 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_zonal.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_zonal.ncl
@@ -20,7 +20,9 @@
 ;   panel_labels:       label individual panels (true, false)
 ;   PanelTop:           manual override for "@gnsPanelTop" used by panel
 ;                       plot(s)
-;   showdiff            calculate and plot differences (default = False)
+;   showdiff:           calculate and plot differences (default = False)
+;   showyears:          add start and end years to the plot titles
+;                       (default = False)
 ;   rel_diff:           if showdiff = True, then plot relative differences (%)
 ;                       (default = False)
 ;   rel_diff_min:       lower cutoff value in case of calculating relative
@@ -45,6 +47,7 @@
 ;   none
 ;
 ; Modification history
+;   20230117-lauer_axel: added support for ICON (code from Manuel)
 ;   20200211-lauer_axel: written.
 ;
 ; ############################################################################
@@ -82,6 +85,7 @@ begin
   set_default_att(diag_script_info, "rel_diff", False)
   set_default_att(diag_script_info, "rel_diff_min", -1.0e19)
   set_default_att(diag_script_info, "showdiff", False)
+  set_default_att(diag_script_info, "showyears", False)
   set_default_att(diag_script_info, "t_test", False)
   set_default_att(diag_script_info, "timemean", "annualclim")
   set_default_att(diag_script_info, "units_to", "")
@@ -327,8 +331,18 @@ begin
 
     ; annotation
 
+    ; if desired, add years to plot title
+    years_str = ""
+    if (diag_script_info@showyears) then
+      years_str = " (" + variable_info[0]@start_year
+      if (variable_info[0]@start_year .ne. variable_info[0]@end_year) then
+        years_str = years_str + "-" + variable_info[0]@end_year
+      end if
+      years_str = years_str + ")"
+    end if
+
 ;    res@tiMainOn             = False
-    res@tiMainString = names(ii)
+    res@tiMainString = names(imod) + years_str
     res@cnLevelSelectionMode = "ExplicitLevels"
     res@cnLinesOn            = False
 
@@ -632,7 +646,7 @@ begin
                                             / tofloat((dim_MOD + 1) / 2)/))
 
 ;      dres@tiMainOn       = False
-      dres@tiMainString = names(ii) + " - " + refname
+      dres@tiMainString = names(imod) + " - " + refname + years_str
 
       dres@cnFillOn       = True      ; color plot desired
       dres@cnLineLabelsOn = False     ; contour lines

From 4b05cc93dd84244306358210d95cc853b102aa39 Mon Sep 17 00:00:00 2001
From: Axel Lauer <axel.lauer@dlr.de>
Date: Wed, 18 Jan 2023 10:31:06 +0100
Subject: [PATCH 2/5] increased flexibility of (old) cloud diagnostics

---
 doc/sphinx/source/recipes/recipe_clouds.rst   |  1 +
 .../diag_scripts/clouds/clouds_bias.ncl       | 31 ++++++++----
 .../diag_scripts/clouds/clouds_ipcc.ncl       | 50 +++++++++++++------
 .../clouds/clouds_lifrac_scatter.ncl          | 49 +++++++++++++++---
 4 files changed, 97 insertions(+), 34 deletions(-)

diff --git a/doc/sphinx/source/recipes/recipe_clouds.rst b/doc/sphinx/source/recipes/recipe_clouds.rst
index 9328c1af3c..3106a616d7 100644
--- a/doc/sphinx/source/recipes/recipe_clouds.rst
+++ b/doc/sphinx/source/recipes/recipe_clouds.rst
@@ -235,6 +235,7 @@ User settings in recipe
    *Optional settings (scripts)*
 
    * explicit_cn_levels: contour levels
+   * highlight_dataset: name of dataset to highlight (default = "MultiModelMean")
    * mask_ts_sea_ice: true = mask T < 272 K as sea ice (only for variable "ts");
      false = no additional grid cells masked for variable "ts"
    * projection: map projection, e.g., Mollweide, Mercator
diff --git a/esmvaltool/diag_scripts/clouds/clouds_bias.ncl b/esmvaltool/diag_scripts/clouds/clouds_bias.ncl
index db8435aaea..af4b04698a 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_bias.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_bias.ncl
@@ -4,7 +4,7 @@
 ; PROJECT-NAME EMBRACE
 ; ############################################################################
 ; Description
-;   Calculates the multi-model mean bias, absolute difference and relative
+;   Calculates the (multi-model mean) bias, absolute difference and relative
 ;   difference of annual mean 2-d cloud variables compared with a
 ;   reference dataset (observations).
 ;
@@ -28,6 +28,7 @@
 ;   none
 ;
 ; Modification history
+;   20230118-lauer_axel: added support to plot just 1 model
 ;   20211006-lauer_axel: removed write_plots
 ;   20190222-lauer_axel: added output of provenance (v2.0)
 ;   20181119-lauer_axel: adapted code to multi-variable capable framework
@@ -126,7 +127,15 @@ begin
   mm_ind = ind(names .eq. "MultiModelMean")
 
   if (ismissing(mm_ind)) then
-    error_msg("f", DIAG_SCRIPT, "", "multi-model mean is missing (required)")
+    mod_ind = ind(names .ne. ref_ind)
+    if (all(ismissing(mod_ind))) then
+      error_msg("f", DIAG_SCRIPT, "", "no dataset besides reference " \
+        + "dataset found. Cannot continue.")
+    end if
+    mm_ind = mod_ind(0)
+    log_info("multi-model mean is missing, using first dataset (" \
+      + names(mm_ind) + ") instead")
+    delete(mod_ind)
   end if
 
   ; basename of diag_script
@@ -233,7 +242,7 @@ begin
 
   diff@res_tiMainFontHeightF  = 0.016
 
-  diff@res_tiMainString = "Multi Model Mean Bias"
+  diff@res_tiMainString = names(mm_ind) + " Bias"
 
   copy_VarMeta(diff, mmdata)
   delete(mmdata@res_cnLevels)
@@ -272,7 +281,7 @@ begin
     end if
   end if
 
-  mmdata@res_tiMainString = "Multi Model Mean"
+  mmdata@res_tiMainString = names(mm_ind)
 
   plotsperline = (/2, 0/)
   plotind = (/0, 1/)  ; mmm and mean bias are always plotted
@@ -298,7 +307,7 @@ begin
     absdiff@res_cnLevels = tmp(1:dimsizes(tmp)-1)
     delete(tmp)
 
-    absdiff@res_tiMainString = "Multi Model Mean of Absolute Error"
+    absdiff@res_tiMainString = names(mm_ind) + " Absolute Error"
 
     iadd = 2
     itmp = array_append_record(plotind, iadd, 0)
@@ -319,7 +328,7 @@ begin
     copy_VarMeta(diff, reldiff)
     delete(reldiff@res_cnLevels)
     reldiff@res_cnLevels = fspan(-90.0, 90.0, 13)
-    reldiff@res_tiMainString = "Multi Model Mean of Relative Error"
+    reldiff@res_tiMainString = names(mm_ind) + " Relative Error"
     reldiff@units = "%"
     reldiff@res_lbTitleString = "(" + reldiff@units + ")"
     if (isvar("pal4")) then
@@ -384,7 +393,7 @@ begin
 
     ; add meta data to plot (for reporting)
 
-    caption = "Multi model values, from top left to bottom right: " \
+    caption = names(mm_ind) + " values, from top left to bottom right: " \
               + "mean, bias"
     if (plot_abs_diff) then
       caption = caption + ", absolute error"
@@ -403,22 +412,22 @@ begin
   nc_filename@existing = "append"
 
   mmdata@var = var0 + "_mean"
-  mmdata@long_name = var0 + " (multi-model mean)"
+  mmdata@long_name = var0 + " " + names(mm_ind) + " (mean)"
   nc_outfile = ncdf_write(mmdata, nc_filename)
 
   diff@var = var0 + "_bias"
-  diff@long_name = var0 + " (multi-model bias)"
+  diff@long_name = var0 + " " + names(mm_ind) + " (bias)"
   nc_outfile = ncdf_write(diff, nc_filename)
 
   if (isvar("absdiff")) then
     absdiff@var = var0 + "_abs_bias"
-    absdiff@long_name = var0 + " (multi-model absolute bias)"
+    absdiff@long_name = var0 + " " + names(mm_ind) + " (absolute bias)"
     nc_outfile = ncdf_write(absdiff, nc_filename)
   end if
 
   if (isvar("reldiff")) then
     reldiff@var = var0 + "_rel_bias"
-    reldiff@long_name = var0 + " (multi-model relative bias)"
+    reldiff@long_name = var0 + " " + names(mm_ind) + " (relative bias)"
     reldiff@units = reldiff@units
     nc_outfile = ncdf_write(reldiff, nc_filename)
   end if
diff --git a/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl b/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
index 5a948a3fe7..0c20ef8290 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
@@ -17,6 +17,7 @@
 ;
 ; Optional diag_script_info attributes (diagnostic specific)
 ;   explicit_cn_levels: contour levels
+;   highlight_dataset:  name of dataset to highlight ("MultiModelMean")
 ;   mask_ts_sea_ice:    - True = mask T < 272 K as sea ice (only for
 ;                         variable "ts")
 ;                       - False = no additional grid cells masked for
@@ -50,6 +51,8 @@
 ;        plot is written to 2 separate netCDFs.
 ;
 ; Modification history
+;   20230118-lauer_axel: added support to highlight other dataset than
+;                        MultiModelMean
 ;   20211006-lauer_axel: removed write_plots
 ;   20190222-lauer_axel: added output of provenance (v2.0)
 ;   20181119-lauer_axel: adapted code to multi-variable capable framework
@@ -127,6 +130,7 @@ begin
   log_info("++++++++++++++++++++++++++++++++++++++++++")
 
   ; Set default values for non-required diag_script_info attributes
+  set_default_att(diag_script_info, "highlight_dataset", "MultiModelMean")
   set_default_att(diag_script_info, "mask_ts_sea_ice", False)
   set_default_att(diag_script_info, "projection", "CylindricalEquidistant")
   set_default_att(diag_script_info, "timemean", "annualclim")
@@ -192,12 +196,22 @@ begin
 
   ; get multi-model mean index
 
-  mm_ind = ind(names .eq. "MultiModelMean")
+  mm_ind = ind(names .eq. diag_script_info@highlight_dataset)
 
   if (ismissing(mm_ind)) then
-    error_msg("f", DIAG_SCRIPT, "", "multi-model mean is missing (required)")
+    mod_ind = ind(names .ne. ref_ind)
+    if (all(ismissing(mod_ind))) then
+      error_msg("f", DIAG_SCRIPT, "", "no dataset besides reference " \
+        + "dataset found. Cannot continue.")
+    end if
+    mm_ind = mod_ind(0)
+    log_info("highlight_dataset (" + diag_script_info@highlight_dataset \
+      + " ) not found, using first dataset (" + names(mm_ind) + ") instead.")
+    delete(mod_ind)
   end if
 
+  highlight_name = names(mm_ind)
+
   mask_ts_sea_ice = diag_script_info@mask_ts_sea_ice
 
   if ((var0 .eq. "ts") .and. (mask_ts_sea_ice))
@@ -357,7 +371,8 @@ begin
 
     delete(tmp)
 
-    ; save maps of multi-model mean and reference data
+    ; save maps of highlight_dataset (default = multi-model mean)
+    ; and reference data
 
     if (imod.eq.mm_ind) then
       mmdata = data
@@ -376,7 +391,8 @@ begin
 
   end do  ; imod
 
-  ; differences between multi-model mean and reference data set
+  ; differences between highlight_dataset (default = multi-model mean)
+  ; and reference data set
 
   diff = mmdata - refdata
   copy_VarMeta(refdata, diff)
@@ -384,17 +400,19 @@ begin
 ;  debugfile->diff = diff
 
   ; we order the zonal mean array in a way so that
-  ; the lines for the multi-model mean and reference model will
-  ; be drawn on top of the lines for the individual models, i.e.:
+  ; the lines for the highligh_dataset (default = multi-model mean)
+  ; and reference model will be drawn on top of the lines for the individual
+  ; models, i.e.:
   ;  (1) individual model(s)
   ;  (2) reference model(s) (observations)
-  ;  (3) multi-model mean
+  ;  (3) highlight_dataset (default = multi-model mean)
 
   dims = dimsizes(zm)
   zonalmean = new(dims, float)
   copy_VarMeta(zm, zonalmean)
 
-  ; model indices with no reference model(s) and no multi-model mean
+  ; model indices with no reference model(s) and no highlight_dataset
+  ; (default = multi-model mean)
 
   model_ind = ispan(0, dim_MOD0 - 1, 1)
   model_ind(ref_ind) = -1
@@ -441,14 +459,14 @@ begin
     zonalmean&model(n) = zm&model(ref_ind2)
   end if
 
-  ; last entry in "zonalmean" = multi-model mean
+  ; last entry in "zonalmean" = highlight_dataset (default = multi-model mean)
 
   n = n + 1
 
   if (numseas.gt.1) then
-    zonalmean(n, :, :) = zm(mm_ind, :, :)  ; multi-model mean
+    zonalmean(n, :, :) = zm(mm_ind, :, :)  ; highlight_dataset
   else
-    zonalmean(n, :) = zm(mm_ind, :)  ; multi-model mean
+    zonalmean(n, :) = zm(mm_ind, :)  ; highlight_dataset
   end if
 
   zonalmean&model(n) = zm&model(mm_ind)
@@ -621,7 +639,7 @@ begin
     n = dimsizes(modelsonly_ind) - 1
 
     ; settings for all models that have been used to calculate the
-    ; multi-model mean (= first entries in "zonalmean")
+    ; highligh_dataset (= first entries in "zonalmean")
 
     linethickness(0:n)   = 1.0
     linecolor(0:n)       = "(/0.5, 0.5, 0.5/)"
@@ -648,9 +666,9 @@ begin
       end do
     else
       linethickness(n+1:dim_MOD0-1) = 4.0      ; reference dataset(s)
-      linethickness(mm_ind)         = 4.0      ; multi-model mean
+      linethickness(mm_ind)         = 4.0      ; highlight_dataset
       linecolor(n+1:dim_MOD0-1)     = "Black"  ; reference data set
-      linecolor(mm_ind)             = "Red"    ; multi-model mean
+      linecolor(mm_ind)             = "Red"    ; highlight_dataset
     end if
 
     res                   = True
@@ -708,7 +726,7 @@ begin
     pres                = True
     pres@gsnPanelCenter = False
     pres@gsnPanelXF     = (/0.075, 0.625/)  ; hor. pos. of sub-plots
-    pres@txString       = season(is)
+    pres@txString       = highlight_name + " (" + season(is) + ")"
 
     outfile = panelling(wks, plots(:, is), 1, 2, pres)
     log_info("Wrote " + wks@fullname)
@@ -745,7 +763,7 @@ begin
   statistics = "clim"
   domain = "global"
   plottype = (/"geo", "zonal"/)
-  caption = "Multi model mean bias (left) and zonal averages (right) " \
+  caption = names(mm_ind) + " bias (left) and zonal averages (right) " \
             + "for variable " + var0 + " (" + allseas \
             + "), reference = " + names(ref_ind) + "."
 
diff --git a/esmvaltool/diag_scripts/clouds/clouds_lifrac_scatter.ncl b/esmvaltool/diag_scripts/clouds/clouds_lifrac_scatter.ncl
index 5149c01bf2..98d2ca9fe2 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_lifrac_scatter.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_lifrac_scatter.ncl
@@ -27,6 +27,7 @@
 ;   none
 ;
 ; Modification history
+;   20230117-lauer_axel: added support for ICON (code from Manuel)
 ;   20210302-lauer_axel: written.
 ;
 ; ############################################################################
@@ -36,13 +37,15 @@ load "$diag_scripts/../interface_scripts/interface.ncl"
 load "$diag_scripts/shared/plot/aux_plotting.ncl"
 load "$diag_scripts/shared/statistics.ncl"
 load "$diag_scripts/shared/plot/style.ncl"
+load "$diag_scripts/shared/dataset_selection.ncl"
 
 begin
 
   enter_msg(DIAG_SCRIPT, "")
 
   diag = "clouds_lifrac_scatter.ncl"
-  variables = metadata_att_as_array(variable_info, "short_name")
+  variables = get_unique_values(metadata_att_as_array(variable_info, \
+                                                      "short_name"))
   dim_VAR = dimsizes(variables)
   refname = new(dim_VAR, string)
 
@@ -73,10 +76,13 @@ begin
   ; is equal for each variable
 
   do i = 0, dim_VAR - 1
-    if (isatt(variable_info[idx(i)], "reference_dataset")) then
-      refname(i) = variable_info[idx(i)]@reference_dataset
+    var = variables(idx(i))
+    var_info = select_metadata_by_name(variable_info, var)
+    var_info := var_info[0]
+    if (isatt(var_info, "reference_dataset")) then
+      refname(i) = var_info@reference_dataset
     end if
-    info = select_metadata_by_name(input_file_info, variables(idx(i)))
+    info = select_metadata_by_name(input_file_info, var)
     if (i .eq. 0) then
       dim_MOD = ListCount(info)
     else
@@ -88,6 +94,8 @@ begin
       end if
     end if
     delete(info)
+    delete(var)
+    delete(var_info)
   end do
 
   delete(idx)
@@ -163,10 +171,24 @@ begin
               + "are not available: " + str)
   end if
 
+  add_dim_MOD = 0
+
   if (mm_mean_median) then
-    add_dim_MOD = 2
-  else
-    add_dim_MOD = 0
+    ; check if enough model datasets are available to calculate multi-model
+    ; mean and median ice/liquid fractions
+
+    ; find all indices of models w/o MultiModelMean/MultiModelMedian
+    ; (if present)
+    idxmod = get_mod(names_clw, projects_clw)
+  
+    if ((idxmod(0) .eq. -1) .or. (dimsizes(idxmod) .le. 2)) then
+      log_info("Not enough model datasets to calculate multi-model " \
+        + "mean and median. Setting of 'mm_mean_median' will be ignored.")
+      mm_mean_median = False
+    else  ; more than one model found
+      add_dim_MOD = 2
+    end if
+    delete(idxmod)
   end if
 
   resulti = new((/dim_MOD + add_dim_MOD, nbins/), float)
@@ -474,6 +496,11 @@ begin
           mmed25i(n) = stat(6)                ; lower quartile
           resulti(index_my_med, n) = stat(8)  ; median
           mmed75i(n) = stat(10)               ; upper quartile
+          ; if number of datasets < 4, quartiles cannot be calculated
+          if (ismissing(mmed25i(n)) .or. ismissing(mmed75i(n))) then
+            mmed25i(n) = 0.0
+            mmed75i(n) = 0.0
+          end if
         else
           resulti(index_my_mm, n) = 0.0
           mmstdi(n) = 0.0
@@ -489,6 +516,11 @@ begin
           mmed25l(n) = stat(6)                ; lower quartile
           resultl(index_my_med, n) = stat(8)  ; median
           mmed75l(n) = stat(10)               ; upper quartile
+          ; if number of datasets < 4, quartiles cannot be calculated
+          if (ismissing(mmed25l(n)) .or. ismissing(mmed75l(n))) then
+            mmed25l(n) = 0.0
+            mmed75l(n) = 0.0
+          end if
         else
           resultl(index_my_mm, n) = 0.0
           mmstdl(n) = 0.0
@@ -504,6 +536,9 @@ begin
     end if
     delete(idx0)
     delete(imod)
+  else
+    index_my_mm = -1
+    index_my_med = -1
   end if
 
 ;  if (mm_ind .ge. 0) then

From 1e1b8f5fe441880999c26179d501c255068f26de Mon Sep 17 00:00:00 2001
From: Axel Lauer <axel.lauer@dlr.de>
Date: Wed, 18 Jan 2023 10:52:48 +0100
Subject: [PATCH 3/5] fixed style issues

---
 esmvaltool/diag_scripts/clouds/clouds_bias.ncl           | 6 +++---
 esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl           | 9 +++++----
 esmvaltool/diag_scripts/clouds/clouds_lifrac_scatter.ncl | 5 +++--
 3 files changed, 11 insertions(+), 9 deletions(-)

diff --git a/esmvaltool/diag_scripts/clouds/clouds_bias.ncl b/esmvaltool/diag_scripts/clouds/clouds_bias.ncl
index af4b04698a..5eeaaaded9 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_bias.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_bias.ncl
@@ -130,11 +130,11 @@ begin
     mod_ind = ind(names .ne. ref_ind)
     if (all(ismissing(mod_ind))) then
       error_msg("f", DIAG_SCRIPT, "", "no dataset besides reference " \
-        + "dataset found. Cannot continue.")
+                + "dataset found. Cannot continue.")
     end if
     mm_ind = mod_ind(0)
     log_info("multi-model mean is missing, using first dataset (" \
-      + names(mm_ind) + ") instead")
+             + names(mm_ind) + ") instead")
     delete(mod_ind)
   end if
 
@@ -394,7 +394,7 @@ begin
     ; add meta data to plot (for reporting)
 
     caption = names(mm_ind) + " values, from top left to bottom right: " \
-              + "mean, bias"
+      + "mean, bias"
     if (plot_abs_diff) then
       caption = caption + ", absolute error"
     end if
diff --git a/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl b/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
index 0c20ef8290..429960c463 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
@@ -202,11 +202,12 @@ begin
     mod_ind = ind(names .ne. ref_ind)
     if (all(ismissing(mod_ind))) then
       error_msg("f", DIAG_SCRIPT, "", "no dataset besides reference " \
-        + "dataset found. Cannot continue.")
+                + "dataset found. Cannot continue.")
     end if
     mm_ind = mod_ind(0)
     log_info("highlight_dataset (" + diag_script_info@highlight_dataset \
-      + " ) not found, using first dataset (" + names(mm_ind) + ") instead.")
+             + " ) not found, using first dataset (" + names(mm_ind) \
+             + \") instead.")
     delete(mod_ind)
   end if
 
@@ -764,8 +765,8 @@ begin
   domain = "global"
   plottype = (/"geo", "zonal"/)
   caption = names(mm_ind) + " bias (left) and zonal averages (right) " \
-            + "for variable " + var0 + " (" + allseas \
-            + "), reference = " + names(ref_ind) + "."
+    + "for variable " + var0 + " (" + allseas \
+    + "), reference = " + names(ref_ind) + "."
 
   do is = 0, numseas - 1
     log_provenance(nc_outfile_bias, plotfile(is), caption, statistics, \
diff --git a/esmvaltool/diag_scripts/clouds/clouds_lifrac_scatter.ncl b/esmvaltool/diag_scripts/clouds/clouds_lifrac_scatter.ncl
index 98d2ca9fe2..f9bf33f30d 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_lifrac_scatter.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_lifrac_scatter.ncl
@@ -180,10 +180,11 @@ begin
     ; find all indices of models w/o MultiModelMean/MultiModelMedian
     ; (if present)
     idxmod = get_mod(names_clw, projects_clw)
-  
+
     if ((idxmod(0) .eq. -1) .or. (dimsizes(idxmod) .le. 2)) then
       log_info("Not enough model datasets to calculate multi-model " \
-        + "mean and median. Setting of 'mm_mean_median' will be ignored.")
+               + "mean and median. Setting of 'mm_mean_median' will" \
+               + " be ignored.")
       mm_mean_median = False
     else  ; more than one model found
       add_dim_MOD = 2

From 88af78363dcb1b62841fccb9ff2e1a29f93d6693 Mon Sep 17 00:00:00 2001
From: Axel Lauer <axel.lauer@dlr.de>
Date: Fri, 20 Jan 2023 11:20:58 +0100
Subject: [PATCH 4/5] fixed syntax error in clouds_ipcc.ncl

---
 esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl b/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
index 429960c463..6f18510c58 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
@@ -207,7 +207,7 @@ begin
     mm_ind = mod_ind(0)
     log_info("highlight_dataset (" + diag_script_info@highlight_dataset \
              + " ) not found, using first dataset (" + names(mm_ind) \
-             + \") instead.")
+             + ") instead.")
     delete(mod_ind)
   end if
 

From 94f3dbc7475e8f352f2d6679ca3234f6f06daca3 Mon Sep 17 00:00:00 2001
From: Axel Lauer <axel.lauer@dlr.de>
Date: Fri, 20 Jan 2023 14:30:21 +0100
Subject: [PATCH 5/5] fixed single dataset bug in clouds_ipcc.ncl

---
 esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl | 18 +++++++++++-------
 1 file changed, 11 insertions(+), 7 deletions(-)

diff --git a/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl b/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
index 6f18510c58..101ca0a080 100644
--- a/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
+++ b/esmvaltool/diag_scripts/clouds/clouds_ipcc.ncl
@@ -425,18 +425,22 @@ begin
   modelsonly_ind = ind(model_ind.ge.0)
   delete(model_ind)
 
-  n = dimsizes(modelsonly_ind) - 1
+  if (.not.all(ismissing(modelsonly_ind))) then
+    n = dimsizes(modelsonly_ind) - 1
 
-  ; first entries in "zonalmean" = individual models
+    ; first entries in "zonalmean" = individual models
 
-  if (numseas.gt.1) then
-    zonalmean(0:n, :, :) = zm(modelsonly_ind, :, :)
+    if (numseas.gt.1) then
+      zonalmean(0:n, :, :) = zm(modelsonly_ind, :, :)
+    else
+      zonalmean(0:n, :) = zm(modelsonly_ind, :)
+    end if
+
+    zonalmean&model(0:n) = zm&model(modelsonly_ind)
   else
-    zonalmean(0:n, :) = zm(modelsonly_ind, :)
+    n = -1
   end if
 
-  zonalmean&model(0:n) = zm&model(modelsonly_ind)
-
   ; observation(s)
 
   n = n + 1