How to use the phoebe.backend.universe function in phoebe

lc['indep_type'] = 'time (hjd)'
   
    #-- then create a BodyBag from WD: we want to make sure the output here
    #   is the same as before
    comp1,comp2,binary = wd.wd_to_phoebe(ps,lc,rv)
    star1,lcdep1,rvdep1 = comp1
    star2,lcdep2,rvdep2 = comp2
    mesh1 = phoebe.ParameterSet(frame='phoebe',context='mesh:marching',delta=0.2,alg='c')
    mesh2 = phoebe.ParameterSet(frame='phoebe',context='mesh:marching',delta=0.2,alg='c')
    star1 = phoebe.BinaryRocheStar(star1,binary,mesh1,pbdep=[lcdep1,rvdep1])
    star2 = phoebe.BinaryRocheStar(star2,binary,mesh2,pbdep=[lcdep2,rvdep2])
    wd_bbag = phoebe.BodyBag([star1,star2])
    
    #-- so write a file to compare the two (that's up to you...)
    curve,params = wd.lc(ps,request='curve',light_curve=lc,rv_curve=rv)
    serial_legacy = universe.serialize(system,color=False)
    serial_wildev = universe.serialize(wd_bbag,color=False)
    with open(phoebe_out_file,'w') as ff:
        for line1,line2 in zip(serial_legacy.split('\n'),serial_wildev.split('\n')):
            ff.write('PH:'+line1+'\n')
            ff.write('WD:'+line2+'\n')
    #============ COMPUTATIONS ===============================================
    #-- get mpi-stuff and details, but only if we're not profiling the code.
    if 'cProfile' in globals():
        mpi = None
    else:
        mpi = get_mpi_parameters()

    #-- compute the system if the light curves haven't been computed before
    if not os.path.isfile(phoebe_lc_file) or recompute:
        #-- compute the system
        #time = time[::5]

def accept_fit(system,fitparams):
    """
    Accept the MCMC fit.
    """
    feedback = fitparams['feedback']
    #-- which parameters were fitted?
    fitted_param_labels = [par.get_unique_label() for par in feedback['parameters']]
    #-- walk through all the parameterSets available:
    walk = utils.traverse(system,list_types=(universe.BodyBag,universe.Body,list,tuple),dict_types=(dict,))
    for parset in walk:
        #-- fore ach parameterSet, walk to all the parameters
        for qual in parset:
            #-- extract those which need to be fitted
            if parset.get_adjust(qual) and parset.has_prior(qual):
                this_param = parset.get_parameter(qual)
                #-- ask a unique ID and update the value of the parameter
                myid = this_param.get_unique_label()
                index = fitted_param_labels.index(myid)
                #-- [iwalker,trace,iparam]
                this_param.set_value(feedback['values'][index])
                logger.info("Set {} = {} from {} fit".format(qual,this_param.as_string(),fitparams['method']))

