diff --git a/MY_RUN/KTEST/003_KW78/004_python/plot_003_KW78.py b/MY_RUN/KTEST/003_KW78/004_python/plot_003_KW78.py
index 9ec2501f9495c586925e8dc649351a883f6db275..139b4a13f35d2d4fdd92d5e94b46716cf7f170ca 100644
--- a/MY_RUN/KTEST/003_KW78/004_python/plot_003_KW78.py
+++ b/MY_RUN/KTEST/003_KW78/004_python/plot_003_KW78.py
@@ -69,7 +69,7 @@ Ltitle = ['Wind at K=2', 'Wind at 3000m', 'Wind at 5000m']
 Lxlab = ['x (m)']*len(Lplot)
 Lylab = ['y (m)']*len(Lplot)
 Llegendval = [10,10,10]
-Lcbarlabel = ['m/s']*len(Lplot)
+Llegendlabel = ['m/s']*len(Lplot)
 Larrowstep = [1]*len(Lplot)
 Lwidth = [0.002]*len(Lplot)
 Lcolor = ['black']*len(Lplot)
@@ -79,7 +79,7 @@ lat = [Dvar['f1']['nj_u'], Dvar['f2']['nj'],  Dvar['f2']['nj'] ]
 Lscale = [200]*len(Lplot)
 fig2 = Panel1.pvector(Lxx=lon, Lyy=lat, Lvar1=Lplot1, Lvar2=Lplot2, Lcarte=[500,23500,500,23500], Llevel=Llvl, 
                       Lxlab=Lxlab, Lylab=Lylab, Ltitle=Ltitle, Lwidth=Lwidth, Larrowstep=Larrowstep, 
-                      Lcolor=Lcolor, Llegendval=Llegendval, Lcbarlabel=Lcbarlabel, Lid_overlap=[4,6,8], ax=fig1.axes, Lscale=Lscale)
+                      Lcolor=Lcolor, Llegendval=Llegendval, Llegendlabel=Llegendlabel, Lid_overlap=[4,6,8], ax=fig1.axes, Lscale=Lscale)
 #  Oblique projection
 i_beg, j_beg = (3,0)
 i_end, j_end = (22,21)
@@ -133,13 +133,13 @@ Lplot1 = [ WIND_sec1]
 Lplot2 = [ WT_sec1]
 Ltitle = ['Wind']
 Llegendval = [25]
-Lcbarlabel = ['m/s']*len(Lplot)
+Llegendlabel = ['m/s']*len(Lplot)
 Larrowstep = [1]*len(Lplot)
 Lwidth = [0.004]*len(Lplot)
 Lscale = [200]*len(Lplot)
 
 fig4 = Panel2.pvector(Lxx=LaxeX, Lyy=LaxeZ, Lvar1=Lplot1, Lvar2=Lplot2, Lxlab=Lxlab, Lylab=Lylab, Ltitle=Ltitle, Lwidth=Lwidth, Larrowstep=Larrowstep, 
-                        Llegendval=Llegendval, Lcbarlabel=Lcbarlabel, Lid_overlap=[0], ax=fig3.axes, Lscale=Lscale)
+                        Llegendval=Llegendval, Llegendlabel=Llegendlabel, Lid_overlap=[0], ax=fig3.axes, Lscale=Lscale)
 
 Lplot = [RRT_sec1]
 LaxeX = [axe_m1]
@@ -195,13 +195,13 @@ Lplot1 = [ Dvar['f1']['VM'][:,:,13]]
 Lplot2 = [ Dvar['f1']['WM'][:,:,13]]
 Ltitle = ['Wind']
 Llegendval = [25]
