Commit 93270bbf authored by Zhang Xin's avatar Zhang Xin
Browse files

Merge branch 'sim_scheduler' of...

Merge branch 'sim_scheduler' of https://csst-tb.bao.ac.cn/code/csst_sim/csst-simulation into sim_scheduler
parents af3e3304 0d1520b4
Loading
Loading
Loading
Loading

tests/PSFInterpTest/readme

deleted100644 → 0
+0 −7
Original line number Original line Diff line number Diff line
unittest of PSF module:
  --test_loadPSFSet.py
    --testdata: S30x30_psfCube_1.h5 config_test.yaml
  --test_PSFInterpModule_coverage.py
    --testdata: S30x30_psfCube_1.h5 S20x20_psfCube_1.h5 config_test.yaml
  --test_Convolve.py
    --testdata: ccd1-w1
+0 −176
Original line number Original line Diff line number Diff line
'''
test galsim.interpolatedImage & galsim.convolve
'''

import os
import unittest

import numpy as np
import matplotlib.pyplot as plt
import galsim
import scipy.io
from scipy import ndimage


def psfEncircle(img, fraction=0.8, psfSampleSizeInMicrons=2.5, focalLengthInMeters=28, cenPix=None):
    #imgMaxPix_x, imgMaxPix_y = findMaxPix(img)
    y,x = ndimage.center_of_mass(img)  #y-rows, x-cols
    imgMaxPix_x = x #int(x)
    imgMaxPix_y = y #int(y)
    if cenPix != None:
        imgMaxPix_x = cenPix[0]
        imgMaxPix_y = cenPix[1]

    im1 = img.copy()
    im1size = im1.shape

    dis = np.zeros_like(img)
    for irow in range(im1size[0]):
        for icol in range(im1size[1]):
            dx = icol - imgMaxPix_x
            dy = irow - imgMaxPix_y
            dis[irow, icol] = np.hypot(dx, dy)

    nn = im1size[1]*im1size[0]
    disX = dis.reshape(nn)
    disXsortId = np.argsort(disX)

    imgX = img.reshape(nn)
    imgY = imgX[disXsortId]
    psfFrac = np.cumsum(imgY)/np.sum(imgY)
    ind = np.where(psfFrac > fraction)[0][0]

    REE80 = np.rad2deg(dis[np.where(img == imgY[ind])]*psfSampleSizeInMicrons*1e-6/focalLengthInMeters)*3600
    return REE80



