Commit d9fd79bc authored by Yan Zhaojun's avatar Yan Zhaojun
Browse files

update

parent aa97d21f
......@@ -821,16 +821,16 @@ def psfMaker_IDW(px, py, PSFMat, cen_col, cen_row, dn=5, IDWindex=3, OnlyNeighbo
dist = np.sqrt((ref_col - cen_col[idx])
** 2 + (ref_row - cen_row[idx])**2)
if IDWindex == 1:
psfWeight[ipsf] = dist
if IDWindex == 2:
psfWeight[ipsf] = dist**2
if IDWindex == 3:
psfWeight[ipsf] = dist**3
if IDWindex == 4:
psfWeight[ipsf] = dist**4
if IDWindex == 5:
psfWeight[ipsf] = dist**5
#if IDWindex == 1:
psfWeight[ipsf] = dist**IDWindex
# if IDWindex == 2:
# psfWeight[ipsf] = dist**2
# if IDWindex == 3:
# psfWeight[ipsf] = dist**3
# if IDWindex == 4:
# psfWeight[ipsf] = dist**4
# if IDWindex == 5:
# psfWeight[ipsf] = dist**5
psfWeight[ipsf] = max(psfWeight[ipsf], minimum_psf_weight)
psfWeight[ipsf] = 1./psfWeight[ipsf]
......@@ -1051,7 +1051,7 @@ class MCIsimulator():
(self.information['ysize'], self.information['xsize']), dtype=float)
return
###############################################################################
##########################################IDWindex#####################################
##############################################################################
def _loadGhostModel(self):
"""
......@@ -1822,11 +1822,11 @@ class MCIsimulator():
# ###### use SED_code to generate star SED #######
if get_file_extension(starcat) == '.fits':
umag = self.star['gmag'][j]
gmag = self.star['gmag'][j]
rmag = self.star['rmag'][j]
imag = self.star['imag'][j]
zmag = self.star['zmag'][j]
umag = self.star['gmag'][j]; gmag = self.star['gmag'][j]; rmag = self.star['rmag'][j]; imag = self.star['imag'][j]; zmag = self.star['zmag'][j];
# SED of j-th star
else:
......@@ -1926,7 +1926,7 @@ class MCIsimulator():
psf[ch] = ndimage.rotate(
temp, theta.deg, order=1, reshape=True)
else:
psf[ch] = temp
psf[ch] = tempIDWindex
conv = psf[ch]
conv = conv/conv.sum()
......@@ -2205,7 +2205,7 @@ class MCIsimulator():
# 'earthshine': earthshine_ifs})
self.zodiacal_wave = wave_mci # in A
self.zodiacal_flux = zodi_mci
self.zodiacal_flux = zodi_mciIDWindex
self.earthshine_wave = wave_mci # A
self.earthshine_flux = earthshine_mci
......@@ -3312,42 +3312,42 @@ class MCIsimulator():
return
##############################################################################
def img_fits_save(self, data, filename):
####
ofd_g = fits.PrimaryHDU()
# World coordinate system and related parameters #####
hdu_g = fits.ImageHDU(data)
hdu_g.header['WCSAXES'] = (np.int16(2), 'number of WCS axes')
hdu_g.header['CRPIX1'] = (
round(float(self.information['CRPIX1']), 1), 'x-coordinate of reference pixel')
hdu_g.header['CRPIX2'] = (
round(float(self.information['CRPIX2']), 1), 'y-coordinate of reference pixel')
hdu_g.header['CRVAL1'] = (
float(self.information['CRVAL1']), 'first axis value at reference pixel')
hdu_g.header['CRVAL2'] = (
float(self.information['CRVAL2']), 'second axis value at reference pixel')
hdu_g.header['CTYPE1'] = (
'RA---TAN', 'the coordinate type for the first axis')
hdu_g.header['CTYPE2'] = (
'DEC--TAN', 'the coordinate type for the second axis')
hdu_g.header['CD1_1'] = (
float(self.information['CD1_1']), 'partial of first axis coordinate w.r.t. x')
hdu_g.header['CD1_2'] = (
float(self.information['CD1_2']), 'partial of first axis coordinate w.r.t. y')
hdu_g.header['CD2_1'] = (
float(self.information['CD2_1']), 'partial of second axis coordinate w.r.t. x')
hdu_g.header['CD2_2'] = (
float(self.information['CD2_2']), 'partial of second axis coordinate w.r.t. y')
hdulist_g = fits.HDUList([ofd_g, hdu_g])
file_g = self.result_path+'/ori_Sky/'+filename + '.fits'
hdulist_g.writeto(file_g, overwrite=True)
return
# def img_fits_save(self, data, filename):
# ####
# ofd_g = fits.PrimaryHDU()
# # World coordinate system and related parameters #####
# hdu_g = fits.ImageHDU(data)
# hdu_g.header['WCSAXES'] = (np.int16(2), 'number of WCS axes')
# hdu_g.header['CRPIX1'] = (
# round(float(self.information['CRPIX1']), 1), 'x-coordinate of reference pixel')
# hdu_g.header['CRPIX2'] = (
# round(float(self.information['CRPIX2']), 1), 'y-coordinate of reference pixel')
# hdu_g.header['CRVAL1'] = (
# float(self.information['CRVAL1']), 'first axis value at reference pixel')
# hdu_g.header['CRVAL2'] = (
# float(self.information['CRVAL2']), 'second axis value at reference pixel')
# hdu_g.header['CTYPE1'] = (
# 'RA---TAN', 'the coordinate type for the first axis')
# hdu_g.header['CTYPE2'] = (
# 'DEC--TAN', 'the coordinate type for the second axis')
# hdu_g.header['CD1_1'] = (
# float(self.information['CD1_1']), 'partial of first axis coordinate w.r.t. x')
# hdu_g.header['CD1_2'] = (
# float(self.information['CD1_2']), 'partial of first axis coordinate w.r.t. y')
# hdu_g.header['CD2_1'] = (
# float(self.information['CD2_1']), 'partial of second axis coordinate w.r.t. x')
# hdu_g.header['CD2_2'] = (
# float(self.information['CD2_2']), 'partial of second axis coordinate w.r.t. y')
# hdulist_g = fits.HDUList([ofd_g, hdu_g])
# file_g = self.result_path+'/ori_Sky/'+filename + '.fits'
# hdulist_g.writeto(file_g, overwrite=True)
# return
###########
########
##############################################################################
......@@ -3370,46 +3370,46 @@ class MCIsimulator():
else:
if direction == 'horizon':
self.information['bleding_direction']='Not_horizon'
# loop over each column, as bleeding is modelled column-wise
for i, column in enumerate(data): # select one solumnn
if column.max() <= self.information['fullwellcapacity']:
continue
sum = 0.
for j, value in enumerate(column):
# first round - from bottom to top (need to half the bleeding)
overload = value - self.information['fullwellcapacity']
if overload > 0.:
overload /= 2.
# self.image[j, i] -= overload
data[i, j] -= overload
sum += overload
elif sum > 0.:
if -overload > sum:
overload = -sum
# self.image[j, i] -= overload
data[i, j] -= overload
sum += overload
######################################################
for i, column in enumerate(data):
if column.max() <= self.information['fullwellcapacity']:
continue
sum = 0.
for j, value in enumerate(column[::-1]):
# second round - from top to bottom (bleeding was half'd already, so now full)
overload = value - self.information['fullwellcapacity']
if overload > 0.:
# self.image[-j-1, i] -= overload
data[i, -j-1] -= overload
sum += overload
elif sum > 0.:
if -overload > sum:
overload = -sum
# self.image[-j-1, i] -= overload
data[i, -j-1,] -= overload
sum += overload
# for i, column in enumerate(data): # select one solumnn
# if column.max() <= self.information['fullwellcapacity']:
# continue
# sum = 0.
# for j, value in enumerate(column):
# # first round - from bottom to top (need to half the bleeding)
# overload = value - self.information['fullwellcapacity']
# if overload > 0.:
# overload /= 2.
# # self.image[j, i] -= overload
# data[i, j] -= overload
# sum += overload
# elif sum > 0.:
# if -overload > sum:
# overload = -sum
# # self.image[j, i] -= overload
# data[i, j] -= overload
# sum += overload
# ######################################################
# for i, column in enumerate(data):
# if column.max() <= self.information['fullwellcapacity']:
# continue
# sum = 0.
# for j, value in enumerate(column[::-1]):
# # second round - from top to bottom (bleeding was half'd already, so now full)
# overload = value - self.information['fullwellcapacity']
# if overload > 0.:
# # self.image[-j-1, i] -= overload
# data[i, -j-1] -= overload
# sum += overload
# elif sum > 0.:
# if -overload > sum:
# overload = -sum
# # self.image[-j-1, i] -= overload
# data[i, -j-1,] -= overload
# sum += overload
else:
......
......@@ -22,12 +22,12 @@ def onOrbitObsPosition(path, input_ra_list, input_dec_list, input_pmra_list, inp
if isinstance(type(input_nstars), type(1)): # type(input_nstars) != type(1):
raise TypeError("Parameter 7 is not int!", input_nstars)
checkInputList(input_ra_list, input_nstars)
checkInputList(input_dec_list, input_nstars)
checkInputList(input_pmra_list, input_nstars)
checkInputList(input_pmdec_list, input_nstars)
checkInputList(input_rv_list, input_nstars)
checkInputList(input_parallax_list, input_nstars)
# checkInputList(input_ra_list, input_nstars)
# checkInputList(input_dec_list, input_nstars)
# checkInputList(input_pmra_list, input_nstars)
# checkInputList(input_pmdec_list, input_nstars)
# checkInputList(input_rv_list, input_nstars)
# checkInputList(input_parallax_list, input_nstars)
# if isinstance(type(input_x), type(1.1)): # type(input_x) != type(1.1):
# raise TypeError("Parameter 8 is not double!", input_x)
......
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