-Lcbarlabel = ['m/s']*len(Lplot)
+Llegendlabel = ['m/s']*len(Lplot)
 Larrowstep = [1]*len(Lplot)
 Lwidth = [0.004]*len(Lplot)
 Lscale = [200]*len(Lplot)
 
 fig7 = Panel3.pvector(Lxx=LaxeX, Lyy=LaxeZ, Lvar1=Lplot1, Lvar2=Lplot2, Lxlab=Lxlab, Lylab=Lylab, Ltitle=Ltitle, Lwidth=Lwidth, Larrowstep=Larrowstep, 
-                        Llegendval=Llegendval, Lcbarlabel=Lcbarlabel, Lid_overlap=[0], ax=fig6.axes, Lscale=Lscale)
+                        Llegendval=Llegendval, Llegendlabel=Llegendlabel, Lid_overlap=[0], ax=fig6.axes, Lscale=Lscale)
 
 
 Lplot = [Dvar['f1']['RRT'][:,:,13]]
diff --git a/MY_RUN/KTEST/004_Reunion/007_python/plot_004_Reunion.py b/MY_RUN/KTEST/004_Reunion/007_python/plot_004_Reunion.py
index ea59e343330e5dfc43a383bb976da141751a6705..5de8f4b897fab5ecd28427472f553f42100b6f99 100644
--- a/MY_RUN/KTEST/004_Reunion/007_python/plot_004_Reunion.py
+++ b/MY_RUN/KTEST/004_Reunion/007_python/plot_004_Reunion.py
@@ -23,7 +23,8 @@ os.system('rm -f tempgraph*')
 LnameFiles = ['REUNI.1.00A20.004dia.nc', 'REUNI.1.00A20.004.nc']
 
 Dvar_input = {
-'f1':['ZS', 'UT', 'VT', 'WT', 'THT', 'ALT_PRESSURE','ALT_U','ALT_V','ALT_THETA','level','ZTOP', 'longitude','latitude','level_w','time'],
+'f1':['ZS', 'UT', 'VT', 'WT', 'THT', 'ALT_PRESSURE','ALT_U','ALT_V','ALT_THETA','level','ZTOP', 'longitude','latitude','level_w','time',
+    'RN', 'H','LE','GFLUX','HU2M','T2M','W10M','CD','CH','Z0','CE','TS','Z0H'],
 'f2':['LSTHM', 'LSVM']}
 
 #  Read the variables in the files
@@ -65,7 +66,7 @@ Ltitle = ['wind vectors at K=2', 'wind vectors at z = 1500m ']
 Lxlab = ['longitude']*len(Lplot1)
 Lylab = ['latitude']*len(Lplot1)
 Llegendval = [25,25]
-Lcbarlabel = ['(m/s)']*len(Lplot1)
+Llegendlabel = ['(m/s)']*len(Lplot1)
 Larrowstep = [4]*len(Lplot1)
 Lwidth = [0.003]*len(Lplot1)
 Lcolor = ['black']*len(Lplot1)
@@ -73,7 +74,7 @@ Lprojection = [ccrs.PlateCarree()]*len(Lplot1)
 Llvl = [0]*len(Lplot1)
 Lscale = [400]*len(Lplot1)
 fig2 = Panel1.pvector(Lxx=LaxeX, Lyy=LaxeY, Llevel=Llvl, Lvar1=Lplot1, Lvar2=Lplot2, Lxlab=Lxlab, Lylab=Lylab, Ltitle=Ltitle, Lwidth=Lwidth, Larrowstep=Larrowstep, 
-                      Llegendval=Llegendval, Lcbarlabel=Lcbarlabel, Lproj=Lprojection, Lid_overlap=[2,6], ax=fig1.axes, Lscale=Lscale)
+                      Llegendval=Llegendval, Llegendlabel=Llegendlabel, Lproj=Lprojection, Lid_overlap=[2,6], ax=fig1.axes, Lscale=Lscale)
 
 ################################################################
 #########          PANEL 2  # Vertical cross-section
@@ -122,7 +123,7 @@ Lplot1 = [ Dvar['f1']['VM'][:,:,i_slice]]
 Lplot2 = [ Dvar['f1']['WM'][:,:,i_slice]]
 Ltitle = ['Wind']
 Llegendval = [15]
-Lcbarlabel = ['m/s']*len(Lplot)
+Llegendlabel = ['m/s']*len(Lplot)
 Lxlab = ['longitude']*len(Lplot)
 Lylab = ['altitude (m)']*len(Lplot)
 Larrowstep = [1]*len(Lplot)
