Python fitbase.FitContribution类代码示例

def makeRecipe(molecule, datname):
    """Create a recipe that uses the DebyePDFGenerator."""

    ## The Profile
    profile = Profile()

    # Load data and add it to the profile
    profile.loadtxt(datname)
    profile.setCalculationRange(xmin=1.2, xmax=8)

    ## The ProfileGenerator
    # Create a DebyePDFGenerator named "G".
    generator = DebyePDFGenerator("G")
    generator.setStructure(molecule)
    # These are metadata needed by the generator
    generator.setQmin(0.68)
    generator.setQmax(22)

    ## The FitContribution
    contribution = FitContribution("bucky")
    contribution.addProfileGenerator(generator)
    contribution.setProfile(profile, xname = "r")

    # Make a FitRecipe.
    recipe = FitRecipe()
    recipe.addContribution(contribution)

    # Specify which parameters we want to refine.
    c60 = generator.phase

    # We're not going to refine the ADPs. However, every atom must have a
    # small, but finite ADP for the PDF calculator to work properly.
    atoms = c60.getScatterers()
    for atom in atoms:
        atom.Uiso.value = 0.001

    # Scale factor. We cannot optimize this efficiently, so we take the value
    # from a previous refinement. In general, care must be taken to properly
    # determine the scale of the profile, or to make sure that the residual is
    # not affected by the scale.
    generator.scale.value = 1.24457360e+4

    # Allow every atom to move.  We define the bounds to be a window of radius
    # 0.1 centered on the current value of the position.
    win = 0.1
    for idx, atom in enumerate(atoms):
        xname, yname, zname = getXYZNames(idx)
        recipe.addVar(atom.x, name = xname).boundWindow(win)
        recipe.addVar(atom.y, name = yname).boundWindow(win)
        recipe.addVar(atom.z, name = zname).boundWindow(win)

    return recipe

def makeRecipe():
    """Create a recipe that uses the GaussianGenerator.

    This will create a FitContribution that uses the GaussianGenerator,
    associate this with a Profile, and use this to define a FitRecipe.

    """

    ## The Profile
    # Create a Profile to hold the experimental and calculated signal.
    profile = Profile()

    # Load data and add it to the profile. This uses the loadtxt function from
    # numpy.
    profile.loadtxt("data/gaussian.dat")

    ## The ProfileGenerator
    # Create a GaussianGenerator named "g". This will be the name we use to
    # refer to the generator from within the FitContribution equation.  
    generator = GaussianGenerator("g")
    
    ## The FitContribution
    # Create a FitContribution that will associate the Profile with the
    # GaussianGenerator.  The GaussianGenerator will be accessible as an
    # attribute of the FitContribution by its name ("g"). Note that this will
    # set the fitting equation to "g", which calls the GaussianGenerator.
    contribution = FitContribution("g1")
    contribution.addProfileGenerator(generator)
    contribution.setProfile(profile)

    ## The FitRecipe
    # Now we create the FitRecipe and add the FitContribution.
    recipe = FitRecipe()
    recipe.addContribution(contribution)

    # Specify which Parameters we want to vary in the fit.  This will add
    # Variables to the FitRecipe that directly modify the Parameters of the
    # FitContribution.
    # 
    # We create a variable for each Parameter of the GaussianGenerator. Note
    # that the Parameters belong to the GaussianGenerator, not the
    # FitContribution as in gaussianrecipe.py. We initialize parameters as in
    # gaussianrecipe.py so we can expect the same output.
    recipe.addVar(generator.A, 1)
    recipe.addVar(generator.x0, 5)
    recipe.addVar(generator.sigma, name = "sig")
    recipe.sig.value = 1

    # Give the recipe away so it can be used!
    return recipe

