Skip to content
GitLab
Commit 364e4c76 authored by xin's avatar xin
Browse files

local wcs sls

parent 9129df8f
Hide whitespace changes
Inline Side-by-side
ObservationSim/MockObject/MockObject.py
View file @ 364e4c76
...@@ -197,29 +197,46 @@ class MockObject(object): ...@@ -197,29 +197,46 @@ class MockObject(object):
star = star.withFlux(nphotons) star = star.withFlux(nphotons)
star = galsim.Convolve(psf, star) star = galsim.Convolve(psf, star)
stamp = star.drawImage(wcs=self.localWCS, method='phot', offset=self.offset, save_photons=True) # stamp = star.drawImage(wcs=self.localWCS, method='phot', offset=self.offset, save_photons=True)
# xmax = max(xmax, stamp.xmax)
# ymax = max(ymax, stamp.ymax)
# photons = stamp.photons
# photons.x += self.x_nominal
# photons.y += self.y_nominal
# photons_list.append(photons)
stamp = star.drawImage(wcs=real_wcs_local, method='phot', offset=offset, save_photons=True)
xmax = max(xmax, stamp.xmax) xmax = max(xmax, stamp.xmax)
ymax = max(ymax, stamp.ymax) ymax = max(ymax, stamp.ymax)
photons = stamp.photons photons = stamp.photons
photons.x += self.x_nominal photons.x += x_nominal
photons.y += self.y_nominal photons.y += y_nominal
photons_list.append(photons) photons_list.append(photons)
stamp = galsim.ImageF(int(xmax * 1.1), int(ymax * 1.1))
stamp.wcs = real_wcs_local
stamp.setCenter(x_nominal, y_nominal)
bounds = stamp.bounds & galsim.BoundsI(0, chip.npix_x-1, 0, chip.npix_y-1)
chip.img.setOrigin(0, 0)
stamp[bounds] = chip.img[bounds]
for i in range(len(photons_list)):
if i == 0:
chip.sensor.accumulate(photons_list[i], stamp)
else:
chip.sensor.accumulate(photons_list[i], stamp, resume=True)
chip.img[bounds] = stamp[bounds]
# stamp = star.drawImage(wcs=real_wcs_local, method='phot', offset=offset, save_photons=True) chip.img.setOrigin(chip.bound.xmin, chip.bound.ymin)
# xmax = max(xmax, stamp.xmax) # Test stamp size
# ymax = max(ymax, stamp.ymax) # print(xmax, ymax)
# photons = stamp.photons
# photons.x += x_nominal
# photons.y += y_nominal
# photons_list.append(photons)
# stamp = galsim.ImageF(int(xmax * 1.1), int(ymax * 1.1)) # stamp = galsim.ImageF(int(xmax*1.1), int(ymax*1.1))
# stamp.wcs = real_wcs_local # stamp.wcs = self.localWCS
# stamp.setCenter(x_nominal, y_nominal) # stamp.setCenter(self.x_nominal, self.y_nominal)
# bounds = stamp.bounds & galsim.BoundsI(0, chip.npix_x-1, 0, chip.npix_y-1) # bounds = stamp.bounds & chip.img.bounds
# chip.img.setOrigin(0, 0)
# stamp[bounds] = chip.img[bounds] # stamp[bounds] = chip.img[bounds]
# for i in range(len(photons_list)): # for i in range(len(photons_list)):
# if i == 0: # if i == 0:
...@@ -228,23 +245,6 @@ class MockObject(object): ...@@ -228,23 +245,6 @@ class MockObject(object):
# chip.sensor.accumulate(photons_list[i], stamp, resume=True) # chip.sensor.accumulate(photons_list[i], stamp, resume=True)
# #
# chip.img[bounds] = stamp[bounds] # chip.img[bounds] = stamp[bounds]
#
# chip.img.setOrigin(chip.bound.xmin, chip.bound.ymin)
# Test stamp size
# print(xmax, ymax)
stamp = galsim.ImageF(int(xmax*1.1), int(ymax*1.1))
stamp.wcs = self.localWCS
stamp.setCenter(self.x_nominal, self.y_nominal)
bounds = stamp.bounds & chip.img.bounds
stamp[bounds] = chip.img[bounds]
for i in range(len(photons_list)):
if i == 0:
chip.sensor.accumulate(photons_list[i], stamp)
else:
chip.sensor.accumulate(photons_list[i], stamp, resume=True)
chip.img[bounds] = stamp[bounds]
# print(chip.img.array.sum()) # print(chip.img.array.sum())
# print("nphotons_sum = ", nphotons_sum) # print("nphotons_sum = ", nphotons_sum)
del photons_list del photons_list
...@@ -252,7 +252,7 @@ class MockObject(object): ...@@ -252,7 +252,7 @@ class MockObject(object):
return True, pos_shear return True, pos_shear
def addSLStoChipImage(self, sdp = None, chip = None, xOrderSigPlus = None): def addSLStoChipImage(self, sdp = None, chip = None, xOrderSigPlus = None, local_wcs=None):
spec_orders = sdp.compute_spec_orders() spec_orders = sdp.compute_spec_orders()
for k, v in spec_orders.items(): for k, v in spec_orders.items():
img_s = v[0] img_s = v[0]
...@@ -267,15 +267,17 @@ class MockObject(object): ...@@ -267,15 +267,17 @@ class MockObject(object):
xlen_imf = int(specImg.xmax - specImg.xmin + 1) xlen_imf = int(specImg.xmax - specImg.xmin + 1)
ylen_imf = int(specImg.ymax - specImg.ymin + 1) ylen_imf = int(specImg.ymax - specImg.ymin + 1)
stamp = galsim.ImageF(xlen_imf, ylen_imf) stamp = galsim.ImageF(xlen_imf, ylen_imf)
stamp.wcs = self.localWCS stamp.wcs = local_wcs
stamp.setOrigin(origin_order_x, origin_order_y) stamp.setOrigin(origin_order_x, origin_order_y)
bounds = stamp.bounds & chip.img.bounds bounds = stamp.bounds & galsim.BoundsI(0, chip.npix_x-1, 0, chip.npix_y-1)
if bounds.area() == 0: if bounds.area() == 0:
continue continue
chip.img.setOrigin(0, 0)
stamp[bounds] = chip.img[bounds] stamp[bounds] = chip.img[bounds]
chip.sensor.accumulate(photons, stamp) chip.sensor.accumulate(photons, stamp)
chip.img[bounds] = stamp[bounds] chip.img[bounds] = stamp[bounds]
chip.img.setOrigin(chip.bound.xmin, chip.bound.ymin)
del stamp del stamp
del spec_orders del spec_orders
...@@ -301,9 +303,20 @@ class MockObject(object): ...@@ -301,9 +303,20 @@ class MockObject(object):
normalSED = Table(np.array([self.sed['WAVELENGTH'], self.sed['FLUX'] * sedNormFactor]).T, normalSED = Table(np.array([self.sed['WAVELENGTH'], self.sed['FLUX'] * sedNormFactor]).T,
names=('WAVELENGTH', 'FLUX')) names=('WAVELENGTH', 'FLUX'))
real_pos = self.getRealPos(chip.img, global_x=self.posImg.x, global_y=self.posImg.y, img_real_wcs=self.real_wcs)
x, y = real_pos.x + 0.5, real_pos.y + 0.5
x_nominal = int(np.floor(x + 0.5))
y_nominal = int(np.floor(y + 0.5))
dx = x - x_nominal
dy = y - y_nominal
offset = galsim.PositionD(dx, dy)
real_wcs_local = self.real_wcs.local(real_pos)
xOrderSigPlus = {'A':1.3909419820029296,'B':1.4760376591236062,'C':4.035447379743442,'D':5.5684364343742825,'E':16.260021029735388} xOrderSigPlus = {'A':1.3909419820029296,'B':1.4760376591236062,'C':4.035447379743442,'D':5.5684364343742825,'E':16.260021029735388}
grating_split_pos_chip = chip.bound.xmin + grating_split_pos grating_split_pos_chip = 0 + grating_split_pos
for i in range(len(bandpass_list)): for i in range(len(bandpass_list)):
bandpass = bandpass_list[i] bandpass = bandpass_list[i]
psf, pos_shear = psf_model.get_PSF(chip=chip, pos_img=pos_img, bandpass=bandpass, folding_threshold=folding_threshold) psf, pos_shear = psf_model.get_PSF(chip=chip, pos_img=pos_img, bandpass=bandpass, folding_threshold=folding_threshold)
...@@ -312,8 +325,8 @@ class MockObject(object): ...@@ -312,8 +325,8 @@ class MockObject(object):
star = galsim.Convolve(psf, star) star = galsim.Convolve(psf, star)
starImg = star.drawImage(nx=100, ny=100, wcs=self.localWCS) starImg = star.drawImage(nx=100, ny=100, wcs=self.localWCS)
origin_star = [self.y_nominal - (starImg.center.y - starImg.ymin), origin_star = [y_nominal - (starImg.center.y - starImg.ymin),
self.x_nominal - (starImg.center.x - starImg.xmin)] x_nominal - (starImg.center.x - starImg.xmin)]
gal_origin = [origin_star[0], origin_star[1]] gal_origin = [origin_star[0], origin_star[1]]
gal_end = [origin_star[0] + starImg.array.shape[0] - 1, origin_star[1] + starImg.array.shape[1] - 1] gal_end = [origin_star[0] + starImg.array.shape[0] - 1, origin_star[1] + starImg.array.shape[1] - 1]
...@@ -325,8 +338,8 @@ class MockObject(object): ...@@ -325,8 +338,8 @@ class MockObject(object):
subImg_p1 = starImg.array[:, 0:subSlitPos] subImg_p1 = starImg.array[:, 0:subSlitPos]
star_p1 = galsim.Image(subImg_p1) star_p1 = galsim.Image(subImg_p1)
origin_p1 = origin_star origin_p1 = origin_star
xcenter_p1 = min(self.x_nominal,grating_split_pos_chip-1) - chip.bound.xmin xcenter_p1 = min(x_nominal,grating_split_pos_chip-1) - 0
ycenter_p1 = self.y_nominal-chip.bound.ymin ycenter_p1 = y_nominal-0
sdp_p1 = SpecDisperser(orig_img=star_p1, xcenter=xcenter_p1, sdp_p1 = SpecDisperser(orig_img=star_p1, xcenter=xcenter_p1,
ycenter=ycenter_p1, origin=origin_p1, ycenter=ycenter_p1, origin=origin_p1,
...@@ -335,13 +348,13 @@ class MockObject(object): ...@@ -335,13 +348,13 @@ class MockObject(object):
conf=chip.sls_conf[0], conf=chip.sls_conf[0],
isAlongY=0) isAlongY=0)
self.addSLStoChipImage(sdp=sdp_p1, chip=chip, xOrderSigPlus = xOrderSigPlus) self.addSLStoChipImage(sdp=sdp_p1, chip=chip, xOrderSigPlus = xOrderSigPlus, local_wcs=real_wcs_local)
subImg_p2 = starImg.array[:, subSlitPos+1:starImg.array.shape[1]] subImg_p2 = starImg.array[:, subSlitPos+1:starImg.array.shape[1]]
star_p2 = galsim.Image(subImg_p2) star_p2 = galsim.Image(subImg_p2)
origin_p2 = [origin_star[0], grating_split_pos_chip] origin_p2 = [origin_star[0], grating_split_pos_chip]
xcenter_p2 = max(self.x_nominal, grating_split_pos_chip - 1) - chip.bound.xmin xcenter_p2 = max(x_nominal, grating_split_pos_chip - 1) - 0
ycenter_p2 = self.y_nominal - chip.bound.ymin ycenter_p2 = y_nominal - 0
sdp_p2 = SpecDisperser(orig_img=star_p2, xcenter=xcenter_p2, sdp_p2 = SpecDisperser(orig_img=star_p2, xcenter=xcenter_p2,
ycenter=ycenter_p2, origin=origin_p2, ycenter=ycenter_p2, origin=origin_p2,
...@@ -350,27 +363,27 @@ class MockObject(object): ...@@ -350,27 +363,27 @@ class MockObject(object):
conf=chip.sls_conf[1], conf=chip.sls_conf[1],
isAlongY=0) isAlongY=0)
self.addSLStoChipImage(sdp=sdp_p2, chip=chip,xOrderSigPlus = xOrderSigPlus) self.addSLStoChipImage(sdp=sdp_p2, chip=chip,xOrderSigPlus = xOrderSigPlus, local_wcs=real_wcs_local)
del sdp_p1 del sdp_p1
del sdp_p2 del sdp_p2
elif grating_split_pos_chip<=gal_origin[1]: elif grating_split_pos_chip<=gal_origin[1]:
sdp = SpecDisperser(orig_img=starImg, xcenter=self.x_nominal - chip.bound.xmin, sdp = SpecDisperser(orig_img=starImg, xcenter=x_nominal - 0,
ycenter=self.y_nominal - chip.bound.ymin, origin=origin_star, ycenter=y_nominal - 0, origin=origin_star,
tar_spec=normalSED, tar_spec=normalSED,
band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10, band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10,
conf=chip.sls_conf[1], conf=chip.sls_conf[1],
isAlongY=0) isAlongY=0)
self.addSLStoChipImage(sdp=sdp, chip=chip, xOrderSigPlus = xOrderSigPlus) self.addSLStoChipImage(sdp=sdp, chip=chip, xOrderSigPlus = xOrderSigPlus, local_wcs=real_wcs_local)
del sdp del sdp
elif grating_split_pos_chip>=gal_end[1]: elif grating_split_pos_chip>=gal_end[1]:
sdp = SpecDisperser(orig_img=starImg, xcenter=self.x_nominal - chip.bound.xmin, sdp = SpecDisperser(orig_img=starImg, xcenter=x_nominal - 0,
ycenter=self.y_nominal - chip.bound.ymin, origin=origin_star, ycenter=y_nominal - 0, origin=origin_star,
tar_spec=normalSED, tar_spec=normalSED,
band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10, band_start=bandpass.blue_limit * 10, band_end=bandpass.red_limit * 10,
conf=chip.sls_conf[0], conf=chip.sls_conf[0],
isAlongY=0) isAlongY=0)
self.addSLStoChipImage(sdp=sdp, chip=chip, xOrderSigPlus = xOrderSigPlus) self.addSLStoChipImage(sdp=sdp, chip=chip, xOrderSigPlus = xOrderSigPlus, local_wcs=real_wcs_local)
del sdp del sdp
del psf del psf
return True, pos_shear return True, pos_shear
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment