def METHOD_NAME(model):
    # Location box
    geo_model = model
    geo_model._rescaling.toggle_axial_anisotropy()
    # gp.compute_model(geo_model, compute_mesh_options={'mask_topography': False})
    geo_model.surface_points.df[['X', 'Y', 'Z']] = geo_model.surface_points.df[['X_c', 'Y_c',
                                                                                'Z_c']]
    geo_model.orientations.df[['X', 'Y', 'Z']] = geo_model.orientations.df[['X_c', 'Y_c',
                                                                            'Z_c']]
    # This is a hack
    geo_model._grid.topography.extent = geo_model._grid.extent_c
    geo_model.set_regular_grid(geo_model._grid.extent_c, [50, 50, 50])
    geo_model.modify_kriging_parameters('range', 0.1)
    geo_model.modify_kriging_parameters('drift equations', [9, 9, 9, 9, 9])
    geo_model.modify_surface_points(
        geo_model.surface_points.df.index,
        smooth=0.001
    )
    gp.set_interpolator(geo_model, aesara_optimizer='fast_run', dtype='float64')
    gp.compute_model(geo_model, compute_mesh_options={'mask_topography': False,
                                                      'masked_marching_cubes': False})
    gp.plot_2d(geo_model,
               section_names=['topography'],
               show_topography=True,
               )
    plt.show()
    gp.plot_3d(geo_model, ve=None, show_topography=False, image=True, show_lith=False,
               kwargs_plot_data={'arrow_size': 10}
               )