Changeset 1164


Ignore:
Timestamp:
Dec 9, 2013 11:29:41 AM (8 years ago)
Author:
vondreele
Message:

Good image calibration in all cases thus far tested

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIimage.py

    r1163 r1164  
    281281    delt = amax-amin
    282282    if len(ring) < 10:             #want more than 10 deg
    283         return [],(delt > 90),0
    284     return ring,(delt > 90),sumI/len(ring)
     283        return [],(delt > 90)
     284    return ring,(delt > 90)
    285285   
    286286def GetEllipse2(tth,dxy,dist,cent,tilt,phi):
     
    525525            break
    526526        ellipse = GetEllipse(dsp,data)
    527         Ring,delt,sumI = makeRing(dsp,ellipse,pixLimit,cutoff,scalex,scaley,ma.array(self.ImageZ,mask=tam))
     527        Ring,delt = makeRing(dsp,ellipse,pixLimit,cutoff,scalex,scaley,ma.array(self.ImageZ,mask=tam))
    528528        if Ring:
    529529            if iH >= skip:
     
    582582       
    583583    #setup 360 points on that ring for "good" fit
    584     Ring,delt,sumI = makeRing(1.0,ellipse,pixLimit,cutoff,scalex,scaley,self.ImageZ)
     584    Ring,delt = makeRing(1.0,ellipse,pixLimit,cutoff,scalex,scaley,self.ImageZ)
    585585    if Ring:
    586586        ellipse,err = FitRing(Ring,delt)
    587         Ring,delt,sumI = makeRing(1.0,ellipse,pixLimit,cutoff,scalex,scaley,self.ImageZ)    #do again
     587        Ring,delt = makeRing(1.0,ellipse,pixLimit,cutoff,scalex,scaley,self.ImageZ)    #do again
    588588        ellipse,err = FitRing(Ring,delt)
    589589    else:
     
    617617    dsp = HKL[0][3]
    618618    tth = 2.0*asind(wave/(2.*dsp))
     619    Ring0 = makeRing(dsp,ellipse,3,cutoff,scalex,scaley,self.ImageZ)[0]
    619620    ttth = nptand(tth)
    620621    stth = npsind(tth)
     
    631632    zdisp = radii[1]*ttth*tand(tilt)
    632633    zdism = radii[1]*ttth*tand(-tilt)
    633     print 'zdis',zdisp,zdism
    634634#cone axis position; 2 choices. Which is right?
    635635    centp = [elcent[0]+zdisp*sind(phi),elcent[1]+zdisp*cosd(phi)]
    636636    centm = [elcent[0]+zdism*sind(phi),elcent[1]+zdism*cosd(phi)]
    637     print 'cent',centp,centm
    638637#check get same ellipse parms either way
    639 #now do next ring; estimate either way & check sum of Imax/Imin in 3x3 block around each point
     638#now do next ring; estimate either way & do a FitDetector each way; best fit is correct one
    640639    dsp = HKL[1][3]
    641640    tth = 2.0*asind(wave/(2.*dsp))
    642641    ellipsep = GetEllipse2(tth,0.,dist,centp,tilt,phi)
    643 #    data['ellipses'].append(copy.deepcopy(ellipsep+('g',)))
    644     Ringp,delt,sumIp = makeRing(dsp,ellipsep,3,cutoff,scalex,scaley,self.ImageZ)
     642    ellipsem = GetEllipse2(tth,0.,dist,centm,-tilt,phi)
     643    Ringp = makeRing(dsp,ellipsep,3,cutoff,scalex,scaley,self.ImageZ)[0]
     644    Ringm = makeRing(dsp,ellipsem,3,cutoff,scalex,scaley,self.ImageZ)[0]
     645    parmDict = {'dist':dist,'det-X':centp[0],'det-Y':centp[1],
     646        'tilt':tilt,'phi':phi,'wave':wave,'dep':0.0}
     647    varyList = ['dist','det-X','det-Y','tilt','phi']
    645648    if len(Ringp) > 10:
    646         outEp,errp = FitRing(Ringp,True)
     649        chip,chisqp = FitDetector(np.array(Ring0+Ringp),varyList,parmDict,False)
    647650    else:
    648         errp = 1e6
    649     print '+',sumIp,errp
    650     ellipsem = GetEllipse2(tth,0.,dist,centm,-tilt,phi)
    651 #    data['ellipses'].append(copy.deepcopy(ellipsem+('r',)))
    652     Ringm,delt,sumIm = makeRing(dsp,ellipsem,3,cutoff,scalex,scaley,self.ImageZ)
     651        chip = 1e6
    653652    if len(Ringm) > 10:
    654         outEm,errm = FitRing(Ringm,True)
     653        parmDict['tilt'] *= -1
     654        chim,chisqm = FitDetector(np.array(Ring0+Ringm),varyList,parmDict,False)
    655655    else:
    656         errm = 1e6
    657     print '-',sumIm,errm
    658 #    if errp < errm:
    659     if sumIp > sumIm:
     656        chim = 1e6
     657    if chip < chim:
    660658        data['tilt'] = tilt
    661659        data['center'] = centp
    662         negTilt = 1
    663660    else:
    664661        data['tilt'] = -tilt
    665662        data['center'] = centm
    666         negTilt = -1
    667663    data['rotation'] = phi
    668664    parmDict = {'dist':data['distance'],'det-X':data['center'][0],'det-Y':data['center'][1],
     
    671667    if data['DetDepthRef']:
    672668        varyList.append('dep')
     669    data['rings'] = []
    673670    for i,H in enumerate(HKL):
    674671        dsp = H[3]
     
    681678        data['ellipses'].append(copy.deepcopy(ellipse+('g',)))
    682679        print fmt%('predicted ellipse:',elcent[0],elcent[1],phi,radii[0],radii[1])
    683         Ring,delt,sumI = makeRing(dsp,ellipse,pixLimit,cutoff,scalex,scaley,self.ImageZ)
     680        Ring,delt = makeRing(dsp,ellipse,pixLimit,cutoff,scalex,scaley,self.ImageZ)
    684681        if Ring:
    685682            data['rings'].append(np.array(Ring))
     
    689686                data['distance'] = parmDict['dist']
    690687                data['center'] = [parmDict['det-X'],parmDict['det-Y']]
    691                 data['rotation'] = np.mod(parmDict['phi'],360.0)
     688                data['rotation'] = np.mod(parmDict['phi']+720.,360.0)
    692689                data['tilt'] = parmDict['tilt']
     690                if data['rotation'] >180.:
     691                    data['rotation'] -= 180.
     692                    data['tilt'] *= -1.
    693693                data['DetDepth'] = parmDict['dep']
    694694                elcent,phi,radii = ellipse = GetEllipse(dsp,data)
     
    707707        data['distance'] = parmDict['dist']
    708708        data['center'] = [parmDict['det-X'],parmDict['det-Y']]
    709         data['rotation'] = np.mod(parmDict['phi'],360.0)
     709        data['rotation'] = np.mod(parmDict['phi']+720.,360.0)
    710710        data['tilt'] = parmDict['tilt']
     711        if data['rotation'] >180.:
     712            data['rotation'] -= 180.
     713            data['tilt'] *= -1.
    711714        data['DetDepth'] = parmDict['dep']
    712715    for H in HKL[:N]:
     
    860863    for ring in StrSta['d-zero']:       #get observed x,y,d points for the d-zeros
    861864        ellipse = GetEllipse(ring['Dset'],StaControls)
    862         Ring,delt,sumI = makeRing(ring['Dset'],ellipse,ring['pixLimit'],ring['cutoff'],scalex,scaley,Image)
     865        Ring,delt = makeRing(ring['Dset'],ellipse,ring['pixLimit'],ring['cutoff'],scalex,scaley,Image)
    863866        Ring = np.array(Ring).T
    864867        ring['ImxyObs'] = np.array(Ring[:2])      #need to apply masks to this to eliminate bad points
Note: See TracChangeset for help on using the changeset viewer.