def makeRecipe(x, y, dy, A, sig, x0):
    """Make a FitRecipe for fitting a Gaussian curve to data.    
    """
    profile = Profile()
    profile.setObservedProfile(x, y, dy)

    contribution = FitContribution("g1")
    contribution.setProfile(profile, xname="x")
    contribution.setEquation("A * exp(-0.5*(x-x0)**2/sigma**2)")

    recipe = FitRecipe()
    recipe.addContribution(contribution)
    recipe.addVar(contribution.A, A)
    recipe.addVar(contribution.x0, x0)
    recipe.addVar(contribution.sigma, sig)
    return recipe

    def __init__(self, name):
        """Create the PDFContribution.

        name        --  The name of the contribution.

        """
        FitContribution.__init__(self, name)
        self._meta = {}
        # Add the profile
        profile = Profile()
        self.setProfile(profile, xname = "r")

        # Need a parameter for the overall scale, in the case that this is a
        # multi-phase fit.
        self.newParameter("scale", 1.0)
        # Profile-related parameters that will be shared between the generators
        self.newParameter("qdamp", 0)
        self.newParameter("qbroad", 0)
        return

def makeRecipe(datname):
    """Create a fitting recipe for ellipsoidal SAS data."""

    ## The Profile
    # This will be used to store the observed and calculated I(Q) data.
    profile = Profile()

    # Load data and add it to the Profile. We use a SASParser to load the data
    # properly and pass the metadata along.
    parser = SASParser()
    parser.parseFile(datname)
    profile.loadParsedData(parser)

    ## The ProfileGenerator
    # The SASGenerator is for configuring and calculating a SAS profile. We use
    # a sans model to configure and serve as the calculation engine of the
    # generator. This allows us to use the full sans model creation
    # capabilities, and tie this into SrFit when we want to fit a model to
    # data. The documentation for the various sans models can be found at
    # http://danse.chem.utk.edu/sansview.html.
    from sans.models.EllipsoidModel import EllipsoidModel
    model = EllipsoidModel()
    generator = SASGenerator("generator", model)

    ## The FitContribution
    # Here we associate the Profile and ProfileGenerator, as has been done
    # before. 
    contribution = FitContribution("ellipsoid")
    contribution.addProfileGenerator(generator)
    contribution.setProfile(profile, xname = "q")

    # We want to fit the log of the signal to the log of the data so that the
    # higher-Q information remains significant. There are no I(Q) uncertainty
    # values with the data, so we do not need to worry about the effect this
    # will have on the estimated parameter uncertainties.
    contribution.setResidualEquation("log(eq) - log(y)")

    ## Make the FitRecipe and add the FitContribution.
    recipe = FitRecipe()
    recipe.addContribution(contribution)

    ## Configure the fit variables
    # The SASGenerator uses the parameters from the params and dispersion
    # attribues of the model. These vary from model to model, but are adopted
    # as SrFit Parameters within the generator. Whereas the dispersion
    # parameters are accessible as, e.g. "radius.width", within the
    # SASGenerator these are named like "radius_width".
    #
    # We want to fit the scale factor, radii and background factors.
    recipe.addVar(generator.scale, 1)
    recipe.addVar(generator.radius_a, 50)
    recipe.addVar(generator.radius_b, 500)
    recipe.addVar(generator.background, 0)
    
    # Give the recipe away so it can be used!
    return recipe

def main():

    p = Profile()
    p.loadtxt("data/gaussian.dat")

    # FitContribution operations
    # "<<"  -   Inject a parameter value
    c = FitContribution("g1")
    c.setProfile(p)
    c.setEquation("A * exp(-0.5*(x-x0)**2/sigma**2)")
    c.A << 0.5
    c.x0 << 5
    c.sigma << 1

    # FitRecipe operations
    # "|="  -   Union of necessary components.
    # "+="  -   Add Parameter or create a new one. Each tuple is a set of
    #           arguments for either setVar or addVar.
    # "*="  -   Constrain a parameter. Think of "*" as a push-pin holding one
    #           parameter's value to that of another.
    # "%="  -   Restrain a parameter or equation. Think of "%" as a rope
    #           loosely tying parameters to a value.
    r = FitRecipe()
    r |= c
    r += (c.A, 0.5), (c.x0, 5), 'sig'
    r *= c.sigma, 'sig'
    r %= c.A, 0.5, 0.5

    from gaussianrecipe import scipyOptimize
    scipyOptimize(r)

    res = FitResults(r)
    # Print the results.
    res.printResults()
    # Plot the results.
    from gaussianrecipe import plotResults
    plotResults(r)

    return

