ObservationSim.py 22.2 KB
Newer Older
1
import os
Fang Yuedong's avatar
Fang Yuedong committed
2
3
4
5
6
import numpy as np
import mpi4py.MPI as MPI
import galsim
import logging
import psutil
Fang Yuedong's avatar
Fang Yuedong committed
7
import gc
8
9
10
from astropy.io import fits
from datetime import datetime

Fang Yuedong's avatar
Fang Yuedong committed
11
12
import traceback

13
14
15
from ObservationSim.Config import config_dir, ChipOutput
from ObservationSim.Config.Header import generatePrimaryHeader, generateExtensionHeader
from ObservationSim.Instrument import Telescope, Filter, FilterParam, FocalPlane, Chip
16
from ObservationSim.Instrument.Chip import Effects
17
from ObservationSim.Straylight import calculateSkyMap_split_g
Zhang Xin's avatar
Zhang Xin committed
18
from ObservationSim.PSF import PSFGauss, FieldDistortion, PSFInterp, PSFInterpSLS
19
from ObservationSim._util import get_shear_field, makeSubDir_PointingList
20
from ObservationSim.Astrometry.Astrometry_util import on_orbit_obs_position
Fang Yuedong's avatar
Fang Yuedong committed
21
22

class Observation(object):
23
    def __init__(self, config, Catalog, work_dir=None, data_dir=None):
24
        self.path_dict = config_dir(config=config, work_dir=work_dir, data_dir=data_dir)
Fang Yuedong's avatar
Fang Yuedong committed
25
        self.config = config
Fang Yuedong's avatar
Fang Yuedong committed
26
        self.tel = Telescope()
27
        self.focal_plane = FocalPlane(survey_type=self.config["obs_setting"]["survey_type"]) 
Fang Yuedong's avatar
Fang Yuedong committed
28
        self.filter_param = FilterParam() 
Fang Yuedong's avatar
Fang Yuedong committed
29
30
        self.chip_list = []
        self.filter_list = []
31
        self.all_filter = []
32
        self.Catalog = Catalog
Fang Yuedong's avatar
Fang Yuedong committed
33
34

        # Construct chips & filters:
Fang Yuedong's avatar
Fang Yuedong committed
35
        for i in range(self.focal_plane.nchips):
Fang Yuedong's avatar
Fang Yuedong committed
36
37
38
            chipID = i + 1

            # Make Chip & Filter lists
39
40
            chip = Chip(
                chipID=chipID, 
Fang Yuedong's avatar
Fang Yuedong committed
41
                config=self.config)
Fang Yuedong's avatar
Fang Yuedong committed
42
            filter_id, filter_type = chip.getChipFilter()
43
44
            filt = Filter(filter_id=filter_id, 
                filter_type=filter_type, 
45
                filter_param=self.filter_param)
46
47
48
49
            if not self.focal_plane.isIgnored(chipID=chipID):
                self.chip_list.append(chip)
                self.filter_list.append(filt)
            self.all_filter.append(filt)
Fang Yuedong's avatar
Fang Yuedong committed
50

51
    def run_one_chip(self, chip, filt, pointing, chip_output, wcs_fp=None, psf_model=None, cat_dir=None, sed_dir=None):
Fang Yuedong's avatar
Fang Yuedong committed
52

Fang Yuedong's avatar
Fang Yuedong committed
53
54
55
56
57
58
59
60
61
        chip_output.Log_info(':::::::::::::::::::Current Pointing Information::::::::::::::::::')
        chip_output.Log_info("RA: %f, DEC; %f" % (pointing.ra, pointing.dec))
        chip_output.Log_info("Time: %s" % datetime.utcfromtimestamp(pointing.timestamp).isoformat())
        chip_output.Log_info("Exposure time: %f" % pointing.exp_time)
        chip_output.Log_info("Satellite Position (x, y, z): (%f, %f, %f)" % (pointing.sat_x, pointing.sat_y, pointing.sat_z))
        chip_output.Log_info("Satellite Velocity (x, y, z): (%f, %f, %f)" % (pointing.sat_vx, pointing.sat_vy, pointing.sat_vz))
        chip_output.Log_info("Position Angle: %f" % pointing.img_pa.deg)
        chip_output.Log_info('Chip : %d' % chip.chipID)
        chip_output.Log_info(':::::::::::::::::::::::::::END:::::::::::::::::::::::::::::::::::')