@@ -133,6 +134,6 @@ Lxlim = [(-21.3,-20.9)]*len(Lplot)
 Lcolor=['lightgray']
 
 fig4 = Panel2.pvector(Lxx=LaxeX, Lyy=LaxeZ, Lvar1=Lplot1, Lvar2=Lplot2, Lxlab=Lxlab, Lylab=Lylab, Ltitle=Ltitle, Lwidth=Lwidth, Larrowstep=Larrowstep, 
-                        Llegendval=Llegendval, Lcbarlabel=Lcbarlabel, Lid_overlap=[6], ax=fig3.axes, Lscale=Lscale, Lylim=Lylim, Lxlim=Lxlim, Lcolor=Lcolor)
+                        Llegendval=Llegendval, Llegendlabel=Llegendlabel, Lid_overlap=[6], ax=fig3.axes, Lscale=Lscale, Lylim=Lylim, Lxlim=Lxlim, Lcolor=Lcolor)
 
-Panel2.save_graph(2,fig4)
\ No newline at end of file
+Panel2.save_graph(2,fig4)
diff --git a/MY_RUN/KTEST/007_16janvier/010_python/plot_007_16janvier.py b/MY_RUN/KTEST/007_16janvier/010_python/plot_007_16janvier.py
index 03bec46158136aacebbe3bb40e9d216a9c16aa85..f5292605c7d9fe3c358ad25b49fc999224725bf9 100644
--- a/MY_RUN/KTEST/007_16janvier/010_python/plot_007_16janvier.py
+++ b/MY_RUN/KTEST/007_16janvier/010_python/plot_007_16janvier.py
@@ -58,7 +58,7 @@ Ltitle = ['Wind at 850hPa', 'Wind at 700hPa', 'Wind at 9000m']
 Lxlab = ['longitude']*len(Lplot1)
 Lylab = ['latitude']*len(Lplot1)
 Llegendval = [20,20,40]
-Lcbarlabel = ['(m/s)']*len(Lplot1)
+Llegendlabel = ['(m/s)']*len(Lplot1)
 Larrowstep = [2]*len(Lplot1)
 Lwidth = [0.002]*len(Lplot1)
 Lcolor = ['black']*len(Lplot1)
@@ -66,7 +66,7 @@ Lprojection = [ccrs.PlateCarree()]*len(Lplot1)
 Llvl = [0]*len(Lplot1)
 fig2 = Panel1.pvector(Lxx=lon, Lyy=lat, Lvar1=Lplot1, Lvar2=Lplot2, Lcarte=[], Llevel=Llvl, Lxlab=Lxlab, Lylab=Lylab, 
                       Ltitle=Ltitle, Lwidth=Lwidth, Larrowstep=Larrowstep, Lproj=Lprojection,
-                      Lcolor=Lcolor, Llegendval=Llegendval, Lcbarlabel=Lcbarlabel, Lid_overlap=[2,4,6], ax=fig1.axes)
+                      Lcolor=Lcolor, Llegendval=Llegendval, Llegendlabel=Llegendlabel, Lid_overlap=[2,4,6], ax=fig1.axes)
 
 Panel1.save_graph(1,fig2)
 
@@ -99,7 +99,7 @@ Ltitle = ['Wind at 850hPa', 'Wind at 700hPa', 'Wind at 9000m']
 Llegendval = [20,20,40]
 Lxlab = ['longitude']*len(Lplot1)
 Lylab = ['latitude']*len(Lplot1)
-Lcbarlabel = ['(m/s)']*len(Lplot1)
+Llegendlabel = ['(m/s)']*len(Lplot1)
 Larrowstep = [2]*len(Lplot1)
 Lwidth = [0.002]*len(Lplot1)
 Lcolor = ['black']*len(Lplot1)