def makeRecipe(ciffile, grdata):
    """Make a recipe to model a crystal-like nanoparticle PDF."""

    # Set up a PDF fit as has been done in other examples.
    pdfprofile = Profile()

    pdfparser = PDFParser()
    pdfparser.parseFile(grdata)
    pdfprofile.loadParsedData(pdfparser)
    pdfprofile.setCalculationRange(xmin = 0.1, xmax = 20)

    pdfcontribution = FitContribution("pdf")
    pdfcontribution.setProfile(pdfprofile, xname = "r")

    pdfgenerator = PDFGenerator("G")
    pdfgenerator.setQmax(30.0)
    stru = CreateCrystalFromCIF(file(ciffile))
    pdfgenerator.setStructure(stru)
    pdfcontribution.addProfileGenerator(pdfgenerator)

    # Register the nanoparticle shape factor.
    from diffpy.srfit.pdf.characteristicfunctions import sphericalCF
    pdfcontribution.registerFunction(sphericalCF, name = "f")

    # Now we set up the fitting equation.
    pdfcontribution.setEquation("f * G")

    # Now make the recipe. Make sure we fit the characteristic function shape
    # parameters, in this case 'psize', which is the diameter of the particle.
    recipe = FitRecipe()
    recipe.addContribution(pdfcontribution)

    phase = pdfgenerator.phase
    for par in phase.sgpars:
        recipe.addVar(par)

    recipe.addVar(pdfcontribution.psize, 20)
    recipe.addVar(pdfgenerator.scale, 1)
    recipe.addVar(pdfgenerator.delta2, 0)
    recipe.B11_0 = 0.1

    return recipe

 def _makeRecipe(self, x, y, dy):
     '''Make a FitRecipe for fitting a Gaussian curve to data.
     '''
     profile = Profile()
     profile.setObservedProfile(x, y, dy)
     contribution = FitContribution("g1")
     contribution.setProfile(profile, xname="x")
     contribution.registerStringFunction(
             '1/sqrt(2 * pi * sig**2)', name='gaussnorm')
     contribution.setEquation(
             "A * gaussnorm * exp(-0.5 * (x - x0)**2/sig**2)")
     recipe = FitRecipe()
     recipe.addContribution(contribution)
     recipe.addVar(contribution.A)
     recipe.addVar(contribution.x0)
     recipe.addVar(contribution.sig)
     recipe.clearFitHooks()
     self.recipe = recipe
     return