Fang Yuedong's avatar
Fang Yuedong committed
62

Fang Yuedong's avatar
Fang Yuedong committed
63
        if self.config["psf_setting"]["psf_model"] == "Gauss":
Fang Yuedong's avatar
Fang Yuedong committed
64
            psf_model = PSFGauss(chip=chip, psfRa=self.config["psf_setting"]["psf_rcont"])
Fang Yuedong's avatar
Fang Yuedong committed
65
        elif self.config["psf_setting"]["psf_model"] == "Interp":
Zhang Xin's avatar
Zhang Xin committed
66
67
68
69
            if chip.survey_type == "spectroscopic":
                psf_model = PSFInterpSLS(chip, filt,PSF_data_prefix=self.path_dict["psf_sls_dir"])
            else:
                psf_model = PSFInterp(chip=chip, npsf=chip.n_psf_samples, PSF_data_file=self.path_dict["psf_dir"])
Fang Yuedong's avatar
Fang Yuedong committed
70
        else:
Fang Yuedong's avatar
Fang Yuedong committed
71
            chip_output.Log_error("unrecognized PSF model type!!", flush=True)
Fang Yuedong's avatar
Fang Yuedong committed
72

Fang Yuedong's avatar
Fang Yuedong committed
73
        # Figure out shear fields
74
        self.g1_field, self.g2_field, self.nshear = get_shear_field(config=self.config)
Fang Yuedong's avatar
Fang Yuedong committed
75

76
77
        # Apply astrometric simulation for pointing
        if self.config["obs_setting"]["enable_astrometric_model"]:
Fang Yuedong's avatar
Fang Yuedong committed
78
            dt = datetime.utcfromtimestamp(pointing.timestamp)
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
            date_str = dt.date().isoformat()
            time_str = dt.time().isoformat()
            ra_cen, dec_cen = on_orbit_obs_position(
                input_ra_list=[pointing.ra],
                input_dec_list=[pointing.dec],
                input_pmra_list=[0.],
                input_pmdec_list=[0.],
                input_rv_list=[0.],
                input_parallax_list=[1e-9],
                input_nstars=1,
                input_x=pointing.sat_x,
                input_y=pointing.sat_y,
                input_z=pointing.sat_z,
                input_vx=pointing.sat_vx,
                input_vy=pointing.sat_vy,
                input_vz=pointing.sat_vz,
Fang Yuedong's avatar
Fang Yuedong committed
95
                input_epoch="J2000",
96
97
98
99
100
101
102
103
                input_date_str=date_str,
                input_time_str=time_str
            )
            ra_cen, dec_cen = ra_cen[0], dec_cen[0]
        else:
            ra_cen = pointing.ra
            dec_cen = pointing.dec

Fang Yuedong's avatar
Fang Yuedong committed
104
105
        # Get WCS for the focal plane
        if wcs_fp == None:
106
            wcs_fp = self.focal_plane.getTanWCS(ra_cen, dec_cen, pointing.img_pa, chip.pix_scale)
Fang Yuedong's avatar
Fang Yuedong committed
107
108
109
110
111

        # Create chip Image
        chip.img = galsim.ImageF(chip.npix_x, chip.npix_y)
        chip.img.setOrigin(chip.bound.xmin, chip.bound.ymin)
        chip.img.wcs = wcs_fp
112