def check_galsimConvolve(path=None, plotImage=True):
    #load psf data
    data=scipy.io.loadmat(path)
    imPSF = data['psf']
    pixSize = np.rad2deg(5.*1e-6/28)*3600
    imPSF = imPSF/np.sum(imPSF)

    #psf -> galsimInterpolatedImage
    img = galsim.ImageF(imPSF, scale=pixSize)
    imgt= galsim.InterpolatedImage(img)

    #imPSFt = imgt.drawImage(nx=256, ny=256, scale=pixSize, method='no_pixel')
    imPSFt = imgt.drawImage(nx=256, ny=256, scale=pixSize)

    ree80 = psfEncircle(imPSF, fraction=0.8, psfSampleSizeInMicrons=5.)
    ree80_pix = ree80/(np.rad2deg((5.*1e-6/28))*3600)

    sliceX = slice(128-int(np.round(ree80_pix[0])), 128+int(np.round(ree80_pix[0]))+1, 1)

    #set a point sorce
    src    = galsim.DeltaFunction(flux=1.0)

    result = galsim.Convolve(src, imgt)

    #drawImage with same pixSize
    #tmp    = result.drawImage(nx=256, ny=256, scale=pixSize, method='no_pixel')
    tmp    = result.drawImage(nx=256, ny=256, scale=pixSize)


    if plotImage != True:
        res = (imPSFt.array - imPSF)/imPSF
        d0 = np.mean(res[sliceX, sliceX].flatten())

        res = (tmp.array - imPSFt.array)/imPSFt.array
        d1 = np.mean(res[sliceX, sliceX].flatten())
        return d0, d1


    #plot images
    fig = plt.figure(figsize=(22, 5))
    ax=plt.subplot(1,3,1)
    plt.imshow(imPSF[128-10:128+10, 128-10:128+10])
    plt.annotate("ORG", [0.1, 0.9],xycoords="axes fraction", fontsize=16, color="w")
    plt.colorbar()

    ax=plt.subplot(1,3,2)
    plt.imshow(imPSFt.array[128-10:128+10, 128-10:128+10])
    plt.annotate("InterpolatedImage", [0.1, 0.9],xycoords="axes fraction", fontsize=16, color="w")
    plt.colorbar()

    ax=plt.subplot(1,3,3)
    plt.imshow(tmp.array[128-10:128+10, 128-10:128+10])
    plt.annotate("ConvolvedImage", [0.1, 0.9],xycoords="axes fraction", fontsize=16, color="w")
    plt.colorbar()
    plt.savefig(OUTPUTPATH+'/fig_check_galsimConvolve_1.pdf')

    fig = plt.figure(figsize=(13, 10))
    ax=plt.subplot(2,2,1)
    res = (imPSFt.array - imPSF)/imPSF
    plt.imshow(res[128-10:128+10, 128-10:128+10])
    plt.annotate("$\Delta_1$", [0.1, 0.9],xycoords="axes fraction", fontsize=16, color="w")
    plt.colorbar()

    ax=plt.subplot(2,2,2)
    plt.hist(res[sliceX, sliceX].flatten(), alpha=0.75, bins=4)
    #plt.annotate("$\Delta_1^{\\rm REE80}$", [0.1, 0.9],xycoords="axes fraction", fontsize=16, color="k")
    plt.xlabel("$\Delta_1^{\\rm REE80}$", fontsize=16)
    plt.ylabel("PDF", fontsize=16)

    ax=plt.subplot(2,2,3)
    res = (tmp.array - imPSFt.array)/imPSFt.array
    plt.imshow(res[128-10:128+10, 128-10:128+10])
    plt.annotate("$\Delta_2$", [0.1, 0.9],xycoords="axes fraction", fontsize=16, color="w")
    plt.colorbar()

    ax=plt.subplot(2,2,4)
    plt.hist(res[sliceX, sliceX].flatten(), alpha=0.75, bins=4)
    #plt.annotate("$\Delta_2^{\\rm REE80}$", [0.1, 0.9],xycoords="axes fraction", fontsize=16, color="k")
    plt.xlabel("$\Delta_2^{\\rm REE80}$", fontsize=16)
    plt.ylabel("PDF", fontsize=16)

    plt.savefig(OUTPUTPATH+'/fig_check_galsimConvolve_2.pdf')


def check_galsimConvolveALL(dataPath):
    d0 = np.zeros(900)
    d1 = np.zeros(900)
    for ipsf in range(1,901,1):
        print("ipsf={:}".format(ipsf), end="\r")
        psfPath = dataPath+"/ccd1-w1/psf_{:}_centroidWgt_BC.mat".format(ipsf)
        t0,t1=check_galsimConvolve(path = psfPath, plotImage=False)
        d0[ipsf-1] = t0
        d1[ipsf-1] = t1
      
    fig = plt.figure(figsize=(12,6))
    ax = plt.subplot(1,2,1)
    #plt.scatter(np.linspace(1,900,900), d0)
    plt.hist(d0, bins=8, alpha=0.75)
    plt.xlabel("mean($\Delta_1^{\\rm REE80}$)", fontsize=16)
    plt.ylabel("PDF", fontsize=16)

    ax = plt.subplot(1,2,2)
    #plt.scatter(np.linspace(1,900,900), d1)
    plt.hist(d1, bins=8, alpha=0.75)
    plt.xlabel("mean($\Delta_2^{\\rm REE80}$)", fontsize=16)
    plt.ylabel("PDF", fontsize=16)

    plt.savefig(OUTPUTPATH+'/fig_check_galsimConvolveALL.pdf')
 

class testConvolve(unittest.TestCase):
    def __init__(self, methodName='runTest'):
        super(testConvolve,self).__init__(methodName)
        self.dataPath = os.path.join(os.getenv('UNIT_TEST_DATA_ROOT'), 'csst_fz_gc1')

        global OUTPUTPATH
        OUTPUTPATH = os.path.join(self.dataPath, 'outputs')

    def test_galsimConvolve(self):
        ipsf = 1
        psfPath = self.dataPath+"/ccd1-w1/psf_{:}_centroidWgt_BC.mat".format(ipsf)
        check_galsimConvolve(path = psfPath)

    def test_galsimConvolveALL(self):
        check_galsimConvolveALL(dataPath=self.dataPath)