class TestWeakBoundMethod(unittest.TestCase):

    def setUp(self):
        self.f = FitContribution('f')
        self.f.setEquation('7')
        self.w = weak_ref(self.f._eq._flush, fallback=_fallback_example)
        return


    def tearDown(self):
        self.f = None
        self.assertTrue(None is self.w._wref())
        obj, args, kw = self.w('any', 'argument', foo=37)
        self.assertTrue(obj is self.w)
        self.assertEqual(('any', 'argument'), args)
        self.assertEqual({'foo' : 37}, kw)
        return


    def test___init__(self):
        """check WeakBoundMethod.__init__()
        """
        self.assertTrue(self.w.fallback is _fallback_example)
        wf = weak_ref(self.f._flush)
        self.assertTrue(None is wf.fallback)
        return


    def test___call__(self):
        """check WeakBoundMethod.__call__()
        """
        f = self.f
        self.assertEqual(7, f.evaluate())
        self.assertEqual(7, f._eq._value)
        # verify f has the same effect as f._eq._flush
        self.w(())
        self.assertTrue(None is f._eq._value)
        # check WeakBoundMethod behavior with no fallback
        x = Parameter('x', value=3)
        wgetx = weak_ref(x.getValue)
        self.assertEqual(3, wgetx())
        del x
        self.assertRaises(ReferenceError, wgetx)
        return


    def test___hash__(self):
        """check WeakBoundMethod.__hash__()
        """
        f1 = FitContribution('f1')
        w1 = weak_ref(f1._flush)
        h0 = hash(w1)
        del f1
        self.assertTrue(None is w1._wref())
        self.assertEqual(h0, hash(w1))
        w1c1 = pickle.loads(pickle.dumps(w1))
        w1c2 = pickle.loads(pickle.dumps(w1))
        self.assertEqual(hash(w1c1), hash(w1c2))
        return


    def test___eq__(self):
        """check WeakBoundMethod.__eq__()
        """
        f1 = FitContribution('f1')
        w1 = weak_ref(f1._flush)
        w2 = weak_ref(f1._flush)
        self.assertEqual(w1, w2)
        w1c = pickle.loads(pickle.dumps(w1))
        # pickle-copied objects should have empty reference
        self.assertTrue(None is w1c._wref())
        self.assertNotEqual(w1, w1c)
        del f1
        self.assertTrue(None is w1._wref())
        self.assertEqual(w1, w1c)
        w1cc = pickle.loads(pickle.dumps(w1c))
        self.assertTrue(None is w1cc._wref())
        self.assertEqual(w1c, w1cc)
        self.assertEqual(w1, w1cc)
        return


    def test_pickling(self):
        """Verify unpickling works when it involves __hash__ call.
        """
        holder = set([self.w])
        objs = [holder, self.f._eq, self.w]
        data = pickle.dumps(objs)
        objs2 = pickle.loads(data)
        h2, feq2, w2 = objs2
        self.assertTrue(w2 in h2)
        self.assertTrue(feq2 is w2._wref())
        return


    def test_observable_deregistration(self):
        """check if Observable drops dead Observer.
        """
        f = self.f
        x = f.newParameter('x', 5)
#.........这里部分代码省略.........

def makeRecipe(niciffile, siciffile, datname):
    """Create a fitting recipe for crystalline PDF data."""

    ## The Profile
    profile = Profile()

    # Load data and add it to the profile
    parser = PDFParser()
    parser.parseFile(datname)
    profile.loadParsedData(parser)
    profile.setCalculationRange(xmax = 20)

    ## The ProfileGenerator
    # In order to fit two phases simultaneously, we must use two PDFGenerators.
    # PDFGenerator is designed to take care of as little information as it
    # must. (Don't do too much, and do it well.) A PDFGenerator can generate
    # the signal from only a single phase at a time. So, we will create one
    # PDFGenerator for each phase and compose them within the same
    # FitContribution. Note that both generators will be associated with the
    # same Profile within the FitContribution, so they will both be
    # automatically configured according to the metadata.
    #
    # The generator for the nickel phase. We call it "G_ni" and will use this
    # name later when we set the fitting equation in the FitContribution.
    generator_ni = PDFGenerator("G_ni")
    stru = CreateCrystalFromCIF(file(niciffile))
    generator_ni.setStructure(stru)
    # The generator for the silicon phase. We call it "G_si".
    generator_si = PDFGenerator("G_si")
    stru = CreateCrystalFromCIF(file(siciffile))
    generator_si.setStructure(stru)

    ## The FitContribution
    # Add both generators to the FitContribution. Add the Profile. This will
    # send the metadata to the generators.
    contribution = FitContribution("nisi")
    contribution.addProfileGenerator(generator_ni)
    contribution.addProfileGenerator(generator_si)
    contribution.setProfile(profile, xname = "r")

    # Write the fitting equation. We want to sum the PDFs from each phase and
    # multiply it by a scaling factor. We also want a certain phase scaling
    # relationship between the PDFs which we will enforce with constraints in
    # the FitRecipe.
    contribution.setEquation("scale * (G_ni +  G_si)")

    # Make the FitRecipe and add the FitContribution.
    recipe = FitRecipe()
    recipe.addContribution(contribution)

    ## Configure the fit variables
    # Start by configuring the scale factor and resolution factors.
    # We want the sum of the phase scale factors to be 1.
    recipe.newVar("scale_ni", 0.1)
    recipe.constrain(generator_ni.scale, "scale_ni")
    recipe.constrain(generator_si.scale, "1 - scale_ni")
    # We also want the resolution factor to be the same on each.
    recipe.newVar("qdamp", 0.03)
    recipe.constrain(generator_ni.qdamp, "qdamp")
    recipe.constrain(generator_si.qdamp, "qdamp")

    # Vary the gloabal scale as well.
    recipe.addVar(contribution.scale, 1)

    # Now we can configure the structural parameters. Since we're using
    # ObjCrystCrystalParSets, the space group constraints are automatically
    # applied to each phase. We must selectively vary the free parameters.
    #
    # First the nickel parameters
    phase_ni = generator_ni.phase
    for par in phase_ni.sgpars:
        recipe.addVar(par, name = par.name + "_ni")
    recipe.addVar(generator_ni.delta2, name = "delta2_ni")
    # Next the silicon parameters
    phase_si = generator_si.phase
    for par in phase_si.sgpars:
        recipe.addVar(par, name = par.name + "_si")
    recipe.addVar(generator_si.delta2, name = "delta2_si")

    # We have prior information from the earlier examples so we'll use it here
    # in the form of restraints.
    #
    # The nickel lattice parameter was measured to be 3.527. The uncertainty
    # values are invalid for that measurement, since the data from which it is
    # derived has no uncertainty. Thus, we will tell the recipe to scale the
    # residual, which means that it will be weighted as much as the average
    # data point during the fit.
    recipe.restrain("a_ni", lb = 3.527, ub = 3.527, scaled = True)
    # Now we do the same with the delta2 and Biso parameters (remember that
    # Biso = 8*pi**2*Uiso)
    recipe.restrain("delta2_ni", lb = 2.22, ub = 2.22, scaled = True)
    recipe.restrain("Biso_0_ni", lb = 0.454, ub = 0.454, scaled = True)
    #
    # We can do the same with the silicon values. We haven't done a thorough
    # job of measuring the uncertainties in the results, so we'll scale these
    # as well.
    recipe.restrain("a_si", lb = 5.430, ub = 5.430, scaled = True)
    recipe.restrain("delta2_si", lb = 3.54, ub = 3.54, scaled = True)
    recipe.restrain("Biso_0_si", lb = 0.645, ub = 0.645, scaled = True)