Zhang Xin's avatar
Zhang Xin committed
113
114
        if self.config["obs_setting"]["enable_straylight_model"]:
            filt.setFilterStrayLightPixel(jtime = pointing.jdt, sat_pos = np.array([pointing.sat_x, pointing.sat_y, pointing.sat_z]), pointing_radec = np.array([pointing.ra,pointing.dec]), sun_pos = np.array([pointing.sun_x,pointing.sun_y,pointing.sun_z]))
115

116
117
        chip_output.Log_info("========================sky pix========================")
        chip_output.Log_info(filt.sky_background)
118

Fang Yuedong's avatar
Fang Yuedong committed
119
120
121
        if chip.survey_type == "photometric":
            sky_map = None
        elif chip.survey_type == "spectroscopic":
xin's avatar
xin committed
122
            # chip.loadSLSFLATCUBE(flat_fn='flat_cube.fits')
123
124
            flat_normal = np.ones_like(chip.img.array)
            if self.config["ins_effects"]["flat_fielding"] == True:
Fang Yuedong's avatar
Fang Yuedong committed
125
                chip_output.Log_info("SLS flat preprocess,CHIP %d : Creating and applying Flat-Fielding"%chip.chipID)
126
                msg = str(chip.img.bounds)
Fang Yuedong's avatar
Fang Yuedong committed
127
                chip_output.Log_info(msg)
128
129
130
131
132
                flat_img = Effects.MakeFlatSmooth(
                    chip.img.bounds,
                    int(self.config["random_seeds"]["seed_flat"]))
                flat_normal = flat_normal * flat_img.array / np.mean(flat_img.array)
            if self.config["ins_effects"]["shutter_effect"] == True:
Fang Yuedong's avatar
Fang Yuedong committed
133
                chip_output.Log_info("SLS flat preprocess,CHIP %d : Apply shutter effect"%chip.chipID)
134
135
136
137
                shuttimg = Effects.ShutterEffectArr(chip.img, t_shutter=1.3, dist_bearing=735,
                                                    dt=1E-3)  # shutter effect normalized image for this chip
                flat_normal = flat_normal*shuttimg
                flat_normal = np.array(flat_normal,dtype='float32')
Fang Yuedong's avatar
Fang Yuedong committed
138
            sky_map = calculateSkyMap_split_g(
139
140
141
142
143
144
                skyMap=flat_normal,
                blueLimit=filt.blue_limit,
                redLimit=filt.red_limit,
                conf=chip.sls_conf,
                pixelSize=chip.pix_scale,
                isAlongY=0,
145
146
                flat_cube=chip.flat_cube, zoldial_spec = filt.zodical_spec)
            sky_map = sky_map+filt.sky_background
147
            del flat_normal
Fang Yuedong's avatar
Fang Yuedong committed
148

149
        if pointing.pointing_type == 'MS':
Fang Yuedong's avatar
Fang Yuedong committed
150
            # Load catalogues and templates
Fang Yuedong's avatar
Fang Yuedong committed
151
            self.cat = self.Catalog(config=self.config, chip=chip, pointing=pointing, cat_dir=cat_dir, sed_dir=sed_dir, chip_output=chip_output, filt=filt)
152
            chip_output.create_output_file()
Fang Yuedong's avatar
Fang Yuedong committed
153
154
            self.nobj = len(self.cat.objs)

Fang Yuedong's avatar
Fang Yuedong committed
155
156
157
158
159
160
161
162
163
            for ifilt in range(len(self.all_filter)):
                temp_filter = self.all_filter[ifilt]
                # Update the limiting magnitude using exposure time in pointing
                temp_filter.update_limit_saturation_mags(exptime=pointing.exp_time, chip=chip)

                # Select cutting band filter for saturation/limiting magnitude
                if temp_filter.filter_type.lower() == self.config["obs_setting"]["cut_in_band"].lower():
                    cut_filter = temp_filter

Fang Yuedong's avatar
Fang Yuedong committed
164
            if self.config["ins_effects"]["field_dist"] == True:
165
                self.fd_model = FieldDistortion(chip=chip, img_rot=pointing.img_pa.deg)