if __name__ == "__main__":
    unittest.main()
+0 −403
Original line number Original line Diff line number Diff line
import unittest

import sys,os,math
from itertools import islice

import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
mpl.use('Agg')

import yaml
from ObservationSim.Instrument import Chip
from ObservationSim.PSF.PSFInterp import PSFInterp


def defineCCD(iccd, config_file):
    with open(config_file, "r") as stream:
        try:
            config = yaml.safe_load(stream)
            #for key, value in config.items():
            #    print (key + " : " + str(value))
        except yaml.YAMLError as exc:
            print(exc)
    chip = Chip(chipID=iccd, config=config)
    #chip = Chip(chipID=iccd, ccdEffCurve_dir=path_dict["ccd_dir"], CRdata_dir=path_dict["CRdata_dir"], normalize_dir=path_dict["normalize_dir"], sls_dir=path_dict['sls_dir'], config=config)
    return chip

def psfSecondMoments(psfMat, cenX, cenY, pixSize=1):
    apr = 0.5 #arcsec, 0.5角秒内测量
    fl  = 28. #meters
    pxs = 2.5 #microns
    apr = np.deg2rad(apr/3600.)*fl*1e6
    apr = apr/pxs
    apr = int(np.ceil(apr))

    I = psfMat
    ncol = I.shape[1]
    nrow = I.shape[0]
    w   = 0.0
    w11 = 0.0
    w12 = 0.0
    w22 = 0.0
    for icol in range(ncol):
        for jrow in range(nrow):
            x = icol*pixSize - cenX
            y = jrow*pixSize - cenY
            rr = np.sqrt(x*x + y*y)
            wgt= 0.0
            if rr <= apr:
                wgt = 1.0
            w   += I[jrow, icol]*wgt
            w11 += x*x*I[jrow, icol]*wgt
            w12 += x*y*I[jrow, icol]*wgt
            w22 += y*y*I[jrow, icol]*wgt
    w11 /= w
    w12 /= w
    w22 /= w
    sz = w11 + w22
    e1 = (w11 - w22)/sz
    e2 = 2.0*w12/sz

    return sz, e1, e2


def test_psfEll(iccd, iwave, psfMat):
    psfMat_iwave = psfMat.psfMat[iwave-1, :,:,:]
    npsf = np.shape(psfMat_iwave)[0]

    imx = np.zeros(npsf)
    imy = np.zeros(npsf)
    psf_e1 = np.zeros(npsf)
    psf_e2 = np.zeros(npsf)
    psf_sz = np.zeros(npsf)

    for ipsf in range(1, npsf+1):
        print('ipsf-{:}'.format(ipsf), end='\r')
        imx[ipsf-1] = psfMat.cen_col[iwave-1, ipsf-1]
        imy[ipsf-1] = psfMat.cen_row[iwave-1, ipsf-1]
        psfMat_iwave_ipsf = psfMat_iwave[ipsf-1, :, :]
        cenX = 256
        cenY = 256
        sz, e1, e2  = psfSecondMoments(psfMat_iwave_ipsf, cenX, cenY, pixSize=1)
        psf_e1[ipsf-1] = e1
        psf_e2[ipsf-1] = e2
        psf_sz[ipsf-1] = sz
        #print('ell======', ipsf, np.sqrt(e1**2 + e2**2))
    #######
    arr = [imx, imy, psf_e1, psf_e2, psf_sz]
    np.save(OUTPUTPATH+'/psfEll{:}_{:}_{:}'.format(int(np.sqrt(npsf)),iccd, iwave), arr)


