diff --git a/csst_ifs_sim/csst_ifs_sim.py b/csst_ifs_sim/csst_ifs_sim.py index bb53d643a83475f8b291c8aabc939789a88c64a8..7ba2cd5c7ee2d64990026c1d421fcd00c3cf5269 100644 --- a/csst_ifs_sim/csst_ifs_sim.py +++ b/csst_ifs_sim/csst_ifs_sim.py @@ -28,7 +28,6 @@ import astropy.coordinates as coord import ctypes import sys - ##sys.path.append('./csst_ifs_sim') conf.auto_max_age = None @@ -1149,6 +1148,25 @@ def beta_angle(x_sat, y_sat, z_sat, vx_sat, vy_sat, vz_sat, ra_obj, dec_obj): return angle ############################################################################### +def find_min(arr): + min_val = arr[0] + min_index = 0 + for i in range(1, len(arr)): + if arr[i] < min_val: + min_val = arr[i] + min_index = i + return min_val, min_index + +################################# +def find_max(arr): + max_val = arr[0] + max_index = 0 + for i in range(1, len(arr)): + if arr[i] > max_val: + max_val = arr[i] + max_index = i + return max_val, max_index + ########################################################## @@ -1534,12 +1552,16 @@ def get_dx_dy_blue(wave): """ # wave is the wavelength in nm; # dx is in dispersion direction, dy is in vertical direction; - dydl = np.array([-423.256, 0.001, 0.00075, - 0.0000078, -0.0000000000007, 0.0, 0.0]) - # 色散方向 - dxdl = 0.2*np.array([-9.1519, -1.00000000e-06, 3.50000000e-08, -5.00000000e-09, - -1.70000000e-11, 4.00949787e-12, -6.16873452e-15]) - # 垂直方向 + # dydl = np.array([-423.256, 0.001, 0.00075, + # 0.0000078, -0.0000000000007, 0.0, 0.0]) + # # 色散方向 + # dxdl = 0.2*np.array([-9.1519, -1.00000000e-06, 3.50000000e-08, -5.00000000e-09, + # -1.70000000e-11, 4.00949787e-12, -6.16873452e-15]) + + ##### update @2024.10.16 + dydl=np.array([ 2447.9, -1/0.141, 0.0000075, 0.00000078, -0.000000000007] ) ; #色散方向 + dxdl=np.array([ 5.46, -1.5e-02, 3.5e-08, -5.0e-09]) #垂直方向 + dx = 0.0 dy = 0.0 for i in range(len(dxdl)): @@ -1567,11 +1589,14 @@ def get_dx_dy_red(wave): """ # wave is the wavelength in nm; - dydl = np.array([-640.0239901372472, 0.0018, 0.00048, - 0.0000028, -0.0000000000007, 0.0, 0.0]) # 色散方向 - dxdl = 0.00325*np.array([-1638.8, 4.0e-2, 5.500e-3, - - 5.2e-10, 1.7000e-10, 7.1e-13, -5.16e-15]) # 垂直方向 - + # dydl = np.array([-640.0239901372472, 0.0018, 0.00048, + # 0.0000028, -0.0000000000007, 0.0, 0.0]) # 色散方向 + # dxdl = 0.00325*np.array([-1638.8, 4.0e-2, 5.500e-3, - + # 5.2e-10, 1.7000e-10, 7.1e-13, -5.16e-15]) # 垂直方向 + ## update @2014.10.17 + dydl=np.array([3519.78622, -1/0.1555, 0.0000048, 0.00000028, -0.0000000000007] ) #色散方向 + dxdl=np.array([-5.6305, 1.0e-2, 5.500e-7, -5.2e-10]) #垂直方向 + dx = 0.0 dy = 0.0 for i in range(len(dxdl)): @@ -1753,6 +1778,9 @@ class IFSsimulator(): self.readoutNoise = self.config.getboolean( self.section, 'readoutnoise') + + self.appbianpai = self.config.getboolean( + self.section, 'appbianpai') self.intscale = True ###################################################################### @@ -1769,7 +1797,8 @@ class IFSsimulator(): bleeding=self.bleeding, sky_noise=self.sky_noise, intscale=self.intscale, - save_cosmicrays=self.save_cosmicrays) + save_cosmicrays=self.save_cosmicrays, + appbianpai=self.appbianpai) ############################################################################ @@ -2134,9 +2163,14 @@ class IFSsimulator(): slice_red['py'][i] = 50+250+randRedpos[i]*4 slice_red['px'][i] = 3.55/0.015*(i-16)+1190.0+118 ####### - - self.slice_blue = slice_blue - self.slice_red = slice_red + ###### flip the fringe up to down,down to up@2024.10.16 + self.slice_blue = dict() + self.slice_red = dict() + self.slice_blue['py'] =2000-slice_blue['py'] + self.slice_blue['px'] = slice_blue['px'] + + self.slice_red['py'] =3000-slice_red['py'] + self.slice_red['px'] =slice_red['px'] ####################################################################### maskSlice = dict() @@ -2223,7 +2257,7 @@ class IFSsimulator(): b_3 = np.sum(temp_b_3, axis=0) temp_b_2 = np.fliplr( - self.image_b[1024:1024+1024, 2048:2048+2048]) # left to right + self.image_b[1024:1024+1024, 2048:2048+2048]) # left to right b_2 = np.sum(temp_b_2, axis=0) b_1 = np.sum(self.image_b[1024:1024+1024, 0:2048], axis=0) @@ -2439,35 +2473,35 @@ class IFSsimulator(): # blue zone 4 self.image_b[0:1024, 0:2048] += self.information['exptime'] * \ - self.information['dark4_b'] + self.information['dark2_b'] ########## zone 1 ################# self.image_b[1024:2048, 0:2048] += self.information['exptime'] * \ - self.information['dark1_b'] + self.information['dark3_b'] ########## zone 3 ################### self.image_b[0:1024, 2048:4096] += self.information['exptime'] * \ - self.information['dark3_b'] + self.information['dark1_b'] # zone 2 self.image_b[1024:2048, 2048:4096] += self.information['exptime'] * \ - self.information['dark2_b'] + self.information['dark4_b'] # red zone 4 self.image_r[0:1536, 0:3072] += self.information['exptime'] * \ - self.information['dark4_r'] + self.information['dark2_r'] ########## zone 1 ################# self.image_r[1536:3712, 0:3072] += self.information['exptime'] * \ - self.information['dark1_r'] + self.information['dark3_r'] ########## zone 3 ################### self.image_r[0:1536, 3072:6144] += self.information['exptime'] * \ - self.information['dark3_r'] + self.information['dark1_r'] # zone 2 self.image_r[1536:3072, 3072:6144] += self.information['exptime'] * \ - self.information['dark2_r'] + self.information['dark4_r'] ############################################################################## @@ -2715,44 +2749,44 @@ class IFSsimulator(): overscan = 320 # zone 4 , OSH #### 1024, 2048 self.image_b[0:1024+overscan, 0:2048+prescan+overscan] += np.random.normal( - loc=0.0, scale=self.information['rn4_b'], size=(1344, 2418)) + loc=0.0, scale=self.information['rn2_b'], size=(1344, 2418)) ########## zone 3, OSG ################# np.random.seed() self.image_b[0:1024+overscan, 2418:2418+2048+prescan+overscan] += np.random.normal( - loc=0.0, scale=self.information['rn3_b'], size=(1344, 2418)) + loc=0.0, scale=self.information['rn1_b'], size=(1344, 2418)) ########## zone 1, OSE ################### np.random.seed() self.image_b[0:1024+overscan, 2418*2:2418*2+2048+prescan+overscan] += np.random.normal( - loc=0.0, scale=self.information['rn1_b'], size=(1344, 2418)) + loc=0.0, scale=self.information['rn3_b'], size=(1344, 2418)) ########## zone 2, OSF ############### np.random.seed() self.image_b[0:1024+overscan, 2418*3:2418*3+2048+prescan+overscan] += np.random.normal( - loc=0.0, scale=self.information['rn2_b'], size=(1344, 2418)) + loc=0.0, scale=self.information['rn4_b'], size=(1344, 2418)) self.log.info('readnoise added in blue channel') # red zone 4 , OSH### 1536* 3072 np.random.seed() self.image_r[0:1536+overscan, 0:3072+prescan+overscan] += np.random.normal( - loc=0.0, scale=self.information['rn4_r'], size=(1856, 3442)) + loc=0.0, scale=self.information['rn2_r'], size=(1856, 3442)) ########## zone 3 ,OSG ################# np.random.seed() self.image_r[0:1536+overscan, 3442:3442+3072+prescan+overscan] += np.random.normal( - loc=0.0, scale=self.information['rn3_r'], size=(1856, 3442)) + loc=0.0, scale=self.information['rn1_r'], size=(1856, 3442)) ########## zone 1 ,OSE ################### np.random.seed() self.image_r[0:1536+overscan, 3442*2:3442*2+3072+prescan+overscan] += np.random.normal( - loc=0.0, scale=self.information['rn1_r'], size=(1856, 3442)) + loc=0.0, scale=self.information['rn3_r'], size=(1856, 3442)) ########## zone 2,OSF ########### np.random.seed() self.image_r[0:1536+overscan, 3442*3:3442*3+3072+prescan+overscan] += np.random.normal( - loc=0.0, scale=self.information['rn2_r'], size=(1856, 3442)) + loc=0.0, scale=self.information['rn4_r'], size=(1856, 3442)) ########################################################################################## @@ -2880,31 +2914,44 @@ class IFSsimulator(): self.log.info( 'Converting from electrons to ADUs using a factor of gain') - # blue zone 4 - self.image_b[0:1344, 0:2418] /= self.information['gain4_b'] + # # blue zone 4 + # self.image_b[0:1344, 0:2418] /= self.information['gain4_b'] - ########## zone 3 ################# - self.image_b[0:1344, 2418:2418*2] /= self.information['gain3_b'] + # ########## zone 3 ################# + # self.image_b[0:1344, 2418:2418*2] /= self.information['gain3_b'] - ########## zone 1 ################### - self.image_b[0:1344, 2418*2:2418*3] /= self.information['gain1_b'] + # ########## zone 1 ################### + # self.image_b[0:1344, 2418*2:2418*3] /= self.information['gain1_b'] - # zone 2 - self.image_b[0:1344, 2418*3:2418*4] /= self.information['gain2_b'] + # # zone 2 + # self.image_b[0:1344, 2418*3:2418*4] /= self.information['gain2_b'] + +################# update @2024.10.18 ### + # first part, menas old blue zone 4 + self.image_b[0:1344, 0:2418] /= self.information['gain1_b'] + ##########second part, means old zone 3 ################# + self.image_b[0:1344, 2418:2418*2] /= self.information['gain2_b'] + + ##########third part, means old zone 1 ################### + self.image_b[0:1344, 2418*2:2418*3] /= self.information['gain3_b'] + + #### fourth part, means old zone 2 + self.image_b[0:1344, 2418*3:2418*4] /= self.information['gain4_b'] ############################################################################ +########################## # red zone 4 - self.image_r[0:1856, 0:3442] /= self.information['gain4_r'] + self.image_r[0:1856, 0:3442] /= self.information['gain1_r'] ########## zone 3 ################# - self.image_r[0:1856, 3442:3442*2] /= self.information['gain3_r'] + self.image_r[0:1856, 3442:3442*2] /= self.information['gain2_r'] ########## zone 1 ################### - self.image_r[0:1856, 3442*2:3442*3] /= self.information['gain1_r'] + self.image_r[0:1856, 3442*2:3442*3] /= self.information['gain3_r'] # zone 2 - self.image_r[0:1856, 3442*3:3442*4] /= self.information['gain2_r'] + self.image_r[0:1856, 3442*3:3442*4] /= self.information['gain4_r'] # 3 @@ -2979,16 +3026,16 @@ class IFSsimulator(): ######################################################################## # blue zone 4 - self.image_b[0:1344, 0:2418] += self.information['bias4_b'] + self.image_b[0:1344, 0:2418] += self.information['bias2_b'] ########## zone 3 ################# - self.image_b[0:1344, 2418:2418*2] += self.information['bias3_b'] + self.image_b[0:1344, 2418:2418*2] += self.information['bias1_b'] ########## zone 1 ################### - self.image_b[0:1344, 2418*2:2418*3] += self.information['bias1_b'] + self.image_b[0:1344, 2418*2:2418*3] += self.information['bias3_b'] # zone 2 - self.image_b[0:1344, 2418*3:2418*4] += self.information['bias2_b'] + self.image_b[0:1344, 2418*3:2418*4] += self.information['bias4_b'] ############################################################################ @@ -2997,16 +3044,16 @@ class IFSsimulator(): ####################################################################### # red zone 4 - self.image_r[0:1856, 0:3442] += self.information['bias4_r'] + self.image_r[0:1856, 0:3442] += self.information['bias2_r'] ########## zone 3 ################# - self.image_r[0:1856, 3442:3442*2] += self.information['bias3_r'] + self.image_r[0:1856, 3442:3442*2] += self.information['bias1_r'] ########## zone 1 ################### - self.image_r[0:1856, 3442*2:3442*3] += self.information['bias1_r'] + self.image_r[0:1856, 3442*2:3442*3] += self.information['bias3_r'] # zone 2 - self.image_r[0:1856, 3442*3:3442*4] += self.information['bias2_r'] + self.image_r[0:1856, 3442*3:3442*4] += self.information['bias4_r'] ####################################################################### @@ -3290,14 +3337,46 @@ class IFSsimulator(): overscan = int(self.information['overscan']) temp = np.zeros((1344, 9672)) - # zone 4, OSH ######0:1024, 50:2048+50 + # # zone 4, OSH ######0:1024, 50:2048+50 + # x1 = 0 + # x2 = x1+1024 + + # y1 = 0+prescan + # y2 = y1+2048 + # temp[x1:x2, y1:y2] = imgb[0:1024, 0:2048] + # # zone 3, OSG , left to right ################# + # # np.fliplr(b2) ## left to right + # # np.flipud(b3) ## down to up + # x1 = 0 + # x2 = x1+1024 + + # y1 = 2418+prescan + # y2 = y1+2048 + # temp[x1:x2, y1:y2] = np.fliplr(imgb[0:1024, 2048:4096]) + # # zone 1, OSE,down to up ################### + # x1 = 0 + # x2 = x1+1024 + + # y1 = 2418*2+prescan + # y2 = y1+2048 + # temp[x1:x2, y1:y2] = np.flipud(imgb[1024:2048, 0:2048]) + # ########## zone 2, OSF down to yp ,left to right ####### + # x1 = 0 + # x2 = x1+1024 + + # y1 = 2418*3+prescan + # y2 = y1+2048 + # temp[x1:x2, y1:y2] = np.flipud(np.fliplr(imgb[1024:2048, 2048:4096])) + + ### update 2024.10.18 + # first part, old OSG part ,## shift: left to right x1 = 0 x2 = x1+1024 y1 = 0+prescan y2 = y1+2048 - temp[x1:x2, y1:y2] = imgb[0:1024, 0:2048] - # zone 3, OSG , left to right ################# + temp[x1:x2, y1:y2] = np.fliplr(imgb[0:1024, 2048:4096]) + # second part, old OSH, no shift ################# # np.fliplr(b2) ## left to right # np.flipud(b3) ## down to up x1 = 0 @@ -3305,21 +3384,22 @@ class IFSsimulator(): y1 = 2418+prescan y2 = y1+2048 - temp[x1:x2, y1:y2] = np.fliplr(imgb[0:1024, 2048:4096]) - # zone 1, OSE,down to up ################### + temp[x1:x2, y1:y2] = imgb[0:1024, 0:2048] + ### third part, old OSE, down to up ################### x1 = 0 x2 = x1+1024 y1 = 2418*2+prescan y2 = y1+2048 temp[x1:x2, y1:y2] = np.flipud(imgb[1024:2048, 0:2048]) - ########## zone 2, OSF down to yp ,left to right ####### + ########## fourth part, old OSF part; down to yp ,left to right ####### x1 = 0 x2 = x1+1024 y1 = 2418*3+prescan y2 = y1+2048 temp[x1:x2, y1:y2] = np.flipud(np.fliplr(imgb[1024:2048, 2048:4096])) + self.image_b = temp @@ -3327,14 +3407,46 @@ class IFSsimulator(): imgr = self.image_r.copy() temp = np.zeros((1856, 13768)) - # zone 4, OSH ######0:1024, 50:2048+50 + # # zone 4, OSH ######0:1024, 50:2048+50 + # x1 = 0 + # x2 = x1+1536 + + # y1 = 0+prescan + # y2 = y1+3072 + # temp[x1:x2, y1:y2] = imgr[0:1536, 0:3072] + # ########## zone 3, OSG , left to right ################# + # # np.fliplr(b2) ## left to right + # # np.flipud(b3) ## down to up + # x1 = 0 + # x2 = x1+1536 + + # y1 = 3442+prescan + # y2 = y1+3072 + # temp[x1:x2, y1:y2] = np.fliplr(imgr[0:1536, 3072:6144]) + # ########## zone 1, OSE,down to up ############################ + # x1 = 0 + # x2 = x1+1536 + + # y1 = 3442*2+prescan + # y2 = y1+3072 + # temp[x1:x2, y1:y2] = np.flipud(imgr[1536:3072, 0:3072]) + # ########## zone 2, OSF down to yp ,left to right ################ + # x1 = 0 + # x2 = x1+1536 + + # y1 = 3442*3+prescan + # y2 = y1+3072 + # temp[x1:x2, y1:y2] = np.flipud(np.fliplr(imgr[1536:3072, 3072:6144])) + + ###### update @2024.10.18 + # readout image ,first part, old OSG, shift: left to right x1 = 0 x2 = x1+1536 y1 = 0+prescan y2 = y1+3072 - temp[x1:x2, y1:y2] = imgr[0:1536, 0:3072] - ########## zone 3, OSG , left to right ################# + temp[x1:x2, y1:y2] = np.fliplr(imgr[0:1536, 3072:6144]) + ########## second part, old OSH ,no change; ################# # np.fliplr(b2) ## left to right # np.flipud(b3) ## down to up x1 = 0 @@ -3342,21 +3454,22 @@ class IFSsimulator(): y1 = 3442+prescan y2 = y1+3072 - temp[x1:x2, y1:y2] = np.fliplr(imgr[0:1536, 3072:6144]) - ########## zone 1, OSE,down to up ############################ + temp[x1:x2, y1:y2] = imgr[0:1536, 0:3072] + ########## third part , old OSE , down to up ############################ x1 = 0 x2 = x1+1536 y1 = 3442*2+prescan y2 = y1+3072 temp[x1:x2, y1:y2] = np.flipud(imgr[1536:3072, 0:3072]) - ########## zone 2, OSF down to yp ,left to right ################ + ########## fourth part, old OSF, down to up ,left to right ################ x1 = 0 x2 = x1+1536 y1 = 3442*3+prescan y2 = y1+3072 temp[x1:x2, y1:y2] = np.flipud(np.fliplr(imgr[1536:3072, 3072:6144])) + self.image_r = temp @@ -3383,6 +3496,18 @@ class IFSsimulator(): Updates header with the input values and flags used during simulation. """ + + ###### + ###### add random pointing error to telescope, the telescope poingt parameter in + ###### fits header have random pointing error + ud_ra = np.random.random() # Choose a random shift in arcsec + self.information['ra_pnt0']=self.information['dec_pnt0']+0.01*ud_ra + + ud_dec= np.random.random() + self.information['dec_pnt0']=self.information['dec_pnt0']+0.01*ud_dec + + + HeaderTest = 'no' sim_ver = str(self.information['sim_ver']) @@ -4760,24 +4885,38 @@ class IFSsimulator(): self.information['ra_obj'] = a[0].header['RA'] # in degree self.information['dec_obj'] = a[0].header['DEC'] # in degree - disRa = (telRa-skyRa)/3600.0 # convert unit of arcsec to degree + disRa = (telRa-skyRa)/3600.0 # convert unit of arcsec to degree disDec = (telDec-skyDec)/3600.0 # convert unit of arcsec to degree # RA_PNT0 =OBJ_RA + DIS_RA/cos(OBJ_DEC), DEC_PNT0= OBJ_DEC DIS_DEC self.information['ra_pnt0'] = a[0].header['RA'] + \ disRa/np.cos(a[0].header['DEC']/180.0*np.pi) + self.information['dec_pnt0'] = a[0].header['DEC'] + disDec ############### calculate the earthshine and zodiacal noise ,new code 2023.11.1 ############ ############### + self.log.info('Real telescope pointing in Ra = %f, Dec = %f' % (self.information['ra_pnt0'], self.information['dec_pnt0'])) + ra = self.information['ra_pnt0'] dec = self.information['dec_pnt0'] time_jd = time2jd(self.dt) + + if self.appbianpai: + + + sn=self.simnumber-1; + x_sat = float(self.bianpai_data['x_sat'][sn*5+self.exptime_start_index]) + y_sat = float(self.bianpai_data['y_sat'][sn*5+self.exptime_start_index]) + z_sat = float(self.bianpai_data['z_sat'][sn*5+self.exptime_start_index]) + + ### + else: - x_sat = float(self.orbit_pars[self.orbit_exp_num, 1]) - y_sat = float(self.orbit_pars[self.orbit_exp_num, 2]) - z_sat = float(self.orbit_pars[self.orbit_exp_num, 3]) + x_sat = float(self.orbit_pars[self.orbit_exp_num, 1]) + y_sat = float(self.orbit_pars[self.orbit_exp_num, 2]) + z_sat = float(self.orbit_pars[self.orbit_exp_num, 3]) wave0, zodi0 = self.zodiacal( ra, dec, self.TianCe_day) # erg/s/cm^2/A/arcsec^2 @@ -4856,14 +4995,26 @@ class IFSsimulator(): width_blue = 0 ################################ ############## doppler effect to photons.wavelength ############# - - # self.orbit_pars - x_sat = float(self.orbit_pars[self.orbit_exp_num, 1]) - y_sat = float(self.orbit_pars[self.orbit_exp_num, 2]) - z_sat = float(self.orbit_pars[self.orbit_exp_num, 3]) - vx_sat = float(self.orbit_pars[self.orbit_exp_num, 4]) - vy_sat = float(self.orbit_pars[self.orbit_exp_num, 5]) - vz_sat = float(self.orbit_pars[self.orbit_exp_num, 6]) + if self.appbianpai: + + sn=self.simnumber-1; + x_sat = float(self.bianpai_data['x_sat'][sn*5+self.exptime_start_index]) + y_sat = float(self.bianpai_data['y_sat'][sn*5+self.exptime_start_index]) + z_sat = float(self.bianpai_data['z_sat'][sn*5+self.exptime_start_index]) + + vx_sat = float(self.bianpai_data['vx_sat'][sn*5+self.exptime_start_index]) + vy_sat = float(self.bianpai_data['vy_sat'][sn*5+self.exptime_start_index]) + vz_sat = float(self.bianpai_data['vz_sat'][sn*5+self.exptime_start_index]) + + else: + + # self.orbit_pars + x_sat = float(self.orbit_pars[self.orbit_exp_num, 1]) + y_sat = float(self.orbit_pars[self.orbit_exp_num, 2]) + z_sat = float(self.orbit_pars[self.orbit_exp_num, 3]) + vx_sat = float(self.orbit_pars[self.orbit_exp_num, 4]) + vy_sat = float(self.orbit_pars[self.orbit_exp_num, 5]) + vz_sat = float(self.orbit_pars[self.orbit_exp_num, 6]) self.information['POSI0_X'] = x_sat self.information['POSI0_Y'] = y_sat @@ -4881,56 +5032,69 @@ class IFSsimulator(): self.information['ra_obj'], self.information['dec_obj'], v1, self.TianCe_day) ################################################# - # exposure end time is t2 ; - t2 = self.dt+timedelta(seconds=self.information['exptime']) - ### data read time is the exposure end time plus readouttime ### - # data read end time is t3 - - t2jd = time2jd(t2) - - if self.orbit_pars[-1, 0] < t2jd: # orbit parameters are not in currenct txt file - self.orbit_file_num = self.orbit_file_num+1 - fn = self.information['dir_path'] + \ - 'IFS_inputdata/TianCe/orbit20160925/' + \ - str(self.orbit_file_num)+'.txt' - self.orbit_pars = np.loadtxt(fn) - self.orbit_exp_num = 0 - - for k in range(self.orbit_exp_num, len(self.orbit_pars), 1): - - if t2jd-self.orbit_pars[k, 0] < 0: - break - if k == 0: - deltaT = jd2time(self.orbit_pars[k, 0])-t2 - p1x = self.orbit_pars[k, 1]-(self.orbit_pars[k+1, 1] - - self.orbit_pars[k, 1])*deltaT.seconds/120 - p1y = self.orbit_pars[k, 2]-(self.orbit_pars[k+1, 2] - - self.orbit_pars[k, 2])*deltaT.seconds/120 - p1z = self.orbit_pars[k, 3]-(self.orbit_pars[k+1, 3] - - self.orbit_pars[k, 3])*deltaT.seconds/120 - - p1vx = self.orbit_pars[k, 4]-(self.orbit_pars[k+1, 4] - - self.orbit_pars[k, 4])*deltaT.seconds/120 - p1vx = self.orbit_pars[k, 5]-(self.orbit_pars[k+1, 5] - - self.orbit_pars[k, 5])*deltaT.seconds/120 - p1vx = self.orbit_pars[k, 6]-(self.orbit_pars[k+1, 6] - - self.orbit_pars[k, 6])*deltaT.seconds/120 - + if self.appbianpai: + sn=self.simnumber-1; + + p1x = float(self.bianpai_data['x_sat'][sn*5+self.exptime_end_index]) + p1y = float(self.bianpai_data['y_sat'][sn*5+self.exptime_end_index]) + p1z = float(self.bianpai_data['z_sat'][sn*5+self.exptime_end_index]) + + p1vx = float(self.bianpai_data['vx_sat'][sn*5+self.exptime_end_index]) + p1vy = float(self.bianpai_data['vy_sat'][sn*5+self.exptime_end_index]) + p1vz = float(self.bianpai_data['vz_sat'][sn*5+self.exptime_end_index]) + else: - deltaT = jd2time(self.orbit_pars[k, 0])-t2 - p1x = self.orbit_pars[k-1, 1]+(self.orbit_pars[k, 1] - - self.orbit_pars[k-1, 1])*deltaT.seconds/120 - p1y = self.orbit_pars[k-1, 2]+(self.orbit_pars[k, 2] - - self.orbit_pars[k-1, 2])*deltaT.seconds/120 - p1z = self.orbit_pars[k-1, 3]+(self.orbit_pars[k, 3] - - self.orbit_pars[k-1, 3])*deltaT.seconds/120 - - p1vx = self.orbit_pars[k-1, 4]+(self.orbit_pars[k, 4] - - self.orbit_pars[k-1, 4])*deltaT.seconds/120 - p1vy = self.orbit_pars[k-1, 5]+(self.orbit_pars[k, 5] - - self.orbit_pars[k-1, 5])*deltaT.seconds/120 - p1vz = self.orbit_pars[k-1, 6]+(self.orbit_pars[k, 6] - - self.orbit_pars[k-1, 6])*deltaT.seconds/120 + + # exposure end time is t2 ; + t2 = self.dt+timedelta(seconds=self.information['exptime']) + ### data read time is the exposure end time plus readouttime ### + # data read end time is t3 + + t2jd = time2jd(t2) + + if self.orbit_pars[-1, 0] < t2jd: # orbit parameters are not in currenct txt file + self.orbit_file_num = self.orbit_file_num+1 + fn = self.information['dir_path'] + \ + 'IFS_inputdata/TianCe/orbit20160925/' + \ + str(self.orbit_file_num)+'.txt' + self.orbit_pars = np.loadtxt(fn) + self.orbit_exp_num = 0 + + for k in range(self.orbit_exp_num, len(self.orbit_pars), 1): + + if t2jd-self.orbit_pars[k, 0] < 0: + break + if k == 0: + deltaT = jd2time(self.orbit_pars[k, 0])-t2 + p1x = self.orbit_pars[k, 1]-(self.orbit_pars[k+1, 1] - + self.orbit_pars[k, 1])*deltaT.seconds/120 + p1y = self.orbit_pars[k, 2]-(self.orbit_pars[k+1, 2] - + self.orbit_pars[k, 2])*deltaT.seconds/120 + p1z = self.orbit_pars[k, 3]-(self.orbit_pars[k+1, 3] - + self.orbit_pars[k, 3])*deltaT.seconds/120 + + p1vx = self.orbit_pars[k, 4]-(self.orbit_pars[k+1, 4] - + self.orbit_pars[k, 4])*deltaT.seconds/120 + p1vy = self.orbit_pars[k, 5]-(self.orbit_pars[k+1, 5] - + self.orbit_pars[k, 5])*deltaT.seconds/120 + p1vz = self.orbit_pars[k, 6]-(self.orbit_pars[k+1, 6] - + self.orbit_pars[k, 6])*deltaT.seconds/120 + + else: + deltaT = jd2time(self.orbit_pars[k, 0])-t2 + p1x = self.orbit_pars[k-1, 1]+(self.orbit_pars[k, 1] - + self.orbit_pars[k-1, 1])*deltaT.seconds/120 + p1y = self.orbit_pars[k-1, 2]+(self.orbit_pars[k, 2] - + self.orbit_pars[k-1, 2])*deltaT.seconds/120 + p1z = self.orbit_pars[k-1, 3]+(self.orbit_pars[k, 3] - + self.orbit_pars[k-1, 3])*deltaT.seconds/120 + + p1vx = self.orbit_pars[k-1, 4]+(self.orbit_pars[k, 4] - + self.orbit_pars[k-1, 4])*deltaT.seconds/120 + p1vy = self.orbit_pars[k-1, 5]+(self.orbit_pars[k, 5] - + self.orbit_pars[k-1, 5])*deltaT.seconds/120 + p1vz = self.orbit_pars[k-1, 6]+(self.orbit_pars[k, 6] - + self.orbit_pars[k-1, 6])*deltaT.seconds/120 ####### self.information['POSI1_X'] = p1x @@ -4961,7 +5125,7 @@ class IFSsimulator(): ilam = klam*alfa - if ilam >= 6000: + if ilam >= 4000: break # print(ilam) @@ -5692,7 +5856,7 @@ class IFSsimulator(): self.debug = self.information['debug'] - self.dt = datetime.utcnow() + if self.information['exptime']>2000 or self.information['exptime']<0: @@ -5718,65 +5882,89 @@ class IFSsimulator(): self.information['sky_fitsin'] print('self.skyfilepath = ', self.skyfilepath) - - # self.earthshine_theta=30.0 # in degree - - ################################################################## - #### load orbit parameters ##### - flag = 0 - - self.dt_num = int( - self.simnumber*(self.information['exptime']+self.information['readouttime']+125)/120) - now_dt = datetime.utcnow() - now_jd = time2jd(now_dt) - - for k in range(1, 50, 1): - - # fn=father_path+'/IFS_inputdata/TianCe/orbit20160925/'+str(k)+'.txt'; - - fn = self.information['dir_path'] + \ - 'IFS_inputdata/TianCe/orbit20160925/'+str(k)+'.txt' - d = np.loadtxt(fn) - - for kk in range(len(d[:, 0])): - if now_jd-d[kk, 0] <= 0: - flag = 1 - break - if flag == 1: - break - - ################################ - - if kk + self.dt_num < len(d[:, 0]): - - exptime_start_jd = d[kk+self.dt_num, 0] - - self.orbit_pars = d - self.orbit_file_num = k - self.orbit_exp_num = kk - + + ################################################################### + if self.appbianpai: + ### load yunxingbianpai csv file + ############ load star data catlog ##################### + starcat=self.information['bianpai_file'] + + ##starcat='selection_20230517_concat.fits' + ################################################### + + self.log.info('Stat catlog file name is %s' % (starcat)) + + ########################################## + df=pd.read_csv(self.information['dir_path']+'IFS_inputdata/TianCe/'+starcat) + + ################################################################### + sn=self.simnumber-1; + arr=np.array(df['time'][sn*5:sn*5+5]); + self.exptime_start_jd,self.exptime_start_index=find_min(arr); + self.exptime_end_jd, self.exptime_end_index =find_max(arr); + + ###self.earthshine_theta=df['earth_angle'][sn*5+index] # in degree + self.dt = jd2time(self.exptime_start_jd); + self.bianpai_data=df; + ################################################################## else: + + #### load orbit parameters ##### + flag = 0 + self.dt = datetime.utcnow() + self.dt_num = int( + self.simnumber*(self.information['exptime']+self.information['readouttime']+125)/120) + now_dt = datetime.utcnow() + now_jd = time2jd(now_dt) + + + + for k in range(1, 50, 1): + + # fn=father_path+'/IFS_inputdata/TianCe/orbit20160925/'+str(k)+'.txt'; + + fn = self.information['dir_path'] + \ + 'IFS_inputdata/TianCe/orbit20160925/'+str(k)+'.txt' + d = np.loadtxt(fn) + + for kk in range(len(d[:, 0])): + if now_jd-d[kk, 0] <= 0: + flag = 1 + break + if flag == 1: + break + + ################################################################## + + if kk + self.dt_num < len(d[:, 0]): + + exptime_start_jd = d[kk+self.dt_num, 0] + + self.orbit_pars = d + self.orbit_file_num = k + self.orbit_exp_num = kk + + else: + + fn = self.information['dir_path'] + \ + 'IFS_inputdata/TianCe/orbit20160925/'+str(k+1)+'.txt' + d = np.loadtxt(fn) + + self.orbit_pars = d + self.orbit_file_num = k+1 + self.orbit_exp_num = self.dt_num + self.dt = jd2time(exptime_start_jd) - fn = self.information['dir_path'] + \ - 'IFS_inputdata/TianCe/orbit20160925/'+str(k+1)+'.txt' - d = np.loadtxt(fn) - - self.orbit_pars = d - self.orbit_file_num = k+1 - self.orbit_exp_num = self.dt_num - - ########################################## - ######################################################################## - #self.dt=julian.from_jd(exptime_start_jd, fmt='jd') - self.dt = jd2time(exptime_start_jd) + ################################################################## ################################################################## - # str(self.dt.year)+'-'+str(self.dt.month)+'-'+str(self.dt.day) + self.TianCe_day = self.dt.strftime("%Y-%m-%d") self.TianCe_exp_start = dt2hmd(self.dt) self.zodiacal_time = self.TianCe_day + ################################################################### skyRa = 0 skyDec = 0 @@ -5787,6 +5975,7 @@ class IFSsimulator(): sky_rot = 2 * (ud-0.5) * 5 dsmax = np.floor(50*np.cos(sky_rot/180*np.pi))-34 + np.random.seed(11*self.simnumber) ud = np.random.random() # Choose a random shift in arcsec telRa = 2 * (ud-0.5) * dsmax * 0.1 # in arcsec @@ -5813,6 +6002,7 @@ class IFSsimulator(): ####################################################################### elif self.source == 'FLAT': + self.dt = datetime.utcnow() self.information['exptime'] = 200+simnumber*100 self.sim_calibration_img(self.information['exptime'], 'FLAT') self.information['sky_rot'] = 0 @@ -5820,6 +6010,7 @@ class IFSsimulator(): ######## elif self.source == 'LAMP': + self.dt = datetime.utcnow() self.information['exptime'] = 200+simnumber*100 self.information['sky_rot'] = 0 @@ -5831,12 +6022,14 @@ class IFSsimulator(): ######### elif self.source == 'DARK': + self.dt = datetime.utcnow() self.information['sky_rot'] = 0 self.information['exptime'] = simnumber*3600*24 self.image_b = np.zeros((2048, 4096)) self.image_r = np.zeros((3072, 6144)) elif self.source == 'BIAS': + self.dt = datetime.utcnow() self.information['sky_rot'] = 0 self.information['exptime'] = 0 self.image_b = np.zeros((2048, 4096))