"""
    if fitparams is None:
        fitparams = parameters.ParameterSet(frame='phoebe',context='fitting:lmfit')
    # We need unique names for the parameters that need to be fitted, we need
    # initial values and identifiers to distinguish parameters with the same
    # name (we'll also use the identifier in the parameter name to make sure
    # they are unique). While we are iterating over all the parameterSets,
    # we'll also have a look at what context they are in. From this, we decide
    # which fitting algorithm to use.
    ids = []
    pars = lmfit.Parameters()
    ppars = [] # Phoebe parameters
    #-- walk through all the parameterSets available. This needs to be via
    #   this utility function because we need to iteratively walk down through
    #   all BodyBags too.
    walk = utils.traverse(system,list_types=(universe.BodyBag,universe.Body,list,tuple),dict_types=(dict,))
    frames = []
    for parset in walk:
        frames.append(parset.frame)
        #-- for each parameterSet, walk through all the parameters
        for qual in parset:
            #-- extract those which need to be fitted
            if parset.get_adjust(qual) and parset.has_prior(qual):
                #-- ask a unique ID and check if this parameter has already
                #   been treated. If so, continue to the next one.
                parameter = parset.get_parameter(qual)
                myid = parameter.get_unique_label()
                if myid in ids: continue
                #-- else, add the name to the list of pnames. Ask for the
                #   prior of this object
                name = '{}_{}'.format(qual,myid).replace('-','_')
                prior = parameter.get_prior()

def _create_system_and_compute_pblums(self, b, compute,
                                          dynamics_method=None,
                                          hier=None,
                                          meshablerefs=None,
                                          datasets=None,
                                          compute_l3=True,
                                          compute_l3_frac=False,
                                          compute_extrinsic=False,
                                          reset=True,
                                          lc_only=True,
                                          **kwargs):

        logger.debug("rank:{}/{} PhoebeBackend._create_system_and_compute_pblums: calling universe.System.from_bundle".format(mpi.myrank, mpi.nprocs))
        system = universe.System.from_bundle(b, compute, datasets=b.datasets, **kwargs)

        if dynamics_method is None:
            computeparams = b.get_compute(compute, force_ps=True)
            dynamics_method = computeparams.get_value(qualifier='dynamics_method', **kwargs)

        if hier is None:
            hier = b.get_hierarchy()

        if meshablerefs is None:
            starrefs  = hier.get_stars()
            meshablerefs = hier.get_meshables()

        t0 = b.get_value(qualifier='t0', context='system', unit=u.d, t0=kwargs.get('t0', None), **_skip_filter_checks)

        if len(meshablerefs) > 1 or hier.get_kind_of(meshablerefs[0])=='envelope':
            logger.debug("rank:{}/{} PhoebeBackend._create_system_and_compute_pblums: computing dynamics at t0".format(mpi.myrank, mpi.nprocs))

mesh1 = parameters.ParameterSet(frame='phoebe', context='mesh:marching',
                                delta=delta, alg='c')
    mesh2 = parameters.ParameterSet(frame='phoebe', context='mesh:marching',
                                delta=delta, alg='c')
    curve, params = wd.lc(ps, request='curve', light_curve=lc, rv_curve=rv)
    lcobs = datasets.LCDataSet(columns=['time','flux','sigma'], time=curve['indeps'],
                             sigma=0.001*curve['lc'],
                             flux=curve['lc'], ref=lcdep1['ref'])
    rvobs1 = datasets.RVDataSet(columns=['time','rv','sigma'], time=curve['indeps'],
                             sigma=0.001*curve['rv1']*100,
                             rv=curve['rv1']*100, ref=rvdep1['ref'])
    rvobs2 = datasets.RVDataSet(columns=['time','rv','sigma'], time=curve['indeps'],
                             sigma=0.001*curve['rv2']*100,
                             rv=curve['rv2']*100, ref=rvdep2['ref'])

    star1 = universe.BinaryRocheStar(star1, orbit, mesh1, pbdep=[lcdep1, rvdep1], obs=[rvobs1])
    star2 = universe.BinaryRocheStar(star2, orbit, mesh2, pbdep=[lcdep2, rvdep2], obs=[rvobs2])
    system = universe.BodyBag([star1, star2], position=position, obs=[lcobs], label='system')
    logger.info("Successfully parsed WD lcin file {}".format(filename))
    return system

current_string = []
            
            
        if item['kind'] == 'Parameter':
            if summary_type == 'full' and not it.get_hidden():
                current_string.append(make_param_label(it, green))
            elif summary_type == 'only_adjust' and it.get_adjust():
                current_string.append(make_param_label(it, green))
            elif summary_type == 'only_adjust':
                pass
            else:
                current_string.append('')
            
            # Get the number of Bodies in the path
            bodies = [body for body in path if isinstance(body, universe.Body)]
            level = len(bodies)
        
        
        # Make sure that the Bodies appear in the right order
        if item['kind'] == 'Body':
            this_bodies = [body for body in path if isinstance(body, universe.Body)]
            this_level = len(this_bodies)
            this_label = make_body_label(this_bodies[this_level-1], this_level, emphasize)
                
            if not this_label in total_string:
                total_string[this_label] = OrderedDict()
                total_obsdep[this_label] = OrderedDict()
    
    # Take care of final thing
    else:
        if current_string:

def walk(mybundle):
    
    for val,path in utils.traverse_memory(mybundle,
                                     list_types=(Bundle, universe.Body, list,tuple),
                                     dict_types=(dict, ),
                                     parset_types=(parameters.ParameterSet, ),
                                     get_label=(universe.Body, ),
                                     get_context=(parameters.ParameterSet, ),
                                     skip=()):
            
            # First one is always root
            path[0] = str(mybundle.__class__.__name__)
            
            
            # All is left is to return it
            
            yield path, val

comp1,comp2,orbit = binary_from_stars(star, be_star, period=(10.,'d'))

    orbit['ecc'] = 0.60
    orbit['per0'] = 20.,'deg'
    orbit['t0'] = conversions.convert('CD','JD',(1914.,8,6.15)) #1914.6
    orbit['long_an'] = long,'deg'
    orbit['incl'] = 100,'deg'
    orbit['c1label'] = 'betacep'
    orbit['c2label'] = 'Be-star_system'
    
    
    betacep = universe.BinaryStar(star, mesh=mesh_betacep, puls=[freq1],
                                magnetic_field=mag_field1, orbit=orbit)
    bestar = universe.Star(be_star, mesh_bestar, magnetic_field=mag_field2)#, circ_spot=[spot])
    mydisk = universe.AccretionDisk(diskpars, mesh=mesh_disk)
    bestar = universe.BodyBag([bestar, mydisk], orbit=orbit, label='Be-star_system')
    system = universe.BodyBag([betacep,bestar], label='system')
    
    return system

# advancing to deeper levels; as long as we don't encounter
                # that reference/label/context, we can't look for the next one
                index = start_index
                
                # Look down the tree structure
                for jlevel, level in enumerate(path):
                    
                    # but only if we still have structure information
                    if index < len(structure_info):
                        
                        # We don't now the name of this level yet, but we'll
                        # figure it out
                        name_of_this_level = None
                        
                        # If this level is a Body, we check it's label
                        if isinstance(level, universe.Body):
                            name_of_this_level = level.get_label()
                            
                        # If it is a ParameterSet, we'll try to match the label,
                        # reference or context
                        elif isinstance(level, parameters.ParameterSet):
                            if 'ref' in level:
                                name_of_this_level = level['ref']
                            elif 'label' in level:
                                name_of_this_level = level['label']
                            if name_of_this_level != structure_info[index]:
                                name_of_this_level = level.get_context()
                        
                        # The walk iterator could also give us 'lcsyn' or something
                        # the like, it apparently doesn't walk over these PSsets
                        # themselves -- but no problem, this works
                        elif isinstance(level, str):

logger.warning('Starting from previous results {}'.format(fitparams['label']))
    else:
        db = 'pickle'
        logger.warning("Starting new chain {}".format(fitparams['label']))
    
    mc = pymc.MCMC(pars,db=db,dbname=dbname)
    #-- store the info in a feedback dictionary
    feedback = dict(parset=fitparams,parameters=[],traces=[],priors=[])
    
    mc.sample(iter=fitparams['iters'],burn=fitparams['burn'],thin=fitparams['thin'])
    if fitparams['incremental']:
        mc.db.close()
    #-- add the posteriors to the parameters
    had = []
    #-- walk through all the parameterSets available:
    walk = utils.traverse(system,list_types=(universe.BodyBag,universe.Body,list,tuple),dict_types=(dict,))
    for parset in walk:
        #-- fore ach parameterSet, walk to all the parameters
        for qual in parset:
            #-- extract those which need to be fitted
            if parset.get_adjust(qual) and parset.has_prior(qual):
                #-- ask a unique ID and add the parameter
                myid = parset.get_parameter(qual).get_unique_label()
                if myid in had: continue
                had.append(myid)
                #-- add the trace, priors and the parameters themselves to
                #   the dictionary
                this_param = parset.get_parameter(qual)
                #-- access all samples, also from previous runs
                trace = mc.trace('{}_{}'.format(qual,myid),chain=None)[:]
                feedback['parameters'].append(this_param)
                feedback['traces'].append(trace)