Page History

def calib_jungfrau_v2(det, evt, cmpars=(7,3,200,10), **kwa):
    """
    v2 - improving performance, reduce time and memory consumption, use peds-offset constants
    Returns calibrated jungfrau data

    - gets constants
    - gets raw data
    - evaluates (code - pedestal - offset)
    - applys common mode correction if turned on
    - apply gain factor

    Parameters

    - det (psana.Detector) - Detector object
    - evt (psana.Event)    - Event object
    - cmpars (tuple) - common mode parameters
        - cmpars[0] - algorithm # 7-for jungfrau
        - cmpars[1] - control bit-word 1-in rows, 2-in columns
        - cmpars[2] - maximal applied correction
    - **kwa - used here and passed to det.mask_v2 or det.mask_comb
      - nda_raw - if not None, substitutes evt.raw()
      - mbits - DEPRECATED parameter of the det.mask_comb(...)
      - mask - user defined mask passed as optional parameter
    """

    t00 = time()

    src = det.source # - src (psana.Source) - Source object

    nda_raw = kwa.get('nda_raw', None)
    loop_segs = kwa.get('loop_segs', False)
    arr = det.raw(evt) if nda_raw is None else nda_raw # shape:(<npanels>, 512, 1024) dtype:uint16
    if arr is None: return None

    t01 = time()

    detname = string_from_source(det.source)

    t02 = time()

    odc = cache.detcache_for_detname(detname)

    t03 = time()

    first_entry = odc is None
    if first_entry:
       odc = cache.add_detcache(det, evt)
       odc.cmps = det.common_mode(evt) if cmpars is None else cmpars
       odc.mask = det.mask_total(evt, **kwa)

    poff = odc.poff # 4d pedestals + offset shape:(3, 1, 512, 1024) dtype:float32
    gfac = odc.gfac # 4d gain factors evaluated form gains
    mask = odc.mask
    outa = odc.outa
    cmps = odc.cmps

    t04 = time()

    if first_entry:
        logger.info('\n  ====================== det.name: %s' % det.name\
                   +'\n  detname from source: %s' % string_from_source(det.source)\
                   +info_ndarr(arr,  '\n  calib_jungfrau arr ')\
                   +info_ndarr(poff, '\n  calib_jungfrau peds+off')\
                   +info_ndarr(gfac, '\n  calib_jungfrau gfac')\
                   +info_ndarr(mask, '\n  calib_jungfrau mask')\
                   +info_ndarr(outa, '\n  calib_jungfrau outa')\
                   +info_ndarr(cmps, '\n  calib_jungfrau common mode parameters ')
                   +'\n    loop over segments: %s' % loop_segs)

    if loop_segs:
      nsegs = arr.shape[0]
      shseg = arr.shape[-2:] # (512, 1024)
      for i in range(nsegs):
        arr1  = arr[i,:]
        mask1 = mask[i,:]
        gfac1 = gfac[:,i,:,:]
        poff1 = poff[:,i,:,:]
        arr1.shape  = (1,) + shseg
        mask1.shape = (1,) + shseg
        gfac1.shape = (3,1,) + shseg
        poff1.shape = (3,1,) + shseg
        out1, times = calib_jungfrau_single_panel(arr1, gfac1, poff1, mask1, cmps)
        outa[i,:] = out1[0,:]

      resp return= outa

    else:
      resp, times = calib_jungfrau_single_panel(arr, gfac, poff, mask, cmps)

      t17 = time()

      times = (t00, t01, t02, t03, t04) + times + (t17,)
      return resp, np.array(times, dtype=np.float64)


def calib_jungfrau_single_panel(arr, gfac, poff, mask, cmps):
    """ Arrays should have a single panel shape, example for 8-panel detector
    arr:  shape:(8, 512, 1024) size:4194304 dtype:uint16 [2906 2945 2813 2861 3093...]
    poff: shape:(3, 8, 512, 1024) size:12582912 dtype:float32 [2922.283 2938.098 2827.207 2855.296 3080.415...]
    gfac: shape:(3, 8, 512, 1024) size:12582912 dtype:float32 [0.02490437 0.02543429 0.02541406 0.02539831 0.02544083...]
    mask: shape:(8, 512, 1024) size:4194304 dtype:uint8 [1 1 1 1 1...]
    cmps: shape:(16,) size:16 dtype:float64 [  7.   1. 100.   0.   0....]
    """

    t05 = time()

    # Define bool arrays of ranges
    gr0 = arr <  BW1              # 490 us
    gr1 =(arr >= BW1) & (arr<BW2) # 714 us
    gr2 = arr >= BW3              # 400 us

    t06 = time()

    factor = np.select((gr0, gr1, gr2), (gfac[0,:], gfac[1,:], gfac[2,:]), default=1) # 2msec

    t07 = time()

    pedoff = np.select((gr0, gr1, gr2), (poff[0,:], poff[1,:], poff[2,:]), default=0)

    t08 = time()

    # Subtract offset-corrected pedestals
    arrf = np.array(arr & MSK, dtype=np.float32)

    t09 = time()

    arrf -= pedoff

    t10 = time()

    if cmps is not None:
      mode, cormax = int(cmps[1]), cmps[2]
      npixmin = cmps[3] if len(cmps)>3 else 10
      if mode>0:
        #arr1 = store.arr1
        #grhg = np.select((gr0,  gr1), (arr1, arr1), default=0)
        logger.debug(info_ndarr(gr0, 'gain group0'))
        logger.debug(info_ndarr(mask, 'mask'))
        t0_sec_cm = time()

        t11 = time()
        gmask = np.bitwise_and(gr0, mask) if mask is not None else gr0

        t12 = time()

        #sh = (nsegs, 512, 1024)
        hrows = 256 #512/2

        for s in range(arrf.shape[0]):

          t13 = time()

          if mode & 4: # in banks: (512/2,1024/16) = (256,64) pixels # 100 ms
            common_mode_2d_hsplit_nbanks(arrf[s,:hrows,:], mask=gmask[s,:hrows,:], nbanks=16, cormax=cormax, npix_min=npixmin)
            common_mode_2d_hsplit_nbanks(arrf[s,hrows:,:], mask=gmask[s,hrows:,:], nbanks=16, cormax=cormax, npix_min=npixmin)

          t14 = time()

          if mode & 1: # in rows per bank: 1024/16 = 64 pixels # 275 ms
            common_mode_rows_hsplit_nbanks(arrf[s,], mask=gmask[s,], nbanks=16, cormax=cormax, npix_min=npixmin)

          t15 = time()

          if mode & 2: # in cols per bank: 512/2 = 256 pixels  # 290 ms
            common_mode_cols(arrf[s,:hrows,:], mask=gmask[s,:hrows,:], cormax=cormax, npix_min=npixmin)
            common_mode_cols(arrf[s,hrows:,:], mask=gmask[s,hrows:,:], cormax=cormax, npix_min=npixmin)

        logger.debug('TIME: common-mode correction time = %.6f sec' % (time()-t0_sec_cm))

    t16 = time()

    arrf = arrf * factor if mask is None else arrf * factor * mask # gain correction

    return arrf, (t05, t06, t07, t08, t09, t10, t11, t12, t13, t14, t15, t16)