def test_psfEllPlot(OVERPLOT=False):
    #if ThisTask == 0:
    if True:
        prefix = 'psfEll30'
        iccd = 1
        iwave= 1
        data = np.load(OUTPUTPATH+'/'+prefix+'_1_1.npy')
        imx= data[0]
        imy= data[1]
        psf_e1 = data[2]
        psf_e2 = data[3]
        print(np.shape(imx))
        npsf = np.shape(imx)[0]

        #######
        plt.cla()
        plt.close("all")

        fig = plt.figure(figsize=(12, 12))
        plt.subplots_adjust(wspace=0.1, hspace=0.1)
        ax = plt.subplot(1, 1, 1)
        for ipsf in range(npsf):
            plt.plot(imx[ipsf], imy[ipsf], 'r.')

            ang = np.arctan2(psf_e2[ipsf], psf_e1[ipsf])/2
            ell = np.sqrt(psf_e1[ipsf]**2 + psf_e2[ipsf]**2)
            ell *= 15
            lcos = ell*np.cos(ang)
            lsin = ell*np.sin(ang)
            plt.plot([imx[ipsf]-lcos, imx[ipsf]+lcos],[imy[ipsf]-lsin, imy[ipsf]+lsin],'r', lw=2)
            ###########
        ang = 0.
        ell = 0.05
        ell*= 15
        lcos = ell*np.cos(ang)
        lsin = ell*np.sin(ang)
        plt.plot([imx[898]-lcos, imx[898]+lcos],[imy[898]+5.-lsin, imy[898]+5.+lsin],'k', lw=2)
        plt.annotate('{:}'.format(ell/15), (imx[898]-2., imy[898]+6.), xycoords='data', fontsize=10)
        plt.xlabel('CCD X (mm)')
        plt.ylabel('CCD Y (mm)')

        if OVERPLOT == True:
            prefix = 'psfEll20'
            data = np.load(OUTPUTPATH+'/'+prefix+'_1_1.npy')
            imx= data[0]
            imy= data[1]
            psf_e1 = data[2]
            psf_e2 = data[3]

            npsf = np.shape(imx)[0]
            for ipsf in range(npsf):
                plt.plot(imx[ipsf], imy[ipsf], 'b.')

                ang = np.arctan2(psf_e2[ipsf], psf_e1[ipsf])/2
                ell = np.sqrt(psf_e1[ipsf]**2 + psf_e2[ipsf]**2)
                ell *= 15
                lcos = ell*np.cos(ang)
                lsin = ell*np.sin(ang)
                plt.plot([imx[ipsf]-lcos, imx[ipsf]+lcos],[imy[ipsf]-lsin, imy[ipsf]+lsin],'b', lw=2)

        plt.gca().set_aspect(1)
        if OVERPLOT == True:
            prefix = 'psfEllOP'
        plt.savefig(OUTPUTPATH+'/'+prefix+'_iccd{:}.pdf'.format(iccd))


def test_psfIDW(iccd, iwave, psfMatA, chip, psfMatB):
    bandpass = iwave-1
    class pos_img():
        def __init__(self,x, y):
            self.x = x*1e3/10.  #in unit of pixels
            self.y = y*1e3/10.

    psfMat_iwave = psfMatA.psfMat[iwave-1, :,:,:]
    npsf = np.shape(psfMat_iwave)[0]
    psf_e1 = np.zeros(npsf)
    psf_e2 = np.zeros(npsf)
    psf_sz = np.zeros(npsf)

    for ipsf in range(1, npsf+1):
        print('ipsf:', ipsf, end='\r', flush=True)
        tpos_img = pos_img(psfMatA.cen_col[iwave-1, ipsf-1], psfMatA.cen_row[iwave-1, ipsf-1])
        psfIDW = psfMatB.get_PSF(chip, tpos_img, bandpass, galsimGSObject=False, findNeighMode='treeFind')
        np.save(OUTPUTPATH+'/psfIDW_{:}_{:}_{:}'.format(iccd, iwave, ipsf), psfIDW)

        cenX = 256
        cenY = 256
        sz, e1, e2  = psfSecondMoments(psfIDW, cenX, cenY, pixSize=1)
        psf_e1[ipsf-1] = e1
        psf_e2[ipsf-1] = e2
        psf_sz[ipsf-1] = sz
    arr = [psf_e1, psf_e2, psf_sz]
    np.save(OUTPUTPATH+'/psfEll20IDW_{:}_{:}'.format(iccd, iwave), arr)