#.........这里部分代码省略.........

def makeRecipe(stru1, stru2, datname):
    """Create a fitting recipe for crystalline PDF data."""

    ## The Profile
    profile = Profile()

    # Load data and add it to the profile
    parser = PDFParser()
    parser.parseFile(datname)
    profile.loadParsedData(parser)
    profile.setCalculationRange(xmin=1.5, xmax = 45, dx = 0.1)

    ## The ProfileGenerator
    # In order to fit the core and shell phases simultaneously, we must use two
    # PDFGenerators.
    #
    # The generator for the CdS core. We call it "G_CdS" and will use this name
    # later when we set the fitting equation in the FitContribution.
    generator_cds = PDFGenerator("G_CdS")
    generator_cds.setStructure(stru1)
    generator_cds.setQmax(26)
    generator_cds.qdamp.value = 0.0396
    # The generator for the ZnS shell. We call it "G_ZnS".
    generator_zns = PDFGenerator("G_ZnS")
    generator_zns.setStructure(stru2)
    generator_zns.setQmax(26)
    generator_zns.qdamp.value = 0.0396

    ## The FitContribution
    # Add both generators and the profile to the FitContribution.
    contribution = FitContribution("cdszns")
    contribution.addProfileGenerator(generator_cds)
    contribution.addProfileGenerator(generator_zns)
    contribution.setProfile(profile, xname = "r")

    # Set up the characteristic functions. We use a spherical CF for the core
    # and a spherical shell CF for the shell. Since this is set up as two
    # phases, we implicitly assume that the core-shell correlations contribute
    # very little to the PDF.
    from diffpy.srfit.pdf.characteristicfunctions import sphericalCF, shellCF
    contribution.registerFunction(sphericalCF, name = "f_CdS")
    contribution.registerFunction(shellCF, name = "f_ZnS")

    # Write the fitting equation. We want to sum the PDFs from each phase and
    # multiply it by a scaling factor.
    contribution.setEquation("scale * (f_CdS * G_CdS +  f_ZnS * G_ZnS)")

    # Make the FitRecipe and add the FitContribution.
    recipe = FitRecipe()
    recipe.addContribution(contribution)

    # Vary the inner radius and thickness of the shell. Constrain the core
    # diameter to twice the shell radius.
    recipe.addVar(contribution.radius, 15)
    recipe.addVar(contribution.thickness, 11)
    recipe.constrain(contribution.psize, "2 * radius")

    ## Configure the fit variables
    # Start by configuring the scale factor and resolution factors.
    # We want the sum of the phase scale factors to be 1.
    recipe.newVar("scale_CdS", 0.7)
    recipe.constrain(generator_cds.scale, "scale_CdS")
    recipe.constrain(generator_zns.scale, "1 - scale_CdS")
    # We also want the resolution factor to be the same on each.

    # Vary the gloabal scale as well.
    recipe.addVar(contribution.scale, 0.3)

    # Now we can configure the structural parameters. We tag the different
    # structural variables so we can easily turn them on and off in the
    # subsequent refinement.
    phase_cds = generator_cds.phase
    for par in phase_cds.sgpars.latpars:
        recipe.addVar(par, name = par.name + "_cds", tag = "lat")
    for par in phase_cds.sgpars.adppars:
        recipe.addVar(par, 1, name = par.name + "_cds", tag = "adp")
    recipe.addVar(phase_cds.sgpars.xyzpars.z_1, name = "z_1_cds", tag = "xyz")
    # Since we know these have stacking disorder, constrain the B33 adps for
    # each atom type.
    recipe.constrain("B33_1_cds", "B33_0_cds")
    recipe.addVar(generator_cds.delta2, name = "delta2_cds", value = 5)

    phase_zns = generator_zns.phase
    for par in phase_zns.sgpars.latpars:
        recipe.addVar(par, name = par.name + "_zns", tag = "lat")
    for par in phase_zns.sgpars.adppars:
        recipe.addVar(par, 1, name = par.name + "_zns", tag = "adp")
    recipe.addVar(phase_zns.sgpars.xyzpars.z_1, name = "z_1_zns", tag = "xyz")
    recipe.constrain("B33_1_zns", "B33_0_zns")
    recipe.addVar(generator_zns.delta2, name = "delta2_zns", value = 2.5)

    # Give the recipe away so it can be used!
    return recipe