@@ -107,6 +107,6 @@ Lprojection = [ccrs.PlateCarree()]*len(Lplot1)
 Llvl = [0]*len(Lplot1)
 fig2 = Panel2.pvector(Lxx=lon, Lyy=lat, Lvar1=Lplot1, Lvar2=Lplot2, Lcarte=[], Llevel=Llvl, Lxlab=Lxlab, Lylab=Lylab, 
                       Ltitle=Ltitle, Lwidth=Lwidth, Larrowstep=Larrowstep, Lproj=Lprojection,
-                      Lcolor=Lcolor, Llegendval=Llegendval, Lcbarlabel=Lcbarlabel, Lid_overlap=[2,4,6], ax=fig1.axes)
+                      Lcolor=Lcolor, Llegendval=Llegendval, Llegendlabel=Llegendlabel, Lid_overlap=[2,4,6], ax=fig1.axes)
 
 Panel2.save_graph(2,fig2)
diff --git a/MY_RUN/KTEST/011_KW78CHEM/004_python/plot_011_KW78CHEM.py b/MY_RUN/KTEST/011_KW78CHEM/004_python/plot_011_KW78CHEM.py
index 3f792a31363a78a8d5db2c9e5163afa1597afc82..9666ee0c3614aea18c5a55add922bbde82015df1 100644
--- a/MY_RUN/KTEST/011_KW78CHEM/004_python/plot_011_KW78CHEM.py
+++ b/MY_RUN/KTEST/011_KW78CHEM/004_python/plot_011_KW78CHEM.py
@@ -68,7 +68,7 @@ Ltitle = ['Wind at K=2', 'Wind at 3000m', 'Wind at 5000m']
 Lxlab = ['x (m)']*len(Lplot)
 Lylab = ['y (m)']*len(Lplot)
 Llegendval = [10,10,10]
-Lcbarlabel = ['m/s']*len(Lplot)
+Llegendlabel = ['m/s']*len(Lplot)
 Larrowstep = [1]*len(Lplot)
 Lwidth = [0.002]*len(Lplot)
 Lcolor = ['black']*len(Lplot)
@@ -78,7 +78,7 @@ lat = [Dvar['f1']['nj_u'], Dvar['f2']['nj'],  Dvar['f2']['nj'] ]
 Lscale = [200]*len(Lplot)
 fig2 = Panel1.pvector(Lxx=lon, Lyy=lat, Lvar1=Lplot1, Lvar2=Lplot2, Lcarte=[500,23500,500,23500], Llevel=Llvl, 
                       Lxlab=Lxlab, Lylab=Lylab, Ltitle=Ltitle, Lwidth=Lwidth, Larrowstep=Larrowstep, 
-                      Lcolor=Lcolor, Llegendval=Llegendval, Lcbarlabel=Lcbarlabel, Lid_overlap=[4,6,8], ax=fig1.axes, Lscale=Lscale)
+                      Lcolor=Lcolor, Llegendval=Llegendval, Llegendlabel=Llegendlabel, Lid_overlap=[4,6,8], ax=fig1.axes, Lscale=Lscale)
 #  Oblique projection
 i_beg, j_beg = (2,0)
 i_end, j_end = (21,22)
@@ -163,13 +163,13 @@ Lplot1 = [ WIND_sec1]
 Lplot2 = [ WT_sec1]
 Ltitle = ['Wind']
 Llegendval = [25]
-Lcbarlabel = ['m/s']*len(Lplot)
+Llegendlabel = ['m/s']*len(Lplot)
 Larrowstep = [1]*len(Lplot)
 Lwidth = [0.004]*len(Lplot)
 Lscale = [200]*len(Lplot)
 
 fig4 = Panel3.pvector(Lxx=LaxeX, Lyy=LaxeZ, Lvar1=Lplot1, Lvar2=Lplot2, Lxlab=Lxlab, Lylab=Lylab, Ltitle=Ltitle, Lwidth=Lwidth, Larrowstep=Larrowstep, 
-                        Llegendval=Llegendval, Lcbarlabel=Lcbarlabel, Lid_overlap=[0], ax=fig3.axes, Lscale=Lscale)
+                        Llegendval=Llegendval, Llegendlabel=Llegendlabel, Lid_overlap=[0], ax=fig3.axes, Lscale=Lscale)
 
 Lplot = [RRT_sec1]
 LaxeX = [axe_m1]