def test_psfResidualPlot(iccd, iwave, ipsf, psfMatA):
    psfMat_iwave = psfMatA.psfMat[iwave-1, :,:,:]

    psfMatORG = psfMat_iwave[ipsf-1, :, :]
    psfMatIDW = np.load(OUTPUTPATH+'/psfIDW_{:}_{:}_{:}.npy'.format(iccd, iwave, ipsf))

    npix = psfMatORG.shape[0]
    pixCutEdge= int(npix/2-15)

    img0 = psfMatORG[pixCutEdge:npix-pixCutEdge, pixCutEdge:npix-pixCutEdge]
    img1 = psfMatIDW[pixCutEdge:npix-pixCutEdge, pixCutEdge:npix-pixCutEdge]
    imgX = (img1 - img0)/img0

    img0 = np.log10(img0)
    img1 = np.log10(img1)
    imgX = np.log10(np.abs(imgX))

    fig = plt.figure(figsize=(18,4))
    ax = plt.subplot(1,3,1)
    plt.imshow(img0, origin='lower', vmin=-7, vmax=-1.3)
    plt.plot([npix/2-pixCutEdge, npix/2-pixCutEdge],[0, (npix/2-pixCutEdge)*2-1],'w--')
    plt.plot([0, (npix/2-pixCutEdge)*2-1],[npix/2-pixCutEdge, npix/2-pixCutEdge],'w--')
    plt.annotate('ORG', [0,(npix/2-pixCutEdge)*2-5], c='w', size=15)
    cticks=[-7, -6, -5, -4, -3, -2, -1]
    cbar = plt.colorbar(ticks=cticks)
    cbar.ax.set_yticklabels(['$10^{-7}$', '$10^{-6}$', '$10^{-5}$','$10^{-4}$','$10^{-3}$','$10^{-2}$', '$10^{-1}$'])
    print(img0.min(), img0.max())


    ax = plt.subplot(1,3,2)
    plt.imshow(img1, origin='lower', vmin=-7, vmax=-1.3)
    plt.plot([npix/2-pixCutEdge, npix/2-pixCutEdge],[0, (npix/2-pixCutEdge)*2-1],'w--')
    plt.plot([0, (npix/2-pixCutEdge)*2-1],[npix/2-pixCutEdge, npix/2-pixCutEdge],'w--')
    plt.annotate('IDW', [0,(npix/2-pixCutEdge)*2-5], c='w', size=15)
    cticks=[-7, -6, -5, -4, -3, -2, -1]
    cbar = plt.colorbar(ticks=cticks)
    cbar.ax.set_yticklabels(['$10^{-7}$', '$10^{-6}$', '$10^{-5}$','$10^{-4}$','$10^{-3}$','$10^{-2}$', '$10^{-1}$'])
    print(img1.min(), img1.max())


    ax = plt.subplot(1,3,3)
    plt.imshow(imgX, origin='lower', vmin =-3, vmax =np.log10(3e-1))
    plt.plot([npix/2-pixCutEdge, npix/2-pixCutEdge],[0, (npix/2-pixCutEdge)*2-1],'w--')
    plt.plot([0, (npix/2-pixCutEdge)*2-1],[npix/2-pixCutEdge, npix/2-pixCutEdge],'w--')
    #plt.annotate('(IDW-ORG)/ORG', [0,(npix/2-pixCutEdge)*2-5], c='w', size=15)
    cticks=[-5, -4, -3, -2, -1]
    cbar = plt.colorbar(ticks=cticks)
    cbar.ax.set_yticklabels(['$10^{-5}$','$10^{-4}$','$10^{-3}$','$10^{-2}$', '$10^{-1}$'])

    print(np.max((psfMatORG-psfMatIDW)))
    plt.savefig(OUTPUTPATH+'/psfResidual_iccd{:}.pdf'.format(iccd))


