diff --git a/pyPlots/plot_colormap.py b/pyPlots/plot_colormap.py
index 2717c046..fe1134c9 100644
--- a/pyPlots/plot_colormap.py
+++ b/pyPlots/plot_colormap.py
@@ -63,7 +63,7 @@ def plot_colormap(filename=None,
                   symmetric=False,
                   pass_vars=None, pass_times=None, pass_full=False,
                   fluxfile=None, fluxdir=None, flux_levels=None,
-                  fluxthick=1.0, fluxlines=1,
+                  fluxthick=1.0, fluxlines=1,fluxlinecolor='k',
                   fsaved=None,
                   nomask=False,
                   Earth=None,
@@ -176,6 +176,7 @@ def plot_colormap(filename=None,
         :kword fluxdir:     Directory in which fluxfunction files can be found
         :kword fluxthick:   Scale fluxfunction line thickness
         :kword fluxlines:   Relative density of fluxfunction contours
+        :kword fluxlinecolor: Set flux line color (default black)
         :kword fsaved:      Overplot locations of fSaved. If keyword is set to a string, that will be the colour used.
         :kword nomask:      Do not mask plotting based on proton density
 
@@ -220,7 +221,6 @@ def exprMA_cust(exprmaps, requestvariables=False):
                 MA = custombulkspeed/va
                 return MA
             plot_colormap(filename=fileLocation, vmin=1 vmax=40, expression=exprMA_cust,lin=True)
-
     '''
 
     # Switch None-keywords to empty lists (this way subsequent calls get correct empty default values
@@ -1017,7 +1017,7 @@ def exprMA_cust(exprmaps, requestvariables=False):
             flux_levels = np.linspace(-10,10,fluxlines*60)
         else:
             pass #This was given, do nothing
-        fluxcont = ax1.contour(XmeshCentres,YmeshCentres,flux_function,flux_levels,colors='k',linestyles='solid',linewidths=0.5*fluxthick,zorder=2)
+        fluxcont = ax1.contour(XmeshCentres,YmeshCentres,flux_function,flux_levels,colors=fluxlinecolor,linestyles='solid',linewidths=0.5*fluxthick,zorder=2)
 
     # add fSaved identifiers
     if fsaved:
diff --git a/pyVlsv/reduction.py b/pyVlsv/reduction.py
index 8db29400..fc95aad1 100644
--- a/pyVlsv/reduction.py
+++ b/pyVlsv/reduction.py
@@ -1148,6 +1148,9 @@ def makelambda(index):
 v5reducers["vg_egradpe_parallel"] =              DataReducerVariable(["vg_e_gradpe", "vg_b_vol"], ParallelVectorComponent, "V/m", 1, latex=r"$E_{\nabla Pe,\parallel}$",latexunits=r"$\mathrm{V}\,\mathrm{m}^{-1}$")
 v5reducers["vg_egradpe_perpendicular"] =         DataReducerVariable(["vg_e_gradpe", "vg_b_vol"], PerpendicularVectorComponent, "V/m", 1, latex=r"$E_{\nabla Pe,\perp}$",latexunits=r"$\mathrm{V}\,\mathrm{m}^{-1}$")
 
+v5reducers["vg_hallterm"] = DataReducerVariable(["vg_e_vol", "vg_v", "vg_b_vol"], Hallterm, "V/m", 3, latex=r"$E_\mathrm{Hall}$", latexunits=r"\mathrm{V}\,\mathrm{m}^{-1}$")
+
+
 v5reducers["vg_pdyn"] =            DataReducerVariable(["vg_v", "vg_rhom"], Pdyn, "Pa", 1, latex=r"$P_\mathrm{dyn}$",latexunits=r"$\mathrm{Pa}$")
 v5reducers["vg_pdynx"] =            DataReducerVariable(["vg_v", "vg_rhom"], Pdynx, "Pa", 1, latex=r"$P_\mathrm{dyn,x}$",latexunits=r"$\mathrm{Pa}$")