ObservationSim.py 22.3 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
18
19
from ObservationSim.PSF import PSFGauss, FieldDistortion, PSFInterp
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, shear_cat_file=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":
Fang Yuedong's avatar
Fang Yuedong committed
66
            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
67
        else:
Fang Yuedong's avatar
Fang Yuedong committed
68
            chip_output.Log_error("unrecognized PSF model type!!", flush=True)
Fang Yuedong's avatar
Fang Yuedong committed
69

Fang Yuedong's avatar
Fang Yuedong committed
70
        # Figure out shear fields
71
        self.g1_field, self.g2_field, self.nshear = get_shear_field(config=self.config, shear_cat_file=shear_cat_file)
Fang Yuedong's avatar
Fang Yuedong committed
72

73
74
        # Apply astrometric simulation for pointing
        if self.config["obs_setting"]["enable_astrometric_model"]:
Fang Yuedong's avatar
Fang Yuedong committed
75
            dt = datetime.utcfromtimestamp(pointing.timestamp)
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
            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
92
                input_epoch="J2000",
93
94
95
96
97
98
99
100
                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
101
102
        # Get WCS for the focal plane
        if wcs_fp == None:
103
            wcs_fp = self.focal_plane.getTanWCS(ra_cen, dec_cen, pointing.img_pa, chip.pix_scale)
Fang Yuedong's avatar
Fang Yuedong committed
104
105
106
107
108

        # 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
109

Zhang Xin's avatar
Zhang Xin committed
110
111
        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]))
112

Zhang Xin's avatar
Zhang Xin committed
113
        print("========================sky pix========================")
114
115
        print(filt.sky_background)

Fang Yuedong's avatar
Fang Yuedong committed
116
117
118
        if chip.survey_type == "photometric":
            sky_map = None
        elif chip.survey_type == "spectroscopic":
xin's avatar
xin committed
119
            # chip.loadSLSFLATCUBE(flat_fn='flat_cube.fits')
120
121
            flat_normal = np.ones_like(chip.img.array)
            if self.config["ins_effects"]["flat_fielding"] == True:
Fang Yuedong's avatar
Fang Yuedong committed
122
                chip_output.Log_info("SLS flat preprocess,CHIP %d : Creating and applying Flat-Fielding"%chip.chipID)
123
                msg = str(chip.img.bounds)
Fang Yuedong's avatar
Fang Yuedong committed
124
                chip_output.Log_info(msg)
125
126
127
128
129
                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
130
                chip_output.Log_info("SLS flat preprocess,CHIP %d : Apply shutter effect"%chip.chipID)
131
132
133
134
                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
135
            sky_map = calculateSkyMap_split_g(
136
137
138
139
140
141
                skyMap=flat_normal,
                blueLimit=filt.blue_limit,
                redLimit=filt.red_limit,
                conf=chip.sls_conf,
                pixelSize=chip.pix_scale,
                isAlongY=0,
142
143
                flat_cube=chip.flat_cube, zoldial_spec = filt.zodical_spec)
            sky_map = sky_map+filt.sky_background
144
            del flat_normal
Fang Yuedong's avatar
Fang Yuedong committed
145

146
        if pointing.pointing_type == 'MS':
Fang Yuedong's avatar
Fang Yuedong committed
147
            # Load catalogues and templates
Fang Yuedong's avatar
Fang Yuedong committed
148
            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)
149
            chip_output.create_output_file()
Fang Yuedong's avatar
Fang Yuedong committed
150
151
            self.nobj = len(self.cat.objs)

Fang Yuedong's avatar
Fang Yuedong committed
152
153
154
155
156
157
158
159
160
            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
161
            if self.config["ins_effects"]["field_dist"] == True:
162
                self.fd_model = FieldDistortion(chip=chip, img_rot=pointing.img_pa.deg)
Fang Yuedong's avatar
Fang Yuedong committed
163
164
165
            else:
                self.fd_model = None

Fang Yuedong's avatar
Fang Yuedong committed
166
167
168
169
            # Loop over objects
            missed_obj = 0
            bright_obj = 0
            dim_obj = 0