def makeRecipe(ciffile, datname):
    """Create a fitting recipe for crystalline PDF data."""

    ## The Profile
    # This will be used to store the observed and calculated PDF profile.
    profile = Profile()

    # Load data and add it to the Profile. As before we use a PDFParser. The
    # metadata is still passed to the PDFGenerator later on. The interaction
    # between the PDFGenerator and the metadata does not depend on type of
    # structure being refined.
    parser = PDFParser()
    parser.parseFile(datname)
    profile.loadParsedData(parser)
    profile.setCalculationRange(xmax = 20)

    ## The ProfileGenerator
    # This time we use the CreateCrystalFromCIF method of pyobjcryst.crystal to
    # create a Crystal object. That object is passed to the PDFGenerator as in
    # the previous example.
    generator = PDFGenerator("G")
    stru = CreateCrystalFromCIF(file(ciffile))
    generator.setStructure(stru)
    generator.setQmax(40.0)
    
    ## The FitContribution
    contribution = FitContribution("nickel")
    contribution.addProfileGenerator(generator)
    contribution.setProfile(profile, xname = "r")

    # Make the FitRecipe and add the FitContribution.
    recipe = FitRecipe()
    recipe.addContribution(contribution)

    ## Configure the fit variables

    # As before, we get a handle to the structure parameter set. In this case,
    # it is a ObjCrystCrystalParSet instance that was created when we called
    # 'setStructure' above. The ObjCrystCrystalParSet has different Parameters
    # and options than the DiffpyStructureParSet used in the last example. See
    # its documentation in diffpy.srfit.structure.objcrystparset.
    phase = generator.phase

    # Here is where we created space group constraints in the previous example.
    # The difference in this example is that the ObjCrystCrystalParSet is aware
    # of space groups, and the DiffpyStructureParSet is not. Constraints are
    # created internally when "sgpars" attribute is called for. These
    # constriants get enforced within the ObjCrystCrystalParSet. Free
    # Parameters are stored within the 'sgpars' member of the
    # ObjCrystCrystalParSet, which is the same as the object returned from
    # 'constrainAsSpaceGroup'.
    #
    # As before, we have one free lattice parameter ('a'). We can simplify
    # things by iterating through all the sgpars.
    for par in phase.sgpars: 
        recipe.addVar(par)
    # set the initial thermal factor to a non-zero value
    assert hasattr(recipe, 'B11_0')
    recipe.B11_0 = 0.1

    # We now select non-structural parameters to refine.
    # This controls the scaling of the PDF.
    recipe.addVar(generator.scale, 1)
    # This is a peak-damping resolution term.
    recipe.addVar(generator.qdamp, 0.01)
    # This is a vibrational correlation term that sharpens peaks at low-r.
    recipe.addVar(generator.delta2, 5)

    # Give the recipe away so it can be used!
    return recipe

