Changeset 3029


Ignore:
Timestamp:
Aug 28, 2017 10:44:37 AM (4 years ago)
Author:
vondreele
Message:

implement correctly errat for protein validation; erratv2 not right

Location:
trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIImath.py

    r3003 r3029  
    25992599    xyzmin = 999.*np.ones(3)
    26002600    xyzmax = -999.*np.ones(3)
    2601     #select residue atoms, skipping main chain C & N; make cartesian
     2601    #select residue atoms,S,Se --> O make cartesian
    26022602    for atom in Atoms:
    26032603        if atom[1] in resNames:
    26042604            if atom[4].strip() in ['S','Se']:
     2605                atom[3] = 'Os'
    26052606                atom[4] = 'O'
    2606 #            if atom[3].strip() in ['C','N']:    #skip main chain C & N atoms
    2607 #                continue
    26082607            cartAtoms.append(atom[:cx+3])
    26092608            cartAtoms[-1][cx:cx+3] = np.inner(Amat,cartAtoms[-1][cx:cx+3])
     2609            if atom[3] in ['N','CA','C','O']:
     2610                cartAtoms[-1].append('B')
     2611            else:
     2612                cartAtoms[-1].append('S')
    26102613    XYZ = np.array([atom[cx:cx+3] for atom in cartAtoms])
    26112614    xyzmin = np.array([np.min(XYZ.T[i]) for i in [0,1,2]])
     
    26442647                newChain = True
    26452648        if atom[0] not in res:  #new residue, get residue no.
    2646             if len(res): print res[-1],resname[-1],intact
    26472649            res.append(atom[0])
    26482650            resname.append('%s-%s%s'%(atom[2],atom[0],atom[1]))
     
    26582660                tgts += list(Boxes[jbox[0],jbox[1],jbox[2]])
    26592661        tgts = list(set(tgts))
    2660         tgts = [tgt for tgt in tgts if atom[1:3] != cartAtoms[tgt][1:3]]    #exclude same residue
    2661         if atom[3].strip() in ['C',]:
    2662             tgts = [tgt for tgt in tgts if cartAtoms[tgt][3].strip() not in ['N',]]  #exclude main chain C-N
    2663         if atom[3].strip() in ['N',]:
    2664             tgts = [tgt for tgt in tgts if cartAtoms[tgt][3].strip() not in ['C',]]  #exclude main chain N-C
    26652662        tgts = [tgt for tgt in tgts if np.sum((XYZ[ia]-XYZ[tgt])**2) <= dsmax]
     2663        tgts = [tgt for tgt in tgts if atom[:3] != cartAtoms[tgt][:3]]    #exclude same residue
     2664        ires = int(atom[0])
     2665        if atom[3].strip() == 'C':
     2666            tgts = [tgt for tgt in tgts if not (cartAtoms[tgt][3].strip() == 'N' and int(cartAtoms[tgt][0]) in [ires-1,ires+1])]
     2667        elif atom[3].strip() == 'N':
     2668            tgts = [tgt for tgt in tgts if not (cartAtoms[tgt][3].strip() in ['C','CA'] and int(cartAtoms[tgt][0]) in [ires-1,ires+1])]
     2669        elif atom[3].strip() == 'CA':
     2670            tgts = [tgt for tgt in tgts if not (cartAtoms[tgt][3].strip() == 'N' and int(cartAtoms[tgt][0]) in [ires-1,ires+1])]
    26662671        for tgt in tgts:
    26672672            dsqt = np.sqrt(np.sum((XYZ[ia]-XYZ[tgt])**2))
     
    26802685        Probs = [0.,0.,0.,0.]
    26812686        for i in range(4,nRes-4):
    2682             print i,IntAct[i]
    26832687            mtrx = np.zeros(5)
    26842688            summ = 0.
     
    26972701                    mtrx[3] += IntAct[j]['NN']
    26982702                    mtrx[4] += IntAct[j]['NO']
    2699             print mtrx,np.sum(mtrx)
    2700             print summ
    27012703            mtrx /= summ
    27022704            if old:
     
    27062708                mtrx -= avg
    27072709                prob = np.inner(np.inner(mtrx,b1),mtrx)
    2708             print i, mtrx
    27092710            Probs.append(prob)
    27102711        Probs += 4*[0.,]
  • trunk/GSASIIphsGUI.py

    r3000 r3029  
    32933293        resNames,Probs = G2mth.validProtein(data)
    32943294        G2plt.PlotAAProb(G2frame,resNames,Probs,Title='Error score for %s'%(data['General']['Name']),
    3295             thresh=[17.190823041860433,11.526684477428809])
     3295            thresh=[8.0,6.0])
     3296#            thresh=[17.190823041860433,11.526684477428809])
    32963297
    32973298    def OnSetAll(event):
Note: See TracChangeset for help on using the changeset viewer.