add optional grid solver output

floenavi/__init__.py

@@ -119,7 +119,8 @@ class FloeNavi(object):
         # --- Update the data set in place ----
         dset.update_floenavi_coords(ice_pos.xc, ice_pos.yc)
 
-    def get_virtual_reference_stations(self, dataset_time_bounds, basestation_setup=None, **refstat_kwargs):
+    def get_virtual_reference_stations(self, dataset_time_bounds, basestation_setup=None,
+                                       return_grid_solution=False, **refstat_kwargs):
         """
         Get two virtual reference stations (origin ref and heading ref) from merging
         the information of all available FloeNavi reference stations. The step is
@@ -131,8 +132,10 @@ class FloeNavi(object):
         position on the the FloeNavi grid.
 
         :param dataset_time_bounds: (datetime list) time range of the dataset
-        :param basestation_setup: (FloeNaviBaseStationSetup instance)
-        :return:
+        :param basestation_setup: (FloeNaviBaseStationSetup instance, optional)
+        :param return_grid_solution: (bool) If set to true, the method will return the refstation
+            and the grid solution
+        :return: icedrift.GeoReferenceStation (+ FloeNaviGrid if return_grid_solution is True)
         """
 
         # Use default base station setup
@@ -145,9 +148,14 @@ class FloeNavi(object):
         floenavi_grid.set_time_bounds(dataset_time_bounds)
 
         # Get Reference station data from two positions
-        pos_ref, pos_hdg = floenavi_grid.get_virtual_reference_station_positions(**refstat_kwargs)
+        pos_ref, pos_hdg, fngrid = floenavi_grid.get_virtual_reference_station_positions(**refstat_kwargs)
         refstat_ais = GeoReferenceStation.from_two_basestations(pos_ref, pos_hdg)
-        return refstat_ais
+
+        # Return value
+        if not return_grid_solution:
+            return refstat_ais
+        else:
+            return refstat_ais, fngrid
 
 
 class FloeNaviConfig(object):
@@ -732,7 +740,7 @@ class FloeNaviGrid(object):
         # Get the virtual origin and heading position data objects
         logger.info("Compute origin and heading station positions")
         orig, hdg = fngrid.get_reference_stations(**refstat_kwargs)
-        return orig, hdg
+        return orig, hdg, fngrid
 
 
 class FloeNaviGridSolver(object):