def makeRecipe():
    """Make a FitRecipe for fitting three double-gaussian curves to data.

    The separation and amplitude ratio of the double peaks follows a specific
    relationship.  The peaks are broadend according to their position and they
    sit on top of a background. We are seeking the absolute locations of the
    peaks as well as their amplitudes.

    The independent variable is t. The relationship between the double
    peaks is
    sin(t2) / l2 = sin(t1) / l1
    amplitude(peak2) = r * amplitude(peak1)
    The values of l1, l2 and r come from experiment. For this example, we
    use l1 = 1.012, l2 = 1.0 and r = 0.23.

    """

    ## The Profile
    # Create a Profile to hold the experimental and calculated signal.
    profile = Profile()
    x, y, dy = profile.loadtxt("data/threedoublepeaks.dat")

    # Create the contribution
    contribution = FitContribution("peaks")
    contribution.setProfile(profile, xname = "t")
    pi = numpy.pi
    exp = numpy.exp

    # This is a building-block of our profile function
    def gaussian(t, mu, sig):
        return 1/(2*pi*sig**2)**0.5 * exp(-0.5 * ((t-mu)/sig)**2)

    contribution.registerFunction(gaussian, name = "peakshape")

    def delta(t, mu):
        """Calculate a delta-function.

        We don't have perfect precision, so we must make this a very thin
        Gaussian.

        """
        sig = t[1] - t[0]
        return gaussian(t, mu, sig)

    contribution.registerFunction(delta)

    # Here is another one
    bkgdstr = "b0 + b1*t + b2*t**2 + b3*t**3 + b4*t**4 + b5*t**5 + b6*t**6"

    contribution.registerStringFunction(bkgdstr, "bkgd")

    # Now define our fitting equation. We will hardcode the peak ratios.
    contribution.setEquation(
        "A1 * ( convolve( delta(t, mu11), peakshape(t, c, sig11) ) \
         + 0.23*convolve( delta(t, mu12), peakshape(t, c, sig12) ) ) + \
         A2 * ( convolve( delta(t, mu21), peakshape(t, c, sig21) ) \
         + 0.23*convolve( delta(t, mu22), peakshape(t, c, sig22) ) ) + \
         A3 * ( convolve( delta(t, mu31), peakshape(t, c, sig31) ) \
         + 0.23*convolve( delta(t, mu32), peakshape(t, c, sig32) ) ) + \
         bkgd")

    # c is the center of the gaussian.
    contribution.c.value =  x[len(x)/2]

    ## The FitRecipe
    # The FitRecipe lets us define what we want to fit. It is where we can
    # create variables, constraints and restraints.
    recipe = FitRecipe()

    # Here we tell the FitRecipe to use our FitContribution. When the FitRecipe
    # calculates its residual function, it will call on the FitContribution to
    # do part of the work.
    recipe.addContribution(contribution)

    # Vary the amplitudes for each double peak
    recipe.addVar(contribution.A1, 100)
    recipe.addVar(contribution.A2, 100)
    recipe.addVar(contribution.A3, 100)

    # Vary the position of the first of the double peaks
    recipe.addVar(contribution.mu11, 13.0)
    recipe.addVar(contribution.mu21, 24.0)
    recipe.addVar(contribution.mu31, 33.0)

    # Constrain the position of the second double peak
    from numpy import sin, arcsin
    def peakloc(mu):
        """Calculate the location of the second peak given the first."""
        l1 = 1.012
        l2 = 1.0
        return 180 / pi * arcsin( pi / 180 * l2 * sin(mu) / l1 )

    recipe.registerFunction(peakloc)
    recipe.constrain(contribution.mu12, "peakloc(mu11)")
    recipe.constrain(contribution.mu22, "peakloc(mu21)")
    recipe.constrain(contribution.mu32, "peakloc(mu31)")

    # Vary the width of the peaks. We know the functional form of the peak
    # broadening.
    sig0 = recipe.newVar("sig0", 0.001)