@@ -287,4 +287,4 @@ fig8 = Panel5.psectionV(Lxx=LaxeX, Lzz=LaxeZ, Lvar=Lplot, Lxlab=Lxlab, Lylab=Lyl
                                 Lstep=Lstep, Lstepticks=Lstepticks, LcolorLine=LcolorLine, Lcbarlabel=Lcbarlabel, Lfacconv=Lfacconv, 
                                 Lpltype=Lpltype, ax=fig7.axes,Lid_overlap=[0, 2, 4, 6, 8, 10, 12, 14, 16],colorbar=False)
 
-Panel5.save_graph(5,fig8)
\ No newline at end of file
+Panel5.save_graph(5,fig8)
diff --git a/MY_RUN/KTEST/012_dust/005_python/plot_012_dust.py b/MY_RUN/KTEST/012_dust/005_python/plot_012_dust.py
index 2116f12329b9163e64761721f006b0e9e06efd1f..34819dabe4071852324a264955e5668dc147fe89 100644
--- a/MY_RUN/KTEST/012_dust/005_python/plot_012_dust.py
+++ b/MY_RUN/KTEST/012_dust/005_python/plot_012_dust.py
@@ -64,14 +64,14 @@ Lplot1 = [ Dvar['f1']['UM']]
 Lplot2 = [ Dvar['f1']['VM']]
 Ltitle = [' vectors at 1st level']
 Llegendval = [20]
-Lcbarlabel = ['m/s']*len(Lplot)
+Llegendlabel = ['m/s']*len(Lplot)
 Larrowstep = [1]*len(Lplot)
 Lwidth = [0.002]*len(Lplot)
 Lcolor = ['black']*len(Lplot)
 Lscale = [200]*len(Lplot)
 
 fig2 = Panel1.pvector(Lxx=lon, Lyy=lat, Lvar1=Lplot1, Llevel=Llvl, Lvar2=Lplot2, Lxlab=Lxlab, Lylab=Lylab, Ltitle=Ltitle, Lwidth=Lwidth, Larrowstep=Larrowstep, 
-                        Llegendval=Llegendval, Lcbarlabel=Lcbarlabel, Lid_overlap=[0], ax=fig1.axes, Lscale=Lscale)
+                        Llegendval=Llegendval, Llegendlabel=Llegendlabel, Lid_overlap=[0], ax=fig1.axes, Lscale=Lscale)
 
 Panel1.save_graph(1,fig2)
 
@@ -121,4 +121,4 @@ fig4 = Panel3.psectionH(lon=lon, lat=lat, Lvar=Lplot, Lxlab=Lxlab, Lylab=Lylab,
                                 Lstep=[], Lstepticks=[], Lcolormap=Lcolormap, Lcbarlabel=Lcbarlabel,Lcarte=Lcarte,
                                 Ltime=Ltime, LaddWhite_cm=LaddWhite, Lproj=Lprojection, ax=[])
 