Fang Yuedong's avatar
Fang Yuedong committed
166
167
168
            else:
                self.fd_model = None

Fang Yuedong's avatar
Fang Yuedong committed
169
170
171
172
            # Loop over objects
            missed_obj = 0
            bright_obj = 0
            dim_obj = 0
Fang Yuedong's avatar
Fang Yuedong committed
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
 
            h_ext = generateExtensionHeader(
                chip=chip,
                xlen=chip.npix_x, 
                ylen=chip.npix_y, 
                ra=pointing.ra, 
                dec=pointing.dec, 
                pa=pointing.img_pa.deg, 
                gain=chip.gain, 
                readout=chip.read_noise, 
                dark=chip.dark_noise, 
                saturation=90000, 
                pixel_scale=chip.pix_scale, 
                pixel_size=chip.pix_size,
                xcen=chip.x_cen,
                ycen=chip.y_cen,
189
190
191
                extName='SCI',
                timestamp = pointing.timestamp,
                exptime = pointing.exp_time,
192
                readoutTime = chip.readout_time)
193
194
            
            chip_wcs = galsim.FitsWCS(header=h_ext)
Fang Yuedong's avatar
Fang Yuedong committed
195

Fang Yuedong's avatar
Fang Yuedong committed
196
            for j in range(self.nobj):
197
198
                
                # (DEBUG)
199
                # if j >= 10:
Xin Zhang's avatar
Xin Zhang committed
200
                #     break
201

Fang Yuedong's avatar
Fang Yuedong committed
202
203
                obj = self.cat.objs[j]

Zhang Xin's avatar
Zhang Xin committed
204
205
206
207
                # (DEBUG)
                # if obj.getMagFilter(filt)>20:
                #     continue

208
                # load and convert SED; also caculate object's magnitude in all CSST bands
209
210
211
                try:
                    sed_data = self.cat.load_sed(obj)
                    norm_filt = self.cat.load_norm_filt(obj)
Fang Yuedong's avatar
Fang Yuedong committed
212
                    obj.sed, obj.param["mag_%s"%filt.filter_type.lower()], obj.param["flux_%s"%filt.filter_type.lower()] = self.cat.convert_sed(
213
214
215
216
217
                        mag=obj.param["mag_use_normal"],
                        sed=sed_data,
                        target_filt=filt, 
                        norm_filt=norm_filt,
                    )
Fang Yuedong's avatar
Fang Yuedong committed
218
                    _, obj.param["mag_%s"%cut_filter.filter_type.lower()], obj.param["flux_%s"%cut_filter.filter_type.lower()] = self.cat.convert_sed(
Fang Yuedong's avatar
Fang Yuedong committed
219
220
221
222
223
                        mag=obj.param["mag_use_normal"],
                        sed=sed_data,
                        target_filt=cut_filter, 
                        norm_filt=norm_filt,
                    )
224
                except Exception as e:
Fang Yuedong's avatar
Fang Yuedong committed
225
                    traceback.print_exc()
Fang Yuedong's avatar
Fang Yuedong committed
226
                    chip_output.Log_error(e)
227
                    continue
228
                
Fang Yuedong's avatar
Fang Yuedong committed
229
                # [TODO] Testing
230
                # chip_output.Log_info("mag_%s = %.3f"%(filt.filter_type.lower(), obj.param["mag_%s"%filt.filter_type.lower()]))
Fang Yuedong's avatar
Fang Yuedong committed
231
232

                # Exclude very bright/dim objects (for now)
233
234
235
                if cut_filter.is_too_bright(
                    mag=obj.param["mag_%s"%self.config["obs_setting"]["cut_in_band"].lower()],
                    margin=self.config["obs_setting"]["mag_sat_margin"]):
Fang Yuedong's avatar
Fang Yuedong committed
236
237
238
239
                    chip_output.Log_info("obj %s too birght!! mag_%s = %.3f"%(obj.id, cut_filter.filter_type, obj.param["mag_%s"%self.config["obs_setting"]["cut_in_band"].lower()]))
                    bright_obj += 1
                    obj.unload_SED()
                    continue
