Changeset 282


Ignore:
Timestamp:
May 5, 2011 2:19:40 PM (11 years ago)
Author:
vondreele
Message:

fix "off by one" error in polymask.for
fix azimuth error in integration & plotting

Location:
trunk
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIimage.py

    r271 r282  
    202202    cphi = cosd(phi)
    203203    sphi = sind(phi)
    204     ring = []
    205204    if azm:
    206         aR = azm[0]-90,azm[1]-90,1
     205        aR = azm[0]-90,azm[1]-90,azm[1]-azm[0]
    207206        if azm[1]-azm[0] > 180:
    208             aR[2] = 2
    209     else:
    210         aR = 0,362,2
    211     for a in range(aR[0],aR[1],aR[2]):
    212         x = radii[0]*cosd(a-phi)
    213         y = radii[1]*sind(a-phi)
    214         X = (cphi*x-sphi*y+cent[0])
    215         Y = (sphi*x+cphi*y+cent[1])
    216         ring.append([X,Y])
    217     return ring
    218    
     207            aR[2] /= 2
     208    else:
     209        aR = 0,362,181
     210       
     211    a = np.linspace(aR[0],aR[1],aR[2])
     212    x = radii[0]*npcosd(a-phi)
     213    y = radii[1]*npsind(a-phi)
     214    X = (cphi*x-sphi*y+cent[0])
     215    Y = (sphi*x+cphi*y+cent[1])
     216    return zip(X,Y)
     217               
    219218def calcDist(radii,tth):
    220219    stth = sind(tth)
     
    511510    scaley = pixelSize[1]/1000.
    512511   
    513     tay,tax = np.mgrid[iLim[0]+0.5:iLim[1]+.5,jLim[0]+.5:jLim[1]+.5]         #bin centers not corners
     512#    tay,tax = np.mgrid[iLim[0]+0.5:iLim[1]+.5,jLim[0]+.5:jLim[1]+.5]         #bin centers not corners
     513    tay,tax = np.mgrid[iLim[0]:iLim[1],jLim[0]:jLim[1]]         #bin corners
    514514    tax = np.asfarray(tax*scalex,dtype=np.float32)
    515515    tay = np.asfarray(tay*scaley,dtype=np.float32)
     
    519519    spots = masks['Points']
    520520    polygons = masks['Polygons']
    521     tam = ma.make_mask_none((iLim[1]-iLim[0],jLim[1]-jLim[0]))
     521    tam = ma.make_mask_none((nI,nJ))
    522522    for polygon in polygons:
    523523        if polygon:
    524524            tamp = ma.make_mask_none((nI*nJ))
    525             tam = ma.mask_or(tam.flatten(),ma.make_mask(pm.polymask(nI*nJ,
    526                 tax.flatten(),tay.flatten(),len(polygon),polygon,tamp)))
     525            tamp = ma.make_mask(pm.polymask(nI*nJ,tax.flatten(),
     526                tay.flatten(),len(polygon),polygon,tamp))
     527            tam = ma.mask_or(tam.flatten(),tamp)
    527528    if tam.shape: tam = np.reshape(tam,(nI,nJ))
    528529    for X,Y,diam in spots:
    529         tam = ma.mask_or(tam,ma.getmask(ma.masked_less((tax-X)**2+(tay-Y)**2,(diam/2.)**2)))
    530     return np.array(GetTthAzm(tax,tay,data)),tam           #2-theta & azimuth arrays & position mask
     530        tamp = ma.getmask(ma.masked_less((tax-X)**2+(tay-Y)**2,(diam/2.)**2))
     531        tam = ma.mask_or(tam,tamp)
     532    TA = np.array(GetTthAzm(tax,tay,data))
     533    TA[1] = np.where(TA[1]<0,TA[1]+360,TA[1])
     534    return np.array(TA),tam           #2-theta & azimuth arrays & position mask
    531535
    532536def Fill2ThetaAzimuthMap(masks,TA,tam,image):
     
    540544        tam = ma.mask_or(tam.flatten(),ma.getmask(ma.masked_inside(tay.flatten(),max(0.01,tth-thick/2.),tth+thick/2.)))
    541545    for tth,azm,thick in arcs:
    542         tam = ma.mask_or(tam.flatten(),ma.getmask(ma.masked_inside(tay.flatten(),max(0.01,tth-thick/2.),tth+thick/2.))* \
    543             ma.getmask(ma.masked_inside(tax.flatten(),azm[0],azm[1])))
     546        tamt = ma.getmask(ma.masked_inside(tay.flatten(),max(0.01,tth-thick/2.),tth+thick/2.))
     547        tama = ma.getmask(ma.masked_inside(tax.flatten(),azm[0],azm[1]))
     548        tam = ma.mask_or(tam.flatten(),tamt*tama)
    544549    taz = ma.masked_outside(image.flatten(),int(Zlim[0]),Zlim[1])
    545550    tam = ma.mask_or(tam.flatten(),ma.getmask(taz))
     
    554559    print 'Begin image integration'
    555560    LUtth = data['IOtth']
    556     if data['fullIntegrate']:
    557         LRazm = [0,360]
    558     else:
    559         LRazm = data['LRazimuth']
     561    LRazm = data['LRazimuth']
    560562    numAzms = data['outAzimuths']
    561563    numChans = data['outChannels']
     
    583585                print 'Process map block:',iBlk,jBlk,' limits:',iBeg,iFin,jBeg,jFin
    584586                TA,tam = Make2ThetaAzimuthMap(data,masks,(iBeg,iFin),(jBeg,jFin))           #2-theta & azimuth arrays & create position mask
     587               
    585588                Nup += 1
    586589                dlg.Update(Nup)
  • trunk/GSASIIplot.py

    r281 r282  
    16481648            x,y = np.hsplit(ring,2)
    16491649            tth,azm = G2img.GetTthAzm(x,y,Data)
     1650            azm = np.where(azm < 0.,azm+360,azm)
    16501651            Plot.plot(tth,azm,'b,')
    16511652    if not newPlot:
  • trunk/fsource/polymask.for

    r100 r282  
    1616      REAL*4       X(0:N-1),Y(0:N-1)
    1717      REAL*8       POLY(0:M-1,0:1)
    18       LOGICAL*1    MASK(N)
     18      LOGICAL*1    MASK(0:N-1)
    1919
    2020      INTEGER*4    I,K
Note: See TracChangeset for help on using the changeset viewer.