Fang Yuedong's avatar
Fang Yuedong committed
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
 
            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,
186
187
188
189
                extName='SCI',
                timestamp = pointing.timestamp,
                exptime = pointing.exp_time,
                readoutTime = 40.)
190
191
            
            chip_wcs = galsim.FitsWCS(header=h_ext)
Fang Yuedong's avatar
Fang Yuedong committed
192

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

Fang Yuedong's avatar
Fang Yuedong committed
199
200
                obj = self.cat.objs[j]

201
                # load and convert SED; also caculate object's magnitude in all CSST bands
202
203
204
                try:
                    sed_data = self.cat.load_sed(obj)
                    norm_filt = self.cat.load_norm_filt(obj)
Fang Yuedong's avatar
Fang Yuedong committed
205
                    obj.sed, obj.param["mag_%s"%filt.filter_type.lower()], obj.param["flux_%s"%filt.filter_type.lower()] = self.cat.convert_sed(
206
207
208
209
210
                        mag=obj.param["mag_use_normal"],
                        sed=sed_data,
                        target_filt=filt, 
                        norm_filt=norm_filt,
                    )
Fang Yuedong's avatar
Fang Yuedong committed
211
                    _, 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
212
213
214
215
216
                        mag=obj.param["mag_use_normal"],
                        sed=sed_data,
                        target_filt=cut_filter, 
                        norm_filt=norm_filt,
                    )
217
                except Exception as e:
Fang Yuedong's avatar
Fang Yuedong committed
218
                    traceback.print_exc()
Fang Yuedong's avatar
Fang Yuedong committed
219
                    chip_output.Log_error(e)
220
                    continue
221
                
Fang Yuedong's avatar
Fang Yuedong committed
222
                # [TODO] Testing
223
                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
224
225

                # Exclude very bright/dim objects (for now)
226
227
228
                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
229
230
231
232
                    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
233
234
235
                if filt.is_too_dim(
                    mag=obj.getMagFilter(filt),
                    margin=self.config["obs_setting"]["mag_lim_margin"]):
Fang Yuedong's avatar
Fang Yuedong committed
236
                    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
237
                    dim_obj += 1
Fang Yuedong's avatar
Fang Yuedong committed
238
239
240
                    obj.unload_SED()
                    continue

Fang Yuedong's avatar
Fang Yuedong committed
241
                # Get corresponding shear values
Fang Yuedong's avatar
Fang Yuedong committed
242
                if self.config["shear_setting"]["shear_type"] == "constant":
Fang Yuedong's avatar
Fang Yuedong committed
243
                    if obj.type == 'star':
244
                        obj.g1, obj.g2 = 0., 0.
Fang Yuedong's avatar
Fang Yuedong committed
245
                    else:
246
                        obj.g1, obj.g2 = self.g1_field, self.g2_field
Fang Yuedong's avatar
Fang Yuedong committed
247
248
                elif self.config["shear_setting"]["shear_type"] == "extra":
                    try:
Fang Yuedong's avatar
Fang Yuedong committed
249
                        # [TODO]: every object with individual shear from input catalog(s)
250
                        obj.g1, obj.g2 = self.g1_field[j], self.g2_field[j]
Fang Yuedong's avatar
Fang Yuedong committed
251
                    except:
Fang Yuedong's avatar
Fang Yuedong committed
252
                        chip_output.Log_error("failed to load external shear.")
253
254
255
                elif self.config["shear_setting"]["shear_type"] == "catalog":
                    pass
                else:
Fang Yuedong's avatar
Fang Yuedong committed
256
                    chip_output.Log_error("Unknown shear input")
257
                    raise ValueError("Unknown shear input")
Fang Yuedong's avatar
Fang Yuedong committed
258
259

                # Get position of object on the focal plane
260
                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
261
262
263
264

                # [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
265
                if pos_img.x == -1 or pos_img.y == -1:
Fang Yuedong's avatar
Fang Yuedong committed
266
267
                    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
268
269
270
                    missed_obj += 1
                    obj.unload_SED()
                    continue