240
241
242
                if filt.is_too_dim(
                    mag=obj.getMagFilter(filt),
                    margin=self.config["obs_setting"]["mag_lim_margin"]):
Fang Yuedong's avatar
Fang Yuedong committed
243
                    chip_output.Log_info("obj %s too dim!! mag_%s = %.3f"%(obj.id, filt.filter_type, obj.getMagFilter(filt)))
Fang Yuedong's avatar
Fang Yuedong committed
244
                    dim_obj += 1
Fang Yuedong's avatar
Fang Yuedong committed
245
246
247
                    obj.unload_SED()
                    continue

Fang Yuedong's avatar
Fang Yuedong committed
248
                # Get corresponding shear values
Fang Yuedong's avatar
Fang Yuedong committed
249
                if self.config["shear_setting"]["shear_type"] == "constant":
Fang Yuedong's avatar
Fang Yuedong committed
250
                    if obj.type == 'star':
251
                        obj.g1, obj.g2 = 0., 0.
Fang Yuedong's avatar
Fang Yuedong committed
252
                    else:
253
                        obj.g1, obj.g2 = self.g1_field, self.g2_field
254
255
256
                elif self.config["shear_setting"]["shear_type"] == "catalog":
                    pass
                else:
Fang Yuedong's avatar
Fang Yuedong committed
257
                    chip_output.Log_error("Unknown shear input")
258
                    raise ValueError("Unknown shear input")
Fang Yuedong's avatar
Fang Yuedong committed
259
260

                # Get position of object on the focal plane
261
                pos_img, offset, local_wcs, real_wcs, fd_shear = obj.getPosImg_Offset_WCS(img=chip.img, fdmodel=self.fd_model, chip=chip, verbose=False, chip_wcs=chip_wcs, img_header=h_ext)
Fang Yuedong's avatar
Fang Yuedong committed
262
263
264
265

                # [TODO] For now, only consider objects which their centers (after field distortion) are projected within the focal plane
                # Otherwise they will be considered missed objects
                # if pos_img.x == -1 or pos_img.y == -1 or (not chip.isContainObj(x_image=pos_img.x, y_image=pos_img.y, margin=0.)):
Fang Yuedong's avatar
Fang Yuedong committed
266
                if pos_img.x == -1 or pos_img.y == -1:
Fang Yuedong's avatar
Fang Yuedong committed
267
268
                    chip_output.Log_info('obj_ra = %.6f, obj_dec = %.6f, obj_ra_orig = %.6f, obj_dec_orig = %.6f'%(obj.ra, obj.dec, obj.ra_orig, obj.dec_orig))
                    chip_output.Log_error("Objected missed: %s"%(obj.id))
Fang Yuedong's avatar
Fang Yuedong committed
269
270
271
                    missed_obj += 1
                    obj.unload_SED()
                    continue
Fang Yuedong's avatar
Fang Yuedong committed
272

Fang Yuedong's avatar
Fang Yuedong committed
273
274
                # Draw object & update output catalog
                try:
Fang Yuedong's avatar
Fang Yuedong committed
275
                    if self.config["run_option"]["out_cat_only"]:
Fang Yuedong's avatar
Fang Yuedong committed
276
                        isUpdated = True
Fang Yuedong's avatar
Fang Yuedong committed
277
                        obj.real_pos = obj.getRealPos(chip.img, global_x=obj.posImg.x, global_y=obj.posImg.y, img_real_wcs=obj.real_wcs)
278
                        pos_shear = 0.
279
                    elif chip.survey_type == "photometric" and not self.config["run_option"]["out_cat_only"]:
