How to use the tractor.ellipses.EllipseE function in tractor

#     plt.clf()
        #     plt.imshow(mod20 - modDev, interpolation='nearest', origin='lower')
        #     plt.colorbar()
        #     plt.title('Sersic n=4 - DEV')
        #     ps.savefig()
        # 
        # self.assertTrue(np.abs(mod20.sum() - 99.95) < 1e-2)
        # target = modDev
        # print('SAD:', np.sum(np.abs(mod20 - target)))
        # self.assertTrue(np.sum(np.abs(mod20 - target)) < 1e-5)


        # A galaxy that will wrap around
        from tractor.ellipses import EllipseE
        mmX = ModelMask(20, 20, 60, 60)
        shapeX = EllipseE(20., 0.7, 0.7)
        tim.psf = pixpsf
        gal6 = DevGalaxy(pos0, bright, shapeX)

        p21 = gal6.getModelPatch(tim, modelMask=mmX)
        mod21 = np.zeros((H,W), np.float32)
        p21.addTo(mod21)
        
        if ps is not None:
            plt.clf()
            plt.imshow(mod21, interpolation='nearest', origin='lower')
            plt.colorbar()
            plt.title('Pixelized PSF')
            ps.savefig()

            
        tim.psf = hybridpsf

wanyband = 'w'

    print('Creating Tractor catalog...')
    cat = []
    for i, t in enumerate(T):
        pos = RaDecPos(t.ra, t.dec)
        flux = NanoMaggies(**{wanyband: 1.})
        if all_ptsrcs:
            cat.append(PointSource(pos, flux))
            continue

        tt = t.type.strip()
        if tt in ['PTSRC', 'STAR', 'S']:
            cat.append(PointSource(pos, flux))
        elif tt in ['EXP', 'E']:
            shape = EllipseE(*t.shapeexp)
            cat.append(ExpGalaxy(pos, flux, shape))
        elif tt in ['DEV', 'D']:
            shape = EllipseE(*t.shapedev)
            cat.append(DevGalaxy(pos, flux, shape))
        elif tt in ['COMP', 'C']:
            eshape = EllipseE(*t.shapeexp)
            dshape = EllipseE(*t.shapedev)
            cat.append(FixedCompositeGalaxy(pos, flux, t.fracdev,
                                            eshape, dshape))
        else:
            print('Did not understand row', i, 'of input catalog:')
            t.about()
            assert(False)

    W = unwise_forcedphot(cat, tiles, roiradecbox=roiradecbox,
                          bands=opt.bands, unwise_dir=opt.unwise_dir,

sources.append(PointSource(PixPos(pos_x, pos_y), Flux(obj['flux'])))

        print("Now you have %d sources" % len(sources))

    elif shape_method is 'decals':
        ## Using DECaLS shapes
        if sources is None:
            sources = []
        for obj in comp_galaxy:
            pos_x, pos_y = obj['x'], obj['y']
            sources.append(
                CompositeGalaxy(
                    PixPos(pos_x, pos_y), Flux(0.4 * obj['flux']),
                    EllipseE(obj['shapeexp_r'], obj['shapeexp_e1'],
                             obj['shapeexp_e2']), Flux(0.6 * obj['flux']),
                    EllipseE(obj['shapedev_r'], obj['shapedev_e1'],
                             obj['shapedev_e2'])))
        for obj in dev_galaxy:
            pos_x, pos_y = obj['x'], obj['y']
            sources.append(
                DevGalaxy(
                    PixPos(pos_x, pos_y), Flux(obj['flux']),
                    EllipseE(obj['shapedev_r'], obj['shapedev_e1'],
                             -obj['shapedev_e2'])))
        for obj in exp_galaxy:
            pos_x, pos_y = obj['x'], obj['y']
            sources.append(
                ExpGalaxy(
                    PixPos(pos_x, pos_y), Flux(obj['flux']),
                    EllipseE(obj['shapeexp_r'], obj['shapeexp_e1'],
                             -obj['shapeexp_e2'])))
        for obj in rex_galaxy:

if all_ptsrcs:
            cat.append(PointSource(pos, flux))
            continue

        tt = t.type.strip()
        if tt in ['PTSRC', 'STAR', 'S']:
            cat.append(PointSource(pos, flux))
        elif tt in ['EXP', 'E']:
            shape = EllipseE(*t.shapeexp)
            cat.append(ExpGalaxy(pos, flux, shape))
        elif tt in ['DEV', 'D']:
            shape = EllipseE(*t.shapedev)
            cat.append(DevGalaxy(pos, flux, shape))
        elif tt in ['COMP', 'C']:
            eshape = EllipseE(*t.shapeexp)
            dshape = EllipseE(*t.shapedev)
            cat.append(FixedCompositeGalaxy(pos, flux, t.fracdev,
                                            eshape, dshape))
        else:
            print('Did not understand row', i, 'of input catalog:')
            t.about()
            assert(False)

    W = unwise_forcedphot(cat, tiles, roiradecbox=roiradecbox,
                          bands=opt.bands, unwise_dir=opt.unwise_dir,
                          use_ceres=opt.ceres, ceres_block=opt.ceresblock,
                          save_fits=opt.save_fits, ps=ps)
    W.writeto(outfn)

def fromCovariance(cov):
        return EllipseESoft.fromEllipseE(EllipseE.fromCovariance(cov))

T = e.getRaDecBasis()
            # print 'T', T
            txy = np.dot(T, xy)
            # print 'txy', txy.shape
            plt.plot(e1 + txy[0, :], e2 + txy[1, :], '-', color=cc, alpha=0.5)
    plt.xlabel('ee1')
    plt.ylabel('ee2')
    plt.axis('scaled')
    plt.title('EllipseESoft')
    ps.savefig()

    plt.clf()
    for re, cc in zip([np.exp(4.), np.exp(5.), np.exp(6.)], 'rgb'):
        for e1, e2 in zip(E1.ravel(), E2.ravel()):
            e = EllipseE(re, e1, e2)
            print(e)
            T = e.getRaDecBasis()
            # print 'T', T
            txy = np.dot(T, xy)
            # print 'txy', txy.shape
            plt.plot(e1 + txy[0, :], e2 + txy[1, :], '-', color=cc, alpha=0.5)
    plt.xlabel('e1')
    plt.ylabel('e2')
    plt.axis('scaled')
    plt.title('EllipseE')
    ps.savefig()

    W, H = 500, 500
    img = np.zeros((H, W), np.float32)
    sig1 = 1.
    pixscale = 1.

for obj in exp_galaxy:
            pos_x, pos_y = obj['x'], obj['y']
            sources.append(
                ExpGalaxy(
                    PixPos(pos_x, pos_y), Flux(obj['flux']),
                    EllipseE(obj['shapeexp_r'], obj['shapeexp_e1'],
                             -obj['shapeexp_e2'])))
        for obj in rex_galaxy:
            #if obj['point_source'] > 0.0:
            #            sources.append(PointSource(PixPos(w.wcs_world2pix([[obj['ra'], obj['dec']]],1)[0]),
            #                                               Flux(obj['flux'])))
            pos_x, pos_y = obj['x'], obj['y']
            sources.append(
                ExpGalaxy(
                    PixPos(pos_x, pos_y), Flux(obj['flux']),
                    EllipseE(obj['shapeexp_r'], obj['shapeexp_e1'],
                             -obj['shapeexp_e2'])))

        for obj in psf_galaxy:
            pos_x, pos_y = obj['x'], obj['y']
            sources.append(PointSource(PixPos(pos_x, pos_y), Flux(obj['flux'])))

        print("Now you have %d sources" % len(sources))
    else:
         raise ValueError('Cannot use this shape method') 
    return sources

cat = []
    for i, t in enumerate(T):
        pos = RaDecPos(t.ra, t.dec)
        flux = NanoMaggies(**{wanyband: 1.})
        if all_ptsrcs:
            cat.append(PointSource(pos, flux))
            continue

        tt = t.type.strip()
        if tt in ['PTSRC', 'STAR', 'S']:
            cat.append(PointSource(pos, flux))
        elif tt in ['EXP', 'E']:
            shape = EllipseE(*t.shapeexp)
            cat.append(ExpGalaxy(pos, flux, shape))
        elif tt in ['DEV', 'D']:
            shape = EllipseE(*t.shapedev)
            cat.append(DevGalaxy(pos, flux, shape))
        elif tt in ['COMP', 'C']:
            eshape = EllipseE(*t.shapeexp)
            dshape = EllipseE(*t.shapedev)
            cat.append(FixedCompositeGalaxy(pos, flux, t.fracdev,
                                            eshape, dshape))
        else:
            print('Did not understand row', i, 'of input catalog:')
            t.about()
            assert(False)

    W = unwise_forcedphot(cat, tiles, roiradecbox=roiradecbox,
                          bands=opt.bands, unwise_dir=opt.unwise_dir,
                          use_ceres=opt.ceres, ceres_block=opt.ceresblock,
                          save_fits=opt.save_fits, ps=ps)
    W.writeto(outfn)

def fromRAbPhi(r, ba, phi):
        ab = 1. / ba
        e = (ab - 1) / (ab + 1)
        angle = math.radians(2. * (-phi))
        e1 = e * math.cos(angle)
        e2 = e * math.sin(angle)
        return EllipseE(r, e1, e2)