def test_psfEllIDWPlot(OVERPLOT=False):
    #if ThisTask == 0:
    if True:
        prefix = 'psfEll20'
        iccd = 1
        iwave= 1
        data = np.load(OUTPUTPATH+'/'+prefix+'_1_1.npy')
        imx= data[0]
        imy= data[1]
        psf_e1 = data[2]
        psf_e2 = data[3]
        print(np.shape(imx))
        npsf = np.shape(imx)[0]

        #######
        plt.cla()
        plt.close("all")

        fig = plt.figure(figsize=(12, 12))
        plt.subplots_adjust(wspace=0.1, hspace=0.1)
        ax = plt.subplot(1, 1, 1)
        for ipsf in range(npsf):
            plt.plot(imx[ipsf], imy[ipsf], 'b.')

            ang = np.arctan2(psf_e2[ipsf], psf_e1[ipsf])/2
            ell = np.sqrt(psf_e1[ipsf]**2 + psf_e2[ipsf]**2)
            ell *= 15
            lcos = ell*np.cos(ang)
            lsin = ell*np.sin(ang)
            plt.plot([imx[ipsf]-lcos, imx[ipsf]+lcos],[imy[ipsf]-lsin, imy[ipsf]+lsin],'b', lw=2)
            ###########
        ang = 0.
        ell = 0.05
        ell*= 15
        lcos = ell*np.cos(ang)
        lsin = ell*np.sin(ang)
        #plt.plot([imx[898]-lcos, imx[898]+lcos],[imy[898]+5.-lsin, imy[898]+5.+lsin],'k', lw=2)
        #plt.annotate('{:}'.format(ell/15), (imx[898]-2., imy[898]+6.), xycoords='data', fontsize=10)
        plt.xlabel('CCD X (mm)')
        plt.ylabel('CCD Y (mm)')

        if OVERPLOT == True:
            prefix = 'psfEll20IDW'
            data = np.load(OUTPUTPATH+'/'+prefix+'_1_1.npy')
            #imx= data[0]
            #imy= data[1]
            psf_e1 = data[0]
            psf_e2 = data[1]

            npsf = np.shape(imx)[0]
            for ipsf in range(npsf):
                #plt.plot(imx[ipsf], imy[ipsf], 'r.')

                ang = np.arctan2(psf_e2[ipsf], psf_e1[ipsf])/2
                ell = np.sqrt(psf_e1[ipsf]**2 + psf_e2[ipsf]**2)
                ell *= 15
                lcos = ell*np.cos(ang)
                lsin = ell*np.sin(ang)
                plt.plot([imx[ipsf]-lcos, imx[ipsf]+lcos],[imy[ipsf]-lsin, imy[ipsf]+lsin],'r', lw=1)

        plt.gca().set_aspect(1)
        if OVERPLOT == True:
            prefix = 'psfEllOPIDW'
        plt.savefig(OUTPUTPATH+'/'+prefix+'_iccd{:}.pdf'.format(iccd))


def test_psfdEllabsPlot(iccd):
    #if ThisTask == 0:
    if True:
        prefix = 'psfEll20'
        #iccd = 1
        #iwave= 1
        data = np.load(OUTPUTPATH+'/'+prefix+'_{:}_1.npy'.format(iccd))
        imx= data[0]
        imy= data[1]
        psf_e1 = data[2]
        psf_e2 = data[3]
        psf_sz = data[4]
        print(np.shape(imx))
        npsf = np.shape(imx)[0]

        ellX = np.sqrt(psf_e1**2 + psf_e2**2)
        angX = np.arctan2(psf_e2, psf_e1)/2
        angX = np.rad2deg(angX)
        szX  = psf_sz

        ##############################
        prefix = 'psfEll20IDW'
        data = np.load(OUTPUTPATH+'/'+prefix+'_{:}_1.npy'.format(iccd))
        #imx= data[0]
        #imy= data[1]
        psf_e1 = data[0]
        psf_e2 = data[1]
        psf_sz = data[2]

        ellY = np.sqrt(psf_e1**2 + psf_e2**2)
        angY = np.arctan2(psf_e2, psf_e1)/2
        angY = np.rad2deg(angY)
        szY  = psf_sz

        ##############################
        fig=plt.figure(figsize=(6, 5))
        grid = plt.GridSpec(3,1,left=0.15, right=0.95, wspace=None, hspace=0.02)
        #plt.subplots_adjust(left=None,bottom=None,right=None,top=None,wspace=None,hspace=0.02)
        ax = plt.subplot(grid[0:2,0])
        plt.plot([0.01,0.1],[0.01,0.1], 'k--', lw=1. )
        plt.scatter(ellX, ellY, color='b', alpha=1., s=3., edgecolors='None')
        #plt.xlim([0.015, 0.085])
        #plt.ylim([0.015, 0.085])
        plt.ylabel('$\epsilon_{\\rm IDW}$')
        plt.gca().axes.get_xaxis().set_visible(False)

        ax = plt.subplot(grid[2,0])
        plt.plot([0.015,0.085],[0.,0.], 'k--', lw=1. )
        plt.scatter(ellX, (ellY-ellX), color='b', s=3., edgecolors='None')
        #plt.xlim([0.015, 0.085])
        #plt.ylim([-0.0018, 0.0018])
        plt.xlabel('$\epsilon_{\\rm ORG}$')
        plt.ylabel('$\Delta$')
        plt.savefig(OUTPUTPATH+'/psfEllOPIDWPDF_{:}.pdf'.format(iccd))

        fig=plt.figure(figsize=(6, 6))
        plt.hist((szY-szX)/szX, bins=20, color='r', alpha=0.5)
        plt.xlabel('$(R_{\\rm IDW}-R_{\\rm ORG})/R_{\\rm ORG}$')
        plt.ylabel('PDF')
        plt.savefig(OUTPUTPATH+'/psfEllOPIDWPDF_dsz_{:}.pdf'.format(iccd))