-Panel3.save_graph(3,fig4)
\ No newline at end of file
+Panel3.save_graph(3,fig4)
diff --git a/src/LIB/Python/Panel_Plot.py b/src/LIB/Python/Panel_Plot.py
index 328d26264ad1f209de65a0ecdc36f563e07dde70..2072d1cf9225c27f14b79849e8b11e5e49885ce8 100644
--- a/src/LIB/Python/Panel_Plot.py
+++ b/src/LIB/Python/Panel_Plot.py
@@ -9,7 +9,7 @@ MNH_LIC for details. version 1.
 @author: 07/2021 Quentin Rodier
 """
 import matplotlib as mpl
-mpl.use('Agg')
+#mpl.use('Agg')
 import matplotlib.pyplot as plt
 from matplotlib import cm
 from matplotlib.colors import ListedColormap
@@ -21,7 +21,7 @@ class PanelPlot():
     
     def __init__(self, nb_l, nb_c, Lfigsize, bigtitle, titlepad=40, minmaxpad=1.03, timepad=-0.06, lateralminmaxpad=0.86, 
                  labelcolorbarpad=6.0, colorbaraspect=20, colorbarpad=0.04, tickspad=0.8,
-                 minmaxTextSize=10, bigtitleSize=13, titleSize=12,
+                 minmaxTextSize=10, bigtitleSize=13, titleSize=12, legendSize=10,
                  xlabelSize=11, ylabelSize=11, timeSize=11, cbTicksLabelSize=11, cbTitleSize=11, xyTicksLabelSize=10, figBoxLinewidth=1,
                  xyTicksWidth=1, xyTicksLength=6):
 
@@ -45,6 +45,7 @@ class PanelPlot():
         self.titleSize = titleSize               #  Graph title fontsize
         self.xlabelSize = xlabelSize             #  X-label fontsize
         self.ylabelSize = ylabelSize             #  Y-label fontsize
+        self.legendSize = legendSize             #  X/Y plot legend fontsize
         self.timeSize = timeSize                 #  Time attribute of the graphs fontsize
         self.cbTicksLabelSize = cbTicksLabelSize #  Colorbar ticks label fontsize
         self.cbTitleSize = cbTitleSize           #  Colorbar title fontsize
@@ -57,7 +58,7 @@ class PanelPlot():
         #  Initialization of the panel plots
         self.fig = plt.figure(figsize=(self.Lfigsize[0],self.Lfigsize[1]))
         self.fig.set_dpi(125)
-        self.fig.suptitle(self.bigtitle,fontsize=16)
+        self.fig.suptitle(self.bigtitle,fontsize=bigtitleSize)
 
     def save_graph(self, iplt, fig, fig_name='tempgraph'):
       """
@@ -221,7 +222,7 @@ class PanelPlot():
                  Lstep=[], Lstepticks=[], Lcolormap=[], Lcbarlabel=[], LcolorLine=[],
                  Lfacconv=[], ax=[], Lid_overlap=[], colorbar=True, orog=[], Lxlim=[], Lylim=[], Ltime=[], Lpltype=[], LaddWhite_cm=[], LwhiteTop=[]):
       """
-        Horizontal cross section plot
+        Vertical cross section plot
         Parameters :
             - Lxx    : List of x or y coordinate variable or time axis
             - Lzz    : List of z coordinates variable
@@ -246,6 +247,7 @@ class PanelPlot():
             - colorbar   : show colorbar or not
             - LaddWhite_cm : List of boolean to add white color to a colormap at the last bottom (low value) tick colorbar
             - LwhiteTop    : List of boolean to add the white color at the first top (high value). If false, the white is added at the bottom if Laddwhite_cm=T
+            - orog         : Orography variable
       """
       self.ax = ax
       firstCall = (len(self.ax) == 0)
@@ -412,13 +414,16 @@ class PanelPlot():
         #  Legend
         #TODO : Handling legend with overlap two axis lines in the same box. For now, placement is by hand
         if not id_overlap: 
-          self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.95))
+          self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.95),fontsize=self.legendSize)
         else:
-          self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.90))
+          self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.90),fontsize=self.legendSize)
 
         #  Title
         if Ltitle: self.set_Title(self.ax, iax, Ltitle[i], id_overlap,Lxlab[i], Lylab[i])
     
+        #  Ticks label 
+        self.ax[iax].tick_params(axis='both', labelsize=self.xyTicksLabelSize, width=self.xyTicksWidth, length=self.xyTicksLength, pad=self.tickspad)
+        
         #  X/Y Axis label
         if id_overlap: 
           self.ax[iax].xaxis.tick_top()
@@ -577,7 +582,9 @@ class PanelPlot():
     
         #  X/Y Axis
         self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i])
-        
+        #  Ticks label 
+        self.ax[iax].tick_params(axis='both', labelsize=self.xyTicksLabelSize, width=self.xyTicksWidth, length=self.xyTicksLength, pad=self.tickspad)
+       
         #  Color label on contour-line
         if Lpltype[i]=='c': #  Contour
           if 'GeoAxes' in str(type(self.ax[self.i])): # cartopy does not like the levels arguments in clabel, known issue
