def get_error(time, x, sets, err_type='block', tool='gmx analyze'):
    """To estimate error using block averaging method

    .. warning::
            To calculate errors by using ``error = 'acf'`` or ``error = 'block'``,
            GROMACS tool ``g_analyze`` or ``gmx analyze`` should be present in ``$PATH``.

    Parameters
    ----------
    time : 1D list or array
        :attr:`DNA.time`
    x   : 2D list or array
        Shape of (nset, nframe); where *nset* is number of set and *nframe* is
        total number of frames. *nframe* should be equal to length of time
        list/array
    sets : int
        Number of sets (*nset*)
    err_type : str
        Error estimation by autocorrelation method ``err_type='acf'`` or
        block averaging method ``err_type='block'``
    tool : str
        GROMACS tool to calculate error. In older versions it is `g_analyze` while in
        newer versions (above 2016) it is `gmx analyze`.

    Returns
    -------
    error : 1D array
        Of length = number of sets (*nset*)

    """
    for i in range(sets):
        if (len(time) != len(x[i])):
            raise ValueError('\nError: number of frame in time {0} mismatched with {1} of x[{2}]!!\n' .format(
                len(time), len(x[i]), i))

    if not((err_type == 'block') or (err_type == 'acf')):
        print('\nWarning: Method {0} is not implemented. Switching to \'acf\'.\n' .format(
            err_type))
        err_type = 'acf'

    error = []
    char_set = string.ascii_lowercase
    name = ''.join(random.sample(string.ascii_lowercase, 10))

    filename = name + '.xvg'
    eefile = 'ee_' + name + '.xvg'
    acfile = 'acf_' + name + '.xvg'

    fout = open(filename, 'w')
    for i in range(len(time)):
        fout.write('{0}' .format(time[i]))
        for j in range(sets):
            fout.write('    {0}' .format(x[j][i]))
        fout.write("\n")
    fout.close()

    command = '{0} -f {1} -ee {2} -ac {3} -fitfn exp' .format(tool, filename, eefile, acfile)

    try:
        p = sub.Popen(command.split(), stdout=sub.PIPE, stderr=sub.PIPE, universal_newlines=True)
        out, outputerror = p.communicate()
    except Exception as e:
        os.remove(filename)
        raise e

    lines = out.split('\n')

    if (err_type == 'block'):
        for line in lines:
            if(re.match('Set', line)):
                temp = line.split()
                error.append(float(temp[3]))

    if (err_type == 'acf'):
        acf_time = []
        for line in lines:
            if re.match('COR: Correlation time', line) is not None:
                temp = line.split('=')
                acf_time.append(abs(float(temp[1].split()[0])))

        total_time = float(time[-1]) - float(time[0])
        dt = total_time / len(time)
        for i in range(sets):
            if(acf_time[i] >= dt):
                n_indp = total_time / acf_time[i]
                tmp_err = np.std(x[i]) / np.sqrt(n_indp)
            else:
                tmp_err = np.std(x[i]) / np.sqrt(len(time))
            error.append(tmp_err)

    os.remove(filename)
    os.remove(eefile)
    os.remove(acfile)
    if os.path.isfile('fitlog.log'):
        os.remove('fitlog.log')

    return np.array(error)
def get_error ( time , x , sets , err_type = 'block' , tool = 'gmx analyze' ) : for i in range ( sets ) : if ( len ( time ) != len ( x [ i ] ) ) : raise ValueError ( '\nError: number of frame in time {0} mismatched with {1} of x[{2}]!!\n' . format ( len ( time ) , len ( x [ i ] ) , i ) ) if not ( ( err_type == 'block' ) or ( err_type == 'acf' ) ) : print ( '\nWarning: Method {0} is not implemented. Switching to \'acf\'.\n' . format ( err_type ) ) err_type = 'acf' error = [ ] char_set = string . ascii_lowercase name = '' . join ( random . sample ( string . ascii_lowercase , 10 ) ) filename = name + '.xvg' eefile = 'ee_' + name + '.xvg' acfile = 'acf_' + name + '.xvg' fout = open ( filename , 'w' ) for i in range ( len ( time ) ) : fout . write ( '{0}' . format ( time [ i ] ) ) for j in range ( sets ) : fout . write ( '    {0}' . format ( x [ j ] [ i ] ) ) fout . write ( "\n" ) fout . close ( ) command = '{0} -f {1} -ee {2} -ac {3} -fitfn exp' . format ( tool , filename , eefile , acfile ) try : p = sub . Popen ( command . split ( ) , stdout = sub . PIPE , stderr = sub . PIPE , universal_newlines = True ) out , outputerror = p . communicate ( ) except Exception as e : os . remove ( filename ) raise e lines = out . split ( '\n' ) if ( err_type == 'block' ) : for line in lines : if ( re . match ( 'Set' , line ) ) : temp = line . split ( ) error . append ( float ( temp [ 3 ] ) ) if ( err_type == 'acf' ) : acf_time = [ ] for line in lines : if re . match ( 'COR: Correlation time' , line ) is not None : temp = line . split ( '=' ) acf_time . append ( abs ( float ( temp [ 1 ] . split ( ) [ 0 ] ) ) ) total_time = float ( time [ - 1 ] ) - float ( time [ 0 ] ) dt = total_time / len ( time ) for i in range ( sets ) : if ( acf_time [ i ] >= dt ) : n_indp = total_time / acf_time [ i ] tmp_err = np . std ( x [ i ] ) / np . sqrt ( n_indp ) else : tmp_err = np . std ( x [ i ] ) / np . sqrt ( len ( time ) ) error . append ( tmp_err ) os . remove ( filename ) os . remove ( eefile ) os . remove ( acfile ) if os . path . isfile ( 'fitlog.log' ) : os . remove ( 'fitlog.log' ) return np . array ( error )
To estimate error using block averaging method

    .. warning::
            To calculate errors by using ``error = 'acf'`` or ``error = 'block'``,
            GROMACS tool ``g_analyze`` or ``gmx analyze`` should be present in ``$PATH``.

    Parameters
    ----------
    time : 1D list or array
        :attr:`DNA.time`
    x   : 2D list or array
        Shape of (nset, nframe); where *nset* is number of set and *nframe* is
        total number of frames. *nframe* should be equal to length of time
        list/array
    sets : int
        Number of sets (*nset*)
    err_type : str
        Error estimation by autocorrelation method ``err_type='acf'`` or
        block averaging method ``err_type='block'``
    tool : str
        GROMACS tool to calculate error. In older versions it is `g_analyze` while in
        newer versions (above 2016) it is `gmx analyze`.

    Returns
    -------
    error : 1D array
        Of length = number of sets (*nset*)