How to use the neuron.h function in NEURON

def test():
            from blenderneuron.quick import bn

            with Blender(keep=False):

                from neuron import h
                h.load_file(test_hoc_file)
                tc = h.TestCell()

                ic = h.IClamp(0.5, sec=tc.soma)
                ic.delay = 1
                ic.dur = 3
                ic.amp = 0.5

                bn.prepare_for_collection()
                h.run()
                bn.to_blender()

                bn.run_command("bpy.data.scenes['Scene'].frame_current = 1;")
                self.assertTrue(bn.run_command("return_value = bpy.data.materials['TestCell[0].dendrites[9][0]'].emit") == 0.0)

                bn.run_command("bpy.data.scenes['Scene'].frame_current = 7;")
                self.assertTrue(bn.run_command("return_value = bpy.data.materials['TestCell[0].dendrites[9][0]'].emit") == 2.0)

def plot_raster(self, color='blue'):
        """ Plot raster with spikes of all cells """
        pyplot.figure()
        pyplot.scatter(self.spkt, self.spkid, marker= "|", s=100, c=color)
        pyplot.xlabel('time (ms)')
        pyplot.ylabel('cell id')
        pyplot.title('Network raster')
        pyplot.show()

# Main code
net = Net(numcells=10)  # create network 
net.connect_cells_ring(syn_weight=0.1, syn_delay=1)  # connect cells in a ring 
net.plot_net()  # plot network cells positions
net.cells[0].add_current_stim(delay=1)  # add stimulation to a cell

h.tstop = 60 # set simulation duration
#h.init()  # initialize sim
#h.run()  # run simulation

def analysis():
    from matplotlib import pyplot

    net.cells[0].plot_voltage()  # plot voltage
    net.cells[4].plot_voltage()  # plot voltage

    # plot raster
    net.plot_raster()

def run(self) :
    h.nrniv_bind_thread(threading.current_thread().ident);
    self.started = True;
    while True:
      self.fun()
      time.sleep(self.interval)

_species_get_all_species = species._get_all_species
_node_get_states = node._get_states
_section1d_transfer_to_legacy = section1d._transfer_to_legacy
_ctypes_c_int = ctypes.c_int
_weakref_ref = weakref.ref

_external_solver = None
_external_solver_initialized = False
_windows_dll_files = []
_windows_dll = []



make_time_ptr = nrn_dll_sym('make_time_ptr')
make_time_ptr.argtypes = [ctypes.py_object, ctypes.py_object]
make_time_ptr(h._ref_dt, h._ref_t)

_double_ptr = ctypes.POINTER(ctypes.c_double)
_int_ptr = ctypes.POINTER(_ctypes_c_int)
_long_ptr = ctypes.POINTER(ctypes.c_long)


fptr_prototype = ctypes.CFUNCTYPE(None)
set_nonvint_block = nrn_dll_sym('set_nonvint_block')
set_nonvint_block(nrn_dll_sym('rxd_nonvint_block'))

set_setup = nrn_dll_sym('set_setup')
set_setup.argtypes = [fptr_prototype]
set_initialize = nrn_dll_sym('set_initialize')
set_initialize.argtypes = [fptr_prototype]

scatter_concentrations = nrn_dll_sym('scatter_concentrations')

def _transfer_to_legacy():
    global  _c_ptr_vector, _c_ptr_vector_storage, _c_ptr_vector_storage_nrn
    global _last_c_ptr_length
    
    size = len(_all_cptrs)
    if _last_c_ptr_length != size:
        if size:
            _c_ptr_vector = h.PtrVector(size)
            _c_ptr_vector.ptr_update_callback(_donothing)
            for i, ptr in enumerate(_all_cptrs):
                _c_ptr_vector.pset(i, ptr)
            _c_ptr_vector_storage_nrn = h.Vector(size)
            _c_ptr_vector_storage = _c_ptr_vector_storage_nrn.as_numpy()
        else:
            _c_ptr_vector = None
        _last_c_ptr_length = size
    if size:
        _c_ptr_vector_storage[:] = node._get_states()[_all_cindices]
        _c_ptr_vector.scatter(_c_ptr_vector_storage_nrn)

def find_parent_seg(join, sdict, objects):
    
    if not join:
        return None
    elif join[0] not in objects:
        pseg = sdict[(join[0]._x0, join[0]._y0, join[0]._z0, join[0]._x1, join[0]._y1, join[0]._z1)]
        # better be all in same cell; so just set root once
        h.distance(0, h.SectionRef(sec=pseg.sec).root(0))
        closest = h.distance(pseg)

    # any other possible instance?

    for item in join:
        if item not in objects:
            s = sdict[(item._x0, item._y0, item._z0, item._x1, item._y1, item._z1)]
            d = h.distance(s)
            if d < closest:
                pseg = s
                closest = d

    return pseg

''' show graphs 
    Parameters
    ----------
    v_vec: h.Vector()
        recorded voltage
    t_vec: h.Vector()
        recorded time
    '''
    dend_plot = pyplot.plot(t_vec, v_vec)
    pyplot.xlabel('time (ms)')
    pyplot.ylabel('mV')

if __name__ == '__main__':
    numofmodel = 8
    cell = cfiber(250, 1, 0, 15000, True, numofmodel)
    for sec in h.allsec():
        h.psection(sec=sec) #show parameters of each section
    branch_vec, t_vec = set_recording_vectors(cell.branch)
    print(cell.numofmodel)
    simulate(cell)
    show_output(branch_vec, t_vec)
    pyplot.show()

def build_subsets(self):
        self.all = h.SectionList()
        for sec in h.allsec():
            self.all.append(sec=sec)    
    def __init__(self):

# create section
            if sectName not in self.secs:
                self.secs[sectName] = Dict()  # create sect dict if doesn't exist
            if 'hObj' not in self.secs[sectName] or self.secs[sectName]['hObj'] in [None, {}, []]: 
                self.secs[sectName]['hObj'] = h.Section(name=sectName, cell=self)  # create h Section object
            sec = self.secs[sectName]  # pointer to section

            # set geometry params 
            if 'geom' in sectParams:
                for geomParamName,geomParamValue in sectParams['geom'].items():  
                    if not type(geomParamValue) in [list, dict]:  # skip any list or dic params
                        setattr(sec['hObj'], geomParamName, geomParamValue)

                # set 3d geometry
                if 'pt3d' in sectParams['geom']:  
                    h.pt3dclear(sec=sec['hObj'])
                    if sim.cfg.pt3dRelativeToCellLocation:
                        x = self.tags['x']
                        y = -self.tags['y'] if sim.cfg.invertedYCoord else self.tags['y'] # Neuron y-axis positive = upwards, so assume pia=0 and cortical depth = neg
                        z = self.tags['z']
                    else:
                        x = y = z = 0
                    for pt3d in sectParams['geom']['pt3d']:
                        h.pt3dadd(x+pt3d[0], y+pt3d[1], z+pt3d[2], pt3d[3], sec=sec['hObj'])

            # add distributed mechanisms 
            if 'mechs' in sectParams:
                for mechName,mechParams in sectParams['mechs'].items(): 
                    if mechName not in sec['mechs']: 
                        sec['mechs'][mechName] = Dict()
                    try:
                        sec['hObj'].insert(mechName)