@@ -594,11 +601,11 @@ class PanelPlot():
       return self.fig
     
     def pvector(self, Lxx=[], Lyy=[], Lvar1=[], Lvar2=[], Lcarte=[], Llevel=[], Lxlab=[], Lylab=[], 
-                Ltitle=[], Lwidth=[], Larrowstep=[], Lcolor=[], Llegendval=[], Lcbarlabel=[], 
+                Ltitle=[], Lwidth=[], Larrowstep=[], Lcolor=[], Llegendval=[], Llegendlabel=[], 
                 Lproj=[], Lfacconv=[], ax=[], coastLines=True, Lid_overlap=[], Ltime=[], Lscale=[],
                 Lylim=[], Lxlim=[]):
       """
-        Horizontal vectors lines
+        Vectors
         Parameters :
             - Lxx    : List of x or y coordinate variable (lat or ni or nm)
             - Lyy    : List of y coordinates variable (lon or level)
@@ -617,7 +624,7 @@ class PanelPlot():
             - Larrowstep : List of sub-sample (frequency) if too much arrows
             - Lcolor : List of colors for the arrows (default: black)
             - Llegendval : List of value for the legend of the default arrow
-            - Lcbarlabel : List of labels for the legend of the default arrow
+            - Llegendlabel : List of labels for the legend of the default arrow
             - Lproj      : List of ccrs cartopy projection
             - Lfacconv   : List of factors for unit conversion of each variables
             - coastLines : Boolean to plot coast lines and grid lines
@@ -703,7 +710,7 @@ class PanelPlot():
         if Lproj: self.draw_Backmap(coastLines, self.ax[iax], Lproj[i])
     
         # Arrow legend key
-        qk = self.ax[iax].quiverkey(cf, 1.0, -0.05, Llegendval[i], str(Llegendval[i]) + Lcbarlabel[i], labelpos='E', color='black')
+        qk = self.ax[iax].quiverkey(cf, 1.0, -0.05, Llegendval[i], str(Llegendval[i]) + Llegendlabel[i], labelpos='E', color='black')
            
       return self.fig
 
diff --git a/src/LIB/Python/read_MNHfile.py b/src/LIB/Python/read_MNHfile.py
index fba716062817184b910fdab3d1bcf2d7bd414082..38072becffe3b1d8c066c6f370683b4988cf7ccf 100644
--- a/src/LIB/Python/read_MNHfile.py
+++ b/src/LIB/Python/read_MNHfile.py
@@ -104,10 +104,10 @@ def read_netcdf(LnameFiles, Dvar_input, path='.', get_data_only=True, del_empty_
                 Dvar[keyFiles] = read_TIMESfiles_55(theFile, Dvar_input[keyFiles], Dvar[keyFiles], get_data_only, del_empty_dim, removeHALO)
         else:
             Dvar[keyFiles]= read_BACKUPfile(theFile, Dvar_input[keyFiles], Dvar[keyFiles], get_data_only, del_empty_dim, removeHALO)
-        theFile.close()
+        #theFile.close()
     return Dvar
 
-def read_var(theFile, Dvar, var_name, get_data_only=True, del_empty_dim=True, removeHALO=True):
+def read_var(theFile, Dvar, var_name, get_data_only, del_empty_dim=True, removeHALO=True):
     """Read a netCDF4 variable
     
     Parameters