class PSFInterpModule_coverage(unittest.TestCase):
    def __init__(self, methodName='runTest'):
        super(PSFInterpModule_coverage,self).__init__(methodName)
        self.dataPath = os.path.join(os.getenv('UNIT_TEST_DATA_ROOT'), 'csst_fz_gc1')
        global OUTPUTPATH
        OUTPUTPATH = os.path.join(self.dataPath, 'outputs')

    def test_psfEll_(self):
        iccd = 1
        iwave= 1

        config_file = os.path.join(self.dataPath, 'config_test.yaml') 
        chip = defineCCD(iccd, config_file)
        print(chip.chipID)
        print(chip.cen_pix_x, chip.cen_pix_y)

        psfMatA = PSFInterp(chip, npsf=400, PSF_data_file=self.dataPath, PSF_data_prefix="S20x20_")
        psfMatB = PSFInterp(chip, npsf=900, PSF_data_file=self.dataPath, PSF_data_prefix="S30x30_")

        test_psfEll(iccd, iwave, psfMatA)
        test_psfEll(iccd, iwave, psfMatB)
        test_psfEllPlot(OVERPLOT=True)

        test_psfIDW(iccd, iwave, psfMatA, chip, psfMatB)
        ipsf = 1
        test_psfResidualPlot(iccd, iwave, ipsf, psfMatA)
        test_psfEllIDWPlot(OVERPLOT=True)
        test_psfdEllabsPlot(iccd)


if __name__ == '__main__':
    unittest.main()
    print('#####haha#####')


+0 −62
Original line number Original line Diff line number Diff line
import unittest

import sys,os,math
from itertools import islice

import numpy as np

import yaml
from ObservationSim.Instrument import Chip
from ObservationSim.PSF.PSFInterp import PSFInterp


def defineCCD(iccd, config_file):
    with open(config_file, "r") as stream:
        try:
            config = yaml.safe_load(stream)
            #for key, value in config.items():
            #    print (key + " : " + str(value))
        except yaml.YAMLError as exc:
            print(exc)
    chip = Chip(chipID=iccd, config=config)
    #chip = Chip(chipID=iccd, ccdEffCurve_dir=path_dict["ccd_dir"], CRdata_dir=path_dict["CRdata_dir"], normalize_dir=path_dict["normalize_dir"], sls_dir=path_dict['sls_dir'], config=config)
    return chip

def loadPSFSet(iccd, dataPath):
    config_file = os.path.join(dataPath, 'config_test.yaml')
    chip = defineCCD(iccd, config_file)
    print(chip.chipID)
    print(chip.cen_pix_x, chip.cen_pix_y)

    psfMat= PSFInterp(chip, npsf=900, PSF_data_file=dataPath, PSF_data_prefix="S30x30_")
    psfSet= psfMat._loadPSF(iccd, dataPath, PSF_data_prefix="S30x30_")

    twave = 0 #[0...3]
    tpsf  = 0 #[0...899]
    field_x= psfSet[twave][tpsf]['field_x']
    field_y= psfSet[twave][tpsf]['field_y']
    image_x= psfSet[twave][tpsf]['image_x']
    image_y= psfSet[twave][tpsf]['image_y']
    centroid_x= psfSet[twave][tpsf]['centroid_x']
    centroid_y= psfSet[twave][tpsf]['centroid_y']
    print("pos_info:", field_x, field_y, image_x, image_y, centroid_x, centroid_y)
    return psfSet



class PSFInterpModule_coverage(unittest.TestCase):
    def __init__(self, methodName='runTest'):
        super(PSFInterpModule_coverage,self).__init__(methodName)
        self.dataPath = os.path.join(os.getenv('UNIT_TEST_DATA_ROOT'), 'csst_fz_gc1')

    def test_psfEll_(self):
        iccd = 1  #[1...30]
        psfSet = loadPSFSet(iccd, dataPath=self.dataPath)


if __name__ == '__main__':
    unittest.main()
    print('#####haha#####')


+0 −153

File deleted.

Preview size limit exceeded, changes collapsed.

Loading