def generate(self,rprs=None, mass=None, radius=None,
                n=2e4, fp_specific=0.01, u1=None, u2=None,
                 starmodel=None,
                Teff=None, logg=None, rbin_width=0.3,
                MAfn=None, lhoodcachefile=None):
        """Generates Population

        All arguments defined in ``__init__``.
        """

        n = int(n)

        if starmodel is None:
            if type(mass) is type((1,)):
                mass = dists.Gaussian_Distribution(*mass)
            if isinstance(mass, dists.Distribution):
                mdist = mass
                mass = mdist.rvs(1e5)

            if type(radius) is type((1,)):
                radius = dists.Gaussian_Distribution(*radius)
            if isinstance(radius, dists.Distribution):
                rdist = radius
                radius = rdist.rvs(1e5)
        else:
            samples = starmodel.random_samples(1e5)
            mass = samples['mass_0_0'].values
            radius = samples['radius_0_0'].values
            Teff = samples['Teff_0_0'].mean()
            logg = samples['logg_0_0'].mean()

        logging.debug('star mass: {}'.format(mass))
        logging.debug('star radius: {}'.format(radius))
        logging.debug('Teff: {}'.format(Teff))
        logging.debug('logg: {}'.format(logg))

        if u1 is None or u2 is None:
            if Teff is None or logg is None:
                logging.warning('Teff, logg not provided; using solar limb darkening')
                u1 = 0.394; u2=0.296
            else:
                u1,u2 = ldcoeffs(Teff, logg)

        #use point estimate of rprs to construct planets in radius bin
        #rp = self.rprs*np.median(radius)
        #rbin_min = (1-rbin_width)*rp
        #rbin_max = (1+rbin_width)*rp

        rprs_bin_min = (1-rbin_width)*self.rprs
        rprs_bin_max = (1+rbin_width)*self.rprs

        radius_p = radius * (np.random.random(int(1e5))*(rprs_bin_max - rprs_bin_min) + rprs_bin_min)
        mass_p = (radius_p*RSUN/REARTH)**2.06 * MEARTH/MSUN #hokey, but doesn't matter

        logging.debug('planet radius: {}'.format(radius_p))

        stars = pd.DataFrame()
        #df_orbpop = pd.DataFrame() #for orbit population

        tot_prob = None; tot_dprob = None; prob_norm = None
        n_adapt = n
        while len(stars) < n:
            n_adapt = int(n_adapt)
            inds = np.random.randint(len(mass), size=n_adapt)

            #calculate eclipses.
            ecl_inds, df, (prob,dprob) = calculate_eclipses(mass[inds], mass_p[inds],
                                                        radius[inds], radius_p[inds],
                                                        15, np.inf, #arbitrary
                                                        u11s=u1, u21s=u2,
                                                        band=self.band,
                                                        period=self.period,
                                                        calc_mininc=True,
                                                        return_indices=True,
                                                        MAfn=MAfn)

            df['mass_A'] = mass[inds][ecl_inds]
            df['mass_B'] = mass_p[inds][ecl_inds]
            df['radius_A'] = radius[inds][ecl_inds]
            df['radius_B'] = radius_p[inds][ecl_inds]
            df['u1'] = u1 * np.ones_like(df['mass_A'])
            df['u2'] = u2 * np.ones_like(df['mass_A'])
            df['P'] = self.period * np.ones_like(df['mass_A'])

            ok = (df['dpri']>0) & (df['T14_pri'] > 0)

            stars = pd.concat((stars, df[ok]))

            logging.info('{} Transiting planet systems generated (target {})'.format(len(stars),n))
            logging.debug('{} nans in stars[dpri]'.format(np.isnan(stars['dpri']).sum()))

            if tot_prob is None:
                prob_norm = (1/dprob**2)
                tot_prob = prob
                tot_dprob = dprob
            else:
                prob_norm = (1/tot_dprob**2 + 1/dprob**2)
                tot_prob = (tot_prob/tot_dprob**2 + prob/dprob**2)/prob_norm
                tot_dprob = 1/np.sqrt(prob_norm)

            n_adapt = min(int(1.2*(n-len(stars)) * n_adapt//len(df)), 5e4)
            n_adapt = max(n_adapt, 100)

        stars = stars.reset_index()
        stars.drop('index', axis=1, inplace=True)
        stars = stars.iloc[:n]

        stars['mass_1'] = stars['mass_A']
        stars['radius_1'] = stars['radius_A']
        stars['mass_2'] = stars['mass_B']
        stars['radius_2'] = stars['radius_B']

        #make OrbitPopulation?

        #finish below.

        if fp_specific is None:
            rp = stars['radius_2'].mean() * RSUN/REARTH
            fp_specific = fp_fressin(rp)

        priorfactors = {'fp_specific':fp_specific}

        self._starmodel = starmodel

        EclipsePopulation.__init__(self, stars=stars,
                                   period=self.period, cadence=self.cadence,
                                   model=self.model,
                                   priorfactors=priorfactors, prob=tot_prob,
                                   lhoodcachefile=lhoodcachefile)