@@ -141,7 +141,6 @@ def read_var(theFile, Dvar, var_name, get_data_only=True, del_empty_dim=True, re
         var_dim_name = theFile.variables[var_name].dimensions
     except:
         raise KeyError("Group and variable name not found in the file, check the group/variable name with ncdump -h YourMNHFile.000.nc. You asked for variable : " + var_name)
-    
     if not get_data_only:
         Dvar[var_name] = theFile.variables[var_name]
     else: 
@@ -180,11 +179,10 @@ def read_var(theFile, Dvar, var_name, get_data_only=True, del_empty_dim=True, re
                 except IndexError:
                     break
             Ldimtosqueeze=tuple(Ldimtosqueeze)
-            Dvar[var_name] = np.squeeze(Dvar[var_name], axis=Ldimtosqueeze) 
-        
+            Dvar[var_name] = np.squeeze(Dvar[var_name], axis=Ldimtosqueeze)
     return Dvar
 
-def read_from_group(theFile, Dvar, group_name, var_name, get_data_only=True, del_empty_dim=True,removeHALO=True):
+def read_from_group(theFile, Dvar, group_name, var_name, get_data_only, del_empty_dim=True,removeHALO=True):
     """Read a variable from a netCDF4 group 
     
     Parameters
@@ -267,7 +265,7 @@ def read_from_group(theFile, Dvar, group_name, var_name, get_data_only=True, del
             Dvar[(group_name,var_name)] = Dvar[(group_name,var_name)].T   
     return Dvar
 
-def read_BACKUPfile(theFile, Dvar_input, Dvar, get_data_only=True, del_empty_dim=True, removeHALO=True):
+def read_BACKUPfile(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=True, removeHALO=True):
     """Read variables from Meso-NH MASDEV >= 5.5.0 synchronous file
     For all variables in Dvar_input of one file, call functions to read the variable of the group+variable
     
@@ -301,19 +299,18 @@ def read_BACKUPfile(theFile, Dvar_input, Dvar, get_data_only=True, del_empty_dim
     #  Reading date since beginning of the model run
     Dvar['time'] = theFile.variables['time'][0]
     Dvar['date'] = nc.num2date(Dvar['time'],units=theFile.variables['time'].units, calendar = theFile.variables['time'].calendar)
-           
     for var in Dvar_input:
         if type(var) == tuple:
             Dvar = read_from_group(theFile, Dvar, var[0], var[1], get_data_only, del_empty_dim, removeHALO)
         else:
             Dvar = read_var(theFile, Dvar, var, get_data_only, del_empty_dim, removeHALO)
-            
             #  For all variables except scalars, change Fill_Value to NaN
-            Dvar[var]= np.where(Dvar[var] != -99999.0, Dvar[var], np.nan)
-            Dvar[var]= np.where(Dvar[var] != 999.0, Dvar[var], np.nan)
+            if get_data_only:
+                Dvar[var]= np.where(Dvar[var] != -99999.0, Dvar[var], np.nan)
+                Dvar[var]= np.where(Dvar[var] != 999.0, Dvar[var], np.nan)
     return Dvar
 
-def read_TIMESfiles_55(theFile, Dvar_input, Dvar, get_data_only=True, del_empty_dim=True, removeHALO=True):
+def read_TIMESfiles_55(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=True, removeHALO=True):
     """Read variables from Meso-NH MASDEV >= 5.5.0 diachronic file
     For all variables in Dvar_input of one file, call functions to read the variable of the group+variable
 
@@ -345,7 +342,6 @@ def read_TIMESfiles_55(theFile, Dvar_input, Dvar, get_data_only=True, del_empty_
     Dvar[ifile][('group_name','var_name')] if the group contains more than one variable
     """  
     for var in Dvar_input:
-        print(var)
         if type(var) == tuple:
             Dvar = read_from_group(theFile, Dvar, var[0], var[1], get_data_only, del_empty_dim, removeHALO)
         else:
@@ -354,6 +350,8 @@ def read_TIMESfiles_55(theFile, Dvar_input, Dvar, get_data_only=True, del_empty_
 
 def removetheHALO(idim, var):
     """Remove a NHALO=1 point [1:-1] at a given dimension idim of a variable var
+    TODO: with NHALO=3 with WENO5
+    TODO: with the use of get_data_only=False (NetCDF4 output)
     
     Parameters
     ----------
diff --git a/src/LIB/Python/readme b/src/LIB/Python/readme
new file mode 100644
index 0000000000000000000000000000000000000000..ddb601ec1ec0ed8b757ed310dda73811edb0840d
--- /dev/null
+++ b/src/LIB/Python/readme
@@ -0,0 +1 @@
+Tutorial given the 17/03/2022 to use the library can be found at https://github.com/QuentinRodier/MNHPy/blob/main/examples/mnhpy.ipynb