Fang Yuedong's avatar
Fang Yuedong committed
271

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

Fang Yuedong's avatar
Fang Yuedong committed
291
                    elif chip.survey_type == "spectroscopic" and not self.config["run_option"]["out_cat_only"]:
Fang Yuedong's avatar
Fang Yuedong committed
292
293
294
295
296
297
298
                        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, 
299
300
                            g1=obj.g1, 
                            g2=obj.g2, 
301
                            exptime=pointing.exp_time,
Fang Yuedong's avatar
Fang Yuedong committed
302
303
304
305
                            normFilter=norm_filt,
                            fd_shear=fd_shear)

                    if isUpdated == 1:
Fang Yuedong's avatar
Fang Yuedong committed
306
                        # TODO: add up stats
307
                        chip_output.cat_add_obj(obj, pos_img, pos_shear)
Fang Yuedong's avatar
Fang Yuedong committed
308
                        pass
Fang Yuedong's avatar
Fang Yuedong committed
309
310
311
                    elif isUpdated == 0:
                        missed_obj += 1
                        chip_output.Log_error("Objected missed: %s"%(obj.id))
Fang Yuedong's avatar
Fang Yuedong committed
312
                    else:
Fang Yuedong's avatar
Fang Yuedong committed
313
                        chip_output.Log_error("Draw error, object omitted: %s"%(obj.id))
Fang Yuedong's avatar
Fang Yuedong committed
314
315
                        continue
                except Exception as e:
Fang Yuedong's avatar
Fang Yuedong committed
316
                    traceback.print_exc()
Fang Yuedong's avatar
Fang Yuedong committed
317
318
319
320
321
322
323
                    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
324
325
326
                # Unload SED:
                obj.unload_SED()
                del obj
Fang Yuedong's avatar
Fang Yuedong committed
327
                gc.collect()
Fang Yuedong's avatar
Fang Yuedong committed
328

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

Fang Yuedong's avatar
Fang Yuedong committed
333
        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
334
335
336

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

        if not self.config["run_option"]["out_cat_only"]:
            chip.img = chip.addEffects(
                config=self.config, 
                img=chip.img, 
                chip_output=chip_output, 
344
                filt=filt,
345
346
347
                ra_cen=pointing.ra, 
                dec_cen=pointing.dec,
                img_rot=pointing.img_pa,
Fang Yuedong's avatar
Fang Yuedong committed
348
                exptime=pointing.exp_time,
349
350
351
352
353
354
355
                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
356
                datetime_obs = datetime.utcfromtimestamp(pointing.timestamp)
357
358
359
360
361
362
363
364
                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
365
                    pixel_scale=chip.pix_scale,
366
367
368
369
370
                    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
371
372
                    sat_vel=[pointing.sat_vx, pointing.sat_vy, pointing.sat_vz],
                    chip_name=str(chip.chipID).rjust(2, '0'))
373
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,
                    readoutTime=40.)
Fang Yuedong's avatar
Fang Yuedong committed
392
                
393
394
                chip.img = galsim.Image(chip.img.array, dtype=np.uint16)
                hdu1 = fits.PrimaryHDU(header=h_prim)
Fang Yuedong's avatar
Fang Yuedong committed
395
                hdu1.add_checksum()
396
397
                hdu1.header.comments['CHECKSUM'] = 'HDU checksum'
                hdu1.header.comments['DATASUM'] = 'data unit checksum'
398
                hdu2 = fits.ImageHDU(chip.img.array, header=h_ext)
Fang Yuedong's avatar
Fang Yuedong committed
399
                hdu2.add_checksum()
400
401
402
                hdu2.header.comments['XTENSION'] = 'extension type'
                hdu2.header.comments['CHECKSUM'] = 'HDU checksum'
                hdu2.header.comments['DATASUM'] = 'data unit checksum'
403
404
405
                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
406
407
408
409

                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
410
411
        del chip.img

Fang Yuedong's avatar
Fang Yuedong committed
412
        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
413

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

Fang Yuedong's avatar
Fang Yuedong committed
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
        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
435

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

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