sens.py

def estimate_sensitivity(msName):
    # find integration time
    time_os = au.timeOnSource(msName)
    print("-------------------------------------------------")
    int_time = time_os[0]['minutes_on_source'] # always in field 0 (from split)
    print "Integration time = ", int_time, " min"

    vm = au.ValueMapping(msName)
    nant = vm.numAntennas # find number of antenna
    npol = vm.nPolarizations # number of polarizations
    spwlist = vm.spwInfo # estimate bandwidth
    bandw = 0.0
    for key in spwlist:
        bandw += spwlist[key]['bandwidth']

    bandw = bandw*1.0e-9 # in GHz
    print "Bandwidth = ", bandw, " GHz"
    print "Number of antenna = ", nant
    print "Number of polarization = ", npol
    
    tb.open(msName + "/ASDM_CALATMOSPHERE")
    tSys = tb.getcol("tSys")
    bnd = tb.getcol("receiverBand")
    tb.close()

    # find band
    bandused = 0
    for b in listband:
        if np.any(bnd == listband[b]):
            bandused = b
            print "Band = ", bandused

    # estimate Tsys
    tsys = tSys.flatten() # 1D array
    tsys0 = tsys[tsys > 0.0] # remove 0's part
    # remove outliers, value > 5 sigma
    meantsys0 = np.mean(tsys0)
    std5 = 5.0*np.std(tsys0)
    tsys5 = tsys0[abs(tsys0 - meantsys0) < std5] 

    meanTsys = np.mean(tsys5)
    print "mean Tsys = ", meanTsys



    #  Assume apropri sensitivity for band 3,4,6,7,8,9,10
    #  Sensitivity units are in mJy and are normalized to:
    #    16 12-m antennas,
    #    8 GHz bandwidth
    #    Dual freq bandwidth
    #    for tsys_nominal given above
    #    Integration time of one minute
    apriori_sensitivity = sensitivities[bandused]

    # sigma = 2*k*Tsys/(aperture_eff * np.sqrt(nant*(nant-1) * npol * bandwidth * int_time))
    # scaling
    scalingTsys = meanTsys/tsys_nominal[bandused]
    scalingInttime_bandwidth_npol = 1.0/np.sqrt(int_time * bandw/8.0 * npol/2.0)
    scalingNant = np.sqrt(240.0/(nant*(nant-1))) # 240 = 16*(16-1)
    scalingfactor = scalingTsys * scalingInttime_bandwidth_npol * scalingNant

    sensitivity = apriori_sensitivity * scalingfactor
    

    return sensitivity # mJy