Fang Yuedong's avatar
Fang Yuedong committed
280
281
282
283
284
285
286
                        isUpdated, pos_shear = obj.drawObj_multiband(
                            tel=self.tel,
                            pos_img=pos_img, 
                            psf_model=psf_model, 
                            bandpass_list=filt.bandpass_sub_list, 
                            filt=filt, 
                            chip=chip, 
287
288
                            g1=obj.g1, 
                            g2=obj.g2, 
Fang Yuedong's avatar
Fang Yuedong committed
289
290
291
                            exptime=pointing.exp_time,
                            fd_shear=fd_shear)

Fang Yuedong's avatar
Fang Yuedong committed
292
                    elif chip.survey_type == "spectroscopic" and not self.config["run_option"]["out_cat_only"]:
Fang Yuedong's avatar
Fang Yuedong committed
293
294
295
296
297
298
299
                        isUpdated, pos_shear = obj.drawObj_slitless(
                            tel=self.tel, 
                            pos_img=pos_img, 
                            psf_model=psf_model, 
                            bandpass_list=filt.bandpass_sub_list, 
                            filt=filt, 
                            chip=chip, 
300
301
                            g1=obj.g1, 
                            g2=obj.g2, 
302
                            exptime=pointing.exp_time,
Fang Yuedong's avatar
Fang Yuedong committed
303
304
305
306
                            normFilter=norm_filt,
                            fd_shear=fd_shear)

                    if isUpdated == 1:
Fang Yuedong's avatar
Fang Yuedong committed
307
                        # TODO: add up stats
308
                        chip_output.cat_add_obj(obj, pos_img, pos_shear)
Fang Yuedong's avatar
Fang Yuedong committed
309
                        pass
Fang Yuedong's avatar
Fang Yuedong committed
310
311
312
                    elif isUpdated == 0:
                        missed_obj += 1
                        chip_output.Log_error("Objected missed: %s"%(obj.id))
Fang Yuedong's avatar
Fang Yuedong committed
313
                    else:
Fang Yuedong's avatar
Fang Yuedong committed
314
                        chip_output.Log_error("Draw error, object omitted: %s"%(obj.id))
Fang Yuedong's avatar
Fang Yuedong committed
315
316
                        continue
                except Exception as e:
Fang Yuedong's avatar
Fang Yuedong committed
317
                    traceback.print_exc()
Fang Yuedong's avatar
Fang Yuedong committed
318
319
320
321
322
323
324
                    chip_output.Log_error(e)
                    
                # # [C6 TEST]
                # chip_output.Log_info("check running:1: pointing-{:} chip-{:} pid-{:} memory-{:6.2}GB".format(pointing.id, chip.chipID, os.getpid(), (psutil.Process(os.getpid()).memory_info().rss / 1024 / 1024 / 1024) ))
                # chip_output.Log_info('draw object %s'%obj.id)
                # chip_output.Log_info('mag = %.3f'%obj.param['mag_use_normal'])

Fang Yuedong's avatar
Fang Yuedong committed
325
326
327
                # Unload SED:
                obj.unload_SED()
                del obj
Fang Yuedong's avatar
Fang Yuedong committed
328
                gc.collect()
Fang Yuedong's avatar
Fang Yuedong committed
329

Fang Yuedong's avatar
Fang Yuedong committed
330
331
            del psf_model
            del self.cat
Fang Yuedong's avatar
Fang Yuedong committed
332
            gc.collect()
Fang Yuedong's avatar
Fang Yuedong committed
333

Fang Yuedong's avatar
Fang Yuedong committed
334
        chip_output.Log_info("check running:1: pointing-%d chip-%d pid-%d memory-%6.2fGB"%(pointing.id, chip.chipID, os.getpid(), (psutil.Process(os.getpid()).memory_info().rss / 1024 / 1024 / 1024) ))
Fang Yuedong's avatar
Fang Yuedong committed
335
336
337

        # Detector Effects
        # ===========================================================
Fang Yuedong's avatar
Fang Yuedong committed
338
        # whether to output zero, dark, flat calibration images.