#.........这里部分代码省略.........

def makeContribution(name, generator, profile):
    """Make a FitContribution and add a generator and profile."""
    contribution = FitContribution(name)
    contribution.addProfileGenerator(generator)
    contribution.setProfile(profile, xname = "r")
    return contribution

def makeRecipe(molecule, datname):
    """Create a recipe that uses the DebyePDFGenerator."""

    ## The Profile
    profile = Profile()

    # Load data and add it to the profile
    profile.loadtxt(datname)
    profile.setCalculationRange(xmin=1.2, xmax=8)

    ## The ProfileGenerator
    # Create a DebyePDFGenerator named "G".
    generator = DebyePDFGenerator("G")
    generator.setStructure(molecule)
    # These are metadata needed by the generator
    generator.setQmin(0.68)
    generator.setQmax(22)

    ## The FitContribution
    contribution = FitContribution("bucky")
    contribution.addProfileGenerator(generator)
    contribution.setProfile(profile, xname = "r")

    # Make a FitRecipe.
    recipe = FitRecipe()
    recipe.addContribution(contribution)

    # Specify which parameters we want to refine. We'll be using the
    # MoleculeParSet within the generator, so let's get a handle to it. See the
    # diffpy.srfit.structure.objcryststructure module for more information
    # about the MoleculeParSet hierarchy.
    c60 = generator.phase

    # First, the isotropic thermal displacement factor.
    Biso = recipe.newVar("Biso")
    for atom in c60.getScatterers():

        # We have defined a 'center' atom that is a dummy, which means that it
        # has no scattering power. It is only used as a reference point for
        # our bond length. We don't want to constrain it.
        if not atom.isDummy():
            recipe.constrain(atom.Biso, Biso)

    # We need to let the molecule expand. If we were modeling it as a crystal,
    # we could let the unit cell expand. For instruction purposes, we use a
    # Molecule to model C60, and molecules have different modeling options than
    # crystals. To make the molecule expand from a central point, we will
    # constrain the distance from each atom to a dummy center atom that was
    # created with the molecule, and allow that distance to vary. (We could
    # also let the nearest-neighbor bond lengths vary, but that would be much
    # more difficult to set up.)
    center = c60.center
    # Create a new Parameter that represents the radius of the molecule. Note
    # that we don't give it an initial value. Since the variable is being
    # directly constrained to further below, its initial value will be inferred
    # from the constraint.
    radius = recipe.newVar("radius")
    for i, atom in enumerate(c60.getScatterers()):

        if atom.isDummy():
            continue

        # This creates a Parameter that moves the second atom according to the
        # bond length. Note that each Parameter needs a unique name.
        par = c60.addBondLengthParameter("rad%i"%i, center, atom)
        recipe.constrain(par, radius)

    # Add the correlation term, scale. The scale is too short to effectively
    # determine qdamp.
    recipe.addVar(generator.delta2, 2)
    recipe.addVar(generator.scale, 1.3e4)

    # Give the recipe away so it can be used!
    return recipe