Python QubicAcquisition.get_hwp_operator方法代碼示例

# 需要導入模塊: from qubic import QubicAcquisition [as 別名]
# 或者: from qubic.QubicAcquisition import get_hwp_operator [as 別名]
def tod2map(tod,pointing,instrument_in,detector_list=False,disp=True,kmax=2,displaytime=False):
    t0=time.time()
    #### Detectors
    mask_packed = np.ones(len(instrument_in.detector.packed), bool)
    if detector_list:
        mask_packed[detector_list] = False
        mask_unpacked = instrument_in.unpack(mask_packed)
        instrument = QubicInstrument('monochromatic', removed=mask_unpacked,nside=instrument_in.sky.nside)
    else:
        instrument = instrument_in
        
    #### Observations
    obs = QubicAcquisition(instrument, pointing)
    projection = obs.get_projection_peak_operator(kmax=kmax)
    coverage = projection.pT1()
    mask = coverage == 0
    projection = pack_projection_inplace(projection, mask)
    hwp = obs.get_hwp_operator()
    polgrid = DenseOperator([[0.5, 0.5, 0],
                             [0.5,-0.5, 0]])
    H = polgrid * hwp * projection
    preconditioner = DiagonalOperator(1/coverage[~mask], broadcast='rightward')
    solution = pcg(H.T * H, H.T(tod), M=preconditioner, disp=disp, tol=1e-3)
    output_map = unpack(solution['x'], mask)
    t1=time.time()
    if displaytime: print(' Map done in {0:.4f} seconds'.format(t1-t0))
    return output_map,coverage

# 需要導入模塊: from qubic import QubicAcquisition [as 別名]
# 或者: from qubic.QubicAcquisition import get_hwp_operator [as 別名]
def map2tod(maps,pointing,instrument_in,detector_list=False,kmax=2):
    #### Detectors
    mask_packed = np.ones(len(instrument_in.detector.packed), bool)
    if detector_list:
        mask_packed[detector_list] = False
        mask_unpacked = instrument_in.unpack(mask_packed)
        instrument = QubicInstrument('monochromatic', removed=mask_unpacked,nside=instrument_in.sky.nside)
    else:
        instrument = instrument_in
        
    #### Observations
    obs = QubicAcquisition(instrument, pointing)
    #C = obs.get_convolution_peak_operator()
    #convmaps=np.transpose(np.array([C(maps[:,0]),C(maps[:,1]),C(maps[:,2])]))
    projection = obs.get_projection_peak_operator(kmax=kmax)
    coverage = projection.pT1()
    mask = coverage == 0
    projection = pack_projection_inplace(projection, mask)
    hwp = obs.get_hwp_operator()
    polgrid = DenseOperator([[0.5, 0.5, 0],
                             [0.5,-0.5, 0]])
    #H = polgrid * hwp * projection * C
    H = polgrid * hwp * projection
    x1 = pack(maps, mask)
    y = H(x1)
    #return y,convmaps
    return y

# 需要導入模塊: from qubic import QubicAcquisition [as 別名]
# 或者: from qubic.QubicAcquisition import get_hwp_operator [as 別名]
    angspeed_psi, maxpsi)

tods = {}
pTxs = {}
pT1s = {}
pTx_pT1s = {}
cbiks = {}
outputs = {}
kinds = ['I', 'QU', 'IQU']
input_maps = [I,
              np.array([Q, U]).T,
              np.array([I, Q, U]).T]
for kind, input_map in zip(kinds, input_maps):
    acq = QubicAcquisition(150, pointings, kind=kind)
    P = acq.get_projection_operator()
    W = acq.get_hwp_operator()
    H = W * P
    coverage = P.pT1()
    tod = H(input_map)
    tods[kind] = tod
    pTxs[kind] = H.T(tod)[coverage > 0]
    if kind != 'QU':
        pTx_pT1 = P.pTx_pT1(tod)
        pTx_pT1s[kind] = pTx_pT1[0][coverage > 0], pTx_pT1[1][coverage > 0]
    cbik = P.canonical_basis_in_kernel()
    mask = coverage > 10
    P = P.restrict(mask, inplace=True)
    unpack = UnpackOperator(mask, broadcast='rightward')
    x0 = unpack.T(input_map)

    M = DiagonalOperator(1 / coverage[mask], broadcast='rightward')

# 需要導入模塊: from qubic import QubicAcquisition [as 別名]
# 或者: from qubic.QubicAcquisition import get_hwp_operator [as 別名]
pointing.angle_hwp = np.random.random_integers(0, 7, pointing.size) * 22.5
ntimes = len(pointing)

# get instrument model with only one detector
idetector = 231
instrument = QubicInstrument('monochromatic')
instrument.plot()
mask_packed = np.ones(len(instrument.detector.packed), bool)
mask_packed[idetector] = False
mask_unpacked = instrument.unpack(mask_packed)
instrument = QubicInstrument('monochromatic', removed=mask_unpacked)

obs = QubicAcquisition(instrument, pointing)
convolution = obs.get_convolution_peak_operator()
projection = obs.get_projection_peak_operator(kmax=0)
hwp = obs.get_hwp_operator()
polarizer = obs.get_polarizer_operator()
coverage = projection.pT1()
mask = coverage > 0
projection.restrict(mask)
pack = PackOperator(mask, broadcast='rightward')
ra,dec=pointings_modbyJC._hor2equ(pointing.azimuth, pointing.elevation, pointing.latitude, pointing.time/3600)
ns_scan = int(delta_az / angspeed / 0.1)
ns_tot = ns_scan * 2
ichunk = (np.arange(len(pointing)) / ns_tot / nsweeps_el).astype(numpy.int64)
isweep = ((np.arange(len(pointing)) / ns_tot).astype(numpy.int64)) % nsweeps_el
clf()
subplot(2,2,1)
plot(pointing.azimuth,pointing.elevation)
e2g = SphericalEquatorial2GalacticOperator(degrees=True)
center = e2g([racenter, deccenter])