339
340
341
342
343
344

        if not self.config["run_option"]["out_cat_only"]:
            chip.img = chip.addEffects(
                config=self.config, 
                img=chip.img, 
                chip_output=chip_output, 
345
                filt=filt,
346
347
348
                ra_cen=pointing.ra, 
                dec_cen=pointing.dec,
                img_rot=pointing.img_pa,
Fang Yuedong's avatar
Fang Yuedong committed
349
                exptime=pointing.exp_time,
350
351
352
353
354
355
356
                pointing_ID=pointing.id,
                timestamp_obs=pointing.timestamp,
                pointing_type=pointing.pointing_type,
                sky_map=sky_map, tel = self.tel,
                logger=chip_output.logger)
            
            if pointing.pointing_type == 'MS':
Fang Yuedong's avatar
Fang Yuedong committed
357
                datetime_obs = datetime.utcfromtimestamp(pointing.timestamp)
358
359
360
361
362
363
364
365
                date_obs = datetime_obs.strftime("%y%m%d")
                time_obs = datetime_obs.strftime("%H%M%S")
                h_prim = generatePrimaryHeader(
                    xlen=chip.npix_x, 
                    ylen=chip.npix_y, 
                    pointNum = str(pointing.id),
                    ra=pointing.ra, 
                    dec=pointing.dec, 
Fang Yuedong's avatar
Fang Yuedong committed
366
                    pixel_scale=chip.pix_scale,
367
368
369
370
371
                    date=date_obs,
                    time_obs=time_obs,
                    exptime=pointing.exp_time,
                    im_type='SCI',
                    sat_pos=[pointing.sat_x, pointing.sat_y, pointing.sat_z],
Fang Yuedong's avatar
Fang Yuedong committed
372
373
                    sat_vel=[pointing.sat_vx, pointing.sat_vy, pointing.sat_vz],
                    chip_name=str(chip.chipID).rjust(2, '0'))
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
                h_ext = generateExtensionHeader(
                    chip=chip,
                    xlen=chip.npix_x,
                    ylen=chip.npix_y,
                    ra=pointing.ra,
                    dec=pointing.dec,
                    pa=pointing.img_pa.deg,
                    gain=chip.gain,
                    readout=chip.read_noise,
                    dark=chip.dark_noise,
                    saturation=90000,
                    pixel_scale=chip.pix_scale,
                    pixel_size=chip.pix_size,
                    xcen=chip.x_cen,
                    ycen=chip.y_cen,
                    extName='SCI',
                    timestamp=pointing.timestamp,
                    exptime=pointing.exp_time,
392
                    readoutTime=chip.readout_time)
Fang Yuedong's avatar
Fang Yuedong committed
393
                
394
395
                chip.img = galsim.Image(chip.img.array, dtype=np.uint16)
                hdu1 = fits.PrimaryHDU(header=h_prim)
Fang Yuedong's avatar
Fang Yuedong committed
396
                hdu1.add_checksum()
397
398
                hdu1.header.comments['CHECKSUM'] = 'HDU checksum'
                hdu1.header.comments['DATASUM'] = 'data unit checksum'
399
                hdu2 = fits.ImageHDU(chip.img.array, header=h_ext)
Fang Yuedong's avatar
Fang Yuedong committed
400
                hdu2.add_checksum()
401
402
403
                hdu2.header.comments['XTENSION'] = 'extension type'
                hdu2.header.comments['CHECKSUM'] = 'HDU checksum'
                hdu2.header.comments['DATASUM'] = 'data unit checksum'
404
405
406
                hdu1 = fits.HDUList([hdu1, hdu2])
                fname = os.path.join(chip_output.subdir, h_prim['FILENAME'] + '.fits')
                hdu1.writeto(fname, output_verify='ignore', overwrite=True)
Fang Yuedong's avatar
Fang Yuedong committed
407
408
409
410

                chip_output.Log_info("# objects that are too bright %d out of %d"%(bright_obj, self.nobj))
                chip_output.Log_info("# objects that are too dim %d out of %d"%(dim_obj, self.nobj))
                chip_output.Log_info("# objects that are missed %d out of %d"%(missed_obj, self.nobj))
Fang Yuedong's avatar
Fang Yuedong committed
411
412
        del chip.img

Fang Yuedong's avatar
Fang Yuedong committed
413
        chip_output.Log_info("check running:2: pointing-%d chip-%d pid-%d memory-%6.2fGB"%(pointing.id, chip.chipID, os.getpid(), (psutil.Process(os.getpid()).memory_info().rss / 1024 / 1024 / 1024) ))
Fang Yuedong's avatar
Fang Yuedong committed
414

415
    def runExposure_MPI_PointingList(self, pointing_list,chips=None, use_mpi=False):
Fang Yuedong's avatar
Fang Yuedong committed
416
417
418
419
        if use_mpi:
            comm = MPI.COMM_WORLD
            ind_thread = comm.Get_rank()
            num_thread = comm.Get_size()
Fang Yuedong's avatar
Fang Yuedong committed
420

Fang Yuedong's avatar
Fang Yuedong committed
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
        if chips is None:
            nchips_per_fp = len(self.chip_list)
            run_chips = self.chip_list
            run_filts = self.filter_list
        else:
            # Only run a particular set of chips
            run_chips = []
            run_filts = []
            nchips_per_fp = len(chips)
            for ichip in range(len(self.chip_list)):
                chip = self.chip_list[ichip]
                filt = self.filter_list[ichip]
                if chip.chipID in chips:
                    run_chips.append(chip)
                    run_filts.append(filt)
Fang Yuedong's avatar
Fang Yuedong committed
436

437
        for ipoint in range(len(pointing_list)):
Fang Yuedong's avatar
Fang Yuedong committed
438
439
            for ichip in range(nchips_per_fp):
                i = ipoint*nchips_per_fp + ichip
440
441
                pointing = pointing_list[ipoint]
                pointing_ID = pointing.id
Fang Yuedong's avatar
Fang Yuedong committed
442
443
444
                if use_mpi:
                    if i % num_thread != ind_thread:
                        continue
Fang Yuedong's avatar
Fang Yuedong committed
445
446
447
448
449

                pid = os.getpid()

                sub_img_dir, prefix = makeSubDir_PointingList(path_dict=self.path_dict, config=self.config, pointing_ID=pointing_ID)

Fang Yuedong's avatar
Fang Yuedong committed
450
451
                chip = run_chips[ichip]
                filt = run_filts[ichip]
Fang Yuedong's avatar
Fang Yuedong committed
452
                # chip_output.Log_info("running pointing#%d, chip#%d, at PID#%d..."%(pointing_ID, chip.chipID, pid))
Fang Yuedong's avatar
Fang Yuedong committed
453
454
455
456
457
                chip_output = ChipOutput(
                    config=self.config, 
                    focal_plane=self.focal_plane, 
                    chip=chip, 
                    filt=filt,  
458
459
                    exptime=pointing.exp_time,
                    pointing_type=pointing.pointing_type,
Fang Yuedong's avatar
Fang Yuedong committed
460
461
462
                    pointing_ID=pointing_ID,  
                    subdir=sub_img_dir,
                    prefix=prefix)
Fang Yuedong's avatar
Fang Yuedong committed
463
                chip_output.Log_info("running pointing#%d, chip#%d, at PID#%d..."%(pointing_ID, chip.chipID, pid))
464
                self.run_one_chip(
Fang Yuedong's avatar
Fang Yuedong committed
465
466
467
                    chip=chip, 
                    filt=filt, 
                    chip_output=chip_output, 
Fang Yuedong's avatar
Fang Yuedong committed
468
469
                    pointing=pointing)
                chip_output.Log_info("finished running chip#%d..."%(chip.chipID))
470
471
                for handler in chip_output.logger.handlers[:]:
                    chip_output.logger.removeHandler(handler)
Fang Yuedong's avatar
Fang Yuedong committed
472
                gc.collect()