Commit 2953a065 authored by Zhang Xin's avatar Zhang Xin
Browse files

Merge branch 'develop' into 'master'

version 1.0.5

See merge request csst_sim/csst-simulation!10
parents 243048a9 ca7fb5db
Loading
Loading
Loading
Loading
+2 −1
Original line number Diff line number Diff line
@@ -111,7 +111,8 @@ def on_orbit_obs_position(input_ra_list, input_dec_list, input_pmra_list, input_
        input_dec = c_double(input_dec_list[i])
        input_pmra = c_double(input_pmra_list[i])
        input_pmdec = c_double(input_pmdec_list[i])
        input_rv = c_double(input_rv_list[i] * 3600.) # Convert from km/s to km/h
        # input_rv = c_double(input_rv_list[i] * 3600.) # Convert from km/s to km/h
        input_rv = c_double(input_rv_list[i])
        input_parallax = c_double(input_parallax_list[i])
        p3 = d3(input_x, input_y, input_z)
        v3 = d3(input_vx, input_vy, input_vz)
+7 −2
Original line number Diff line number Diff line
@@ -68,8 +68,13 @@ class ChipOutput(object):
            self.cat.write(self.hdr)

    def cat_add_obj(self, obj, pos_img, pos_shear):
        ximg = pos_img.x - self.chip.bound.xmin + 1.0
        yimg = pos_img.y - self.chip.bound.ymin + 1.0
        # ximg = pos_img.x - self.chip.bound.xmin + 1.0
        # yimg = pos_img.y - self.chip.bound.ymin + 1.0
        # print('-------------debug-----------------')
        # print('from global',ximg, yimg)
        ximg = obj.real_pos.x + 1.0
        yimg = obj.real_pos.y + 1.0
        # print('from loacl',ximg, yimg)
        
        line = self.fmt%(
            obj.id, int(self.chipLabel), self.filt.filter_type, ximg, yimg, obj.ra, obj.dec, obj.ra_orig, obj.dec_orig, obj.z, obj.getMagFilter(self.filt), obj.type, 
+71 −3
Original line number Diff line number Diff line
@@ -7,7 +7,7 @@ from astropy.io import fits
import astropy.wcs as pywcs
from collections import OrderedDict

from scipy import math
# from scipy import math
import random

import os
@@ -87,6 +87,17 @@ def rotate_CD_matrix(cd, pa_aper):
    cd_rot = np.dot(mat, cd)
    return cd_rot

def calcaluteSLSRotSkyCoor(pix_xy = None,rot_angle = 1, xlen = 9216, ylen = 9232, w = None):
    rad = np.deg2rad(rot_angle)
    mat = np.zeros((2,2))
    mat[0,:] = np.array([np.cos(rad),-np.sin(rad)])
    mat[1,:] = np.array([np.sin(rad),np.cos(rad)])
    center = np.array([xlen/2, ylen/2])
    rot_pix = np.dot(mat, pix_xy-center) + center
    skyCoor = w.wcs_pix2world(np.array([rot_pix]), 1)

    return skyCoor


# def Header_extention(xlen = 9216, ylen = 9232, gain = 1.0, readout = 5.0, dark = 0.02,saturation=90000, row_num = 1, col_num = 1):
#
@@ -153,7 +164,7 @@ def rotate_CD_matrix(cd, pa_aper):


##9232 9216  898 534 1309 60 -40  -23.4333
def WCS_def(xlen = 9216, ylen = 9232, gapy = 898.0, gapx1 = 534, gapx2 = 1309, ra = 60, dec = -40, pa = -23.433,psize = 0.074, row_num = 1, col_num = 1):
def WCS_def(xlen = 9216, ylen = 9232, gapy = 898.0, gapx1 = 534, gapx2 = 1309, ra = 60, dec = -40, pa = -23.433,psize = 0.074, row_num = 1, col_num = 1, filter = 'GI'):

    """ Creat a wcs frame for CCST with multiple extensions

@@ -228,6 +239,63 @@ def WCS_def(xlen = 9216, ylen = 9232, gapy = 898.0, gapx1 = 534, gapx2 = 1309, r
        r_dat['CD1_2'] = cd_rot[0,1]
        r_dat['CD2_1'] = cd_rot[1,0]
        r_dat['CD2_2'] = cd_rot[1,1]

        if filter in ['GU', 'GV', 'GI']:
            from astropy import wcs

            w = wcs.WCS(naxis=2)
            w.wcs.crpix = [r_dat['CRPIX1'], r_dat['CRPIX2']]
            w.wcs.cd = cd_rot
            w.wcs.crval = [ra_ref, dec_ref]
            w.wcs.ctype = [r_dat['CTYPE1'], r_dat['CTYPE2']]

            # test_center_o = w.wcs_pix2world(np.array([[xlen / 2, ylen / 2]]), 1)

            sls_rot = 1
            if i > 2:
                sls_rot = -sls_rot

            sn_x = 30
            sn_y = 30
            x_pixs = np.zeros(sn_y * sn_x)
            y_pixs = np.zeros(sn_y * sn_x)
            xpixs_line = np.linspace(1, xlen, sn_x)
            ypixs_line = np.linspace(1, ylen, sn_y)

            sky_coors = []

            for n1, y in enumerate(ypixs_line):
                for n2, x in enumerate(xpixs_line):
                    i_pix = n1 * sn_x + n2
                    x_pixs[i_pix] = x
                    y_pixs[i_pix] = y

                    pix_coor = np.array([x, y])
                    sc1 = calcaluteSLSRotSkyCoor(pix_xy=pix_coor, rot_angle=sls_rot, w=w)
                    # print(sc1[0,0],sc1[0,1])
                    sky_coors.append((sc1[0, 0], sc1[0, 1]))

            from astropy.coordinates import SkyCoord
            from astropy.wcs.utils import fit_wcs_from_points

            wcs_new = fit_wcs_from_points(xy=np.array([x_pixs, y_pixs]),
                                          world_coords=SkyCoord(sky_coors, frame="icrs", unit="deg"), projection='TAN')

            # print(wcs_new)
            # test_center = wcs_new.wcs_pix2world(np.array([[xlen / 2, ylen / 2]]), 1)
            #
            # print(test_center - test_center_o)

            r_dat['CD1_1'] = wcs_new.wcs.cd[0, 0]
            r_dat['CD1_2'] = wcs_new.wcs.cd[0, 1]
            r_dat['CD2_1'] = wcs_new.wcs.cd[1, 0]
            r_dat['CD2_2'] = wcs_new.wcs.cd[1, 1]
            r_dat['CRPIX1'] = wcs_new.wcs.crpix[0]
            r_dat['CRPIX2'] = wcs_new.wcs.crpix[1]

            r_dat['CRVAL1'] = wcs_new.wcs.crval[0]
            r_dat['CRVAL2'] = wcs_new.wcs.crval[1]

        return r_dat


@@ -347,7 +415,7 @@ def generateExtensionHeader(xlen = 9216, ylen = 9232,ra = 60, dec = -40, pa = -2
    h_ext['CCDCHIP'] =  'ccd' + CCDID[k-1].rjust(2,'0')
    h_ext['POS_ANG'] = pa
    header_wcs = WCS_def(xlen=xlen, ylen=ylen, gapy=898.0, gapx1=534, gapx2=1309, ra=ra, dec=dec, pa=pa, psize=psize,
                         row_num=row_num, col_num=col_num)
                         row_num=row_num, col_num=col_num, filter = h_ext['FILTER'])

    h_ext['CRPIX1'] = header_wcs['CRPIX1']
    h_ext['CRPIX2'] = header_wcs['CRPIX2']
+20 −1
Original line number Diff line number Diff line
@@ -82,6 +82,8 @@ class Chip(FocalPlane):
        else:
            self.sensor = galsim.Sensor()

        self.flat_cube = None # for spectroscopic flat field cube simulation

    # def _getChipRowCol(self):
    #     self.rowID = (self.chipID - 1) // 5 + 1
    #     self.colID = (self.chipID - 1) % 5 + 1
@@ -719,7 +721,7 @@ class Chip(FocalPlane):
            if config["ins_effects"]["add_bias"] == True:
                biaslevel = self.bias_level
                overscan = biaslevel-2
            elif config["ins_effects"]["add_bias"] == True:
            elif config["ins_effects"]["add_bias"] == False:
                biaslevel = 0
                overscan = 0
            for i in range(NDark):
@@ -846,3 +848,20 @@ class Chip(FocalPlane):
        #             sub_img.write("%s/%s_%s.fits" % (chip_output.subdir, img_name_root, rowcoltag))
        #     del sub_img
        return img

    def loadSLSFLATCUBE(self, flat_fn='flat_cube.fits'):
        from astropy.io import fits
        try:
            with pkg_resources.files('ObservationSim.Instrument.data').joinpath(flat_fn) as data_path:
                flat_fits = fits.open(data_path, ignore_missing_simple=True)
        except AttributeError:
            with pkg_resources.path('ObservationSim.Instrument.data', flat_fn) as data_path:
                flat_fits = fits.open(data_path, ignore_missing_simple=True)

        fl = len(flat_fits)
        fl_sh = flat_fits[0].data.shape
        assert fl == 4, 'FLAT Field Cube is Not 4 layess!!!!!!!'
        self.flat_cube = np.zeros([fl, fl_sh[0], fl_sh[1]])
        for i in np.arange(0, fl, 1):
            self.flat_cube[i] = flat_fits[i].data
+2.54 GiB

File added.

No diff preview for this file type.

Loading