diff --git a/pyVlsv/reduction.py b/pyVlsv/reduction.py
index 7ee6d9ab..f47d9502 100644
--- a/pyVlsv/reduction.py
+++ b/pyVlsv/reduction.py
@@ -686,6 +686,14 @@ def gyrophase_relstddev( variables, velocity_cell_data, velocity_coordinates ):
    histo = pl.hist(gyrophase_data[0].data, weights=gyrophase_data[1].data, bins=36, range=[-180.0,180.0], log=False, normed=1)
    return np.std(histo[0])/np.mean(histo[0])
 
+def vspace_dummy( variables, velocity_cell_data, velocity_coordinates ):
+   # This reducer needs to be verified
+   rho = variables[0]
+
+   print('velocity_cell_data.shape = '.format(velocity_cell_data.shape))
+   print('velocity_coordinates.shape = '.format(velocity_coordinates.shape))
+   return np.sum(velocity_cell_data)
+
 def Dng( variables ):
    # This reducer needs to be verified
    # This routine is still very slow due to for-loops
@@ -1194,6 +1202,9 @@ def makelambda(index):
 v5reducers["vg_jacobian_bper"] =          DataReducerVariable(["vg_dperbxvoldx","vg_dperbxvoldy","vg_dperbxvoldz","vg_dperbyvoldx","vg_dperbyvoldy","vg_dperbyvoldz","vg_dperbzvoldx","vg_dperbzvoldy","vg_dperbzvoldz"], TensorFromScalars, "T/m", 1, latex=r"$\vec{J}$",latexunits=r"$\mathrm{A}\,\mathrm{m}^{-2}$")
 v5reducers["vg_j"] =                     DataReducerVariable(["vg_jacobian_bper"], J, "A/m^2", 1, latex=r"$\vec{J}$",latexunits=r"$\mathrm{A}\,\mathrm{m}^{-2}$")
 
+#reducers with useVspace
+v5reducers["vg_vspace_dummy"] =    DataReducerVariable(["vg_rho"], vspace_dummy, "", 1, useVspace=True) # I think this has vector length 1?
+
 #multipopv5reducers
 multipopv5reducers = {}
 multipopv5reducers["pop/vg_rhom"] =                   DataReducerVariable(["pop/vg_rho"], rhom, "kg/m3", 1, latex=r"$\rho_{m,\mathrm{REPLACEPOP}}$",latexunits=r"$\mathrm{kg}\,\mathrm{m}^{-3}$")
diff --git a/pyVlsv/vlsvreader.py b/pyVlsv/vlsvreader.py
index 6d2b8f27..19f61df7 100644
--- a/pyVlsv/vlsvreader.py
+++ b/pyVlsv/vlsvreader.py
@@ -853,7 +853,7 @@ def read_attribute(self, name="", mesh="", attribute="", tag=""):
       raise ValueError("Variable or attribute not found")
 
 
-   def read(self, name="", tag="", mesh="", operator="pass", cellids=-1):
+   def read(self, name="", tag="", mesh="", operator="pass", cellids=-1, pop='proton'):
       ''' Read data from the open vlsv file. 
       
       :param name: Name of the data array
@@ -1030,7 +1030,7 @@ def read(self, name="", tag="", mesh="", operator="pass", cellids=-1):
             output = np.zeros(len(actualcellids))
             index = 0
             for singlecellid in actualcellids:
-               velocity_cell_data = self.read_velocity_cells(singlecellid)
+               velocity_cell_data = self.read_velocity_cells(singlecellid, pop=pop)
                # Get cells:
                vcellids = list(velocity_cell_data.keys())
                # Get coordinates: