Changeset 298


Ignore:
Timestamp:
Jun 9, 2011 4:05:07 PM (10 years ago)
Author:
vondreele
Message:

GSASIIimage.py - new print routines for ellipse fitting

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIimage.py

    r292 r298  
    126126    from scipy.optimize import leastsq
    127127       
     128    def CalibPrint(ValSig):
     129        print 'Image Parameters:'
     130        ptlbls = 'names :'
     131        ptstr =  'values:'
     132        sigstr = 'esds  :'
     133        for name,fmt,value,sig in ValSig:
     134            ptlbls += "%s" % (name.rjust(12))
     135            if name == 'rotate':
     136                ptstr += fmt % (value-90.)      #kluge to get rotation from vertical - see GSASIIimgGUI
     137            else:
     138                ptstr += fmt % (value)
     139            if sig:
     140                sigstr += fmt % (sig)
     141            else:
     142                sigstr += 12*' '
     143        print ptlbls
     144        print ptstr
     145        print sigstr       
     146       
    128147    def ellipseCalcD(B,xyd,wave):
    129148        x = xyd[0]
     
    144163        YR = X*npsind(phi)+Y*npcosd(phi)
    145164        return (XR/R0)**2+(YR/R1)**2-1
     165       
    146166    def ellipseCalcW(C,xyd):
    147167        dist,x0,y0,phi,tilt,wave = C
    148168        B = dist,x0,y0,phi,tilt
    149169        return ellipseCalcD(B,xyd,wave)
    150     result = leastsq(ellipseCalcD,p0,args=(rings.T,wave))
     170       
     171    names = ['distance','det-X0','det-Y0','rotate','tilt','wavelength']
     172    fmt = ['%12.2f','%12.2f','%12.2f','%12.2f','%12.2f','%12.5f']   
     173    result = leastsq(ellipseCalcD,p0,args=(rings.T,wave),full_output=True)
     174    result[0][3] = np.mod(result[0][3],360.0)               #remove random full circles
     175    vals = list(result[0])
     176    vals.append(wave)
     177    chi = np.sqrt(np.sum(ellipseCalcD(result[0],rings.T,wave)**2))
     178    sig = list(chi*np.sqrt(np.diag(result[1])))
     179    sig.append(0.)
     180    ValSig = zip(names,fmt,vals,sig)
     181    CalibPrint(ValSig)
    151182    if len(rings) > 1:
     183        print 'Trial refinement of wavelength - not used for calibration'
    152184        p0 = result[0]
    153185        p0 = np.append(p0,wave)
    154         resultW = leastsq(ellipseCalcW,p0,args=(rings.T,))
     186        resultW = leastsq(ellipseCalcW,p0,args=(rings.T,),full_output=True)
     187        resultW[0][3] = np.mod(result[0][3],360.0)          #remove random full circles
     188        sig = np.sqrt(np.sum(ellipseCalcW(resultW[0],rings.T)**2))
     189        ValSig = zip(names,fmt,resultW[0],sig*np.sqrt(np.diag(resultW[1])))
     190        CalibPrint(ValSig)
    155191        return result[0],resultW[0][-1]
    156192    else:
     
    344380    import copy
    345381    import ImageCalibrants as calFile
    346     print 'image calibrate'
     382    print 'Image calibration:'
    347383    time0 = time.time()
    348384    ring = data['ring']
     
    360396    outE = FitRing(ring,True)
    361397    if outE:
    362         print 'start ellipse:',outE
     398#        print 'start ellipse:',outE
    363399        ellipse = outE
    364400    else:
     
    374410        print '1st ring not sufficiently complete to proceed'
    375411        return False
    376     print 'inner ring:',ellipse
     412#    print 'inner ring:',ellipse
    377413    data['center'] = copy.copy(ellipse[0])           #not right!! (but useful for now)
    378414    data['ellipses'].append(ellipse[:]+('r',))
     
    441477            distR = 1.-dist/data['distance']
    442478            if abs(distR) > 0.1:
    443                 print dsp,dist,data['distance'],distR,len(Ring),delt
     479#                print dsp,dist,data['distance'],distR,len(Ring),delt
    444480                break
    445481            if distR > 0.001:
     
    471507                    print 'Bad fit for ring # %i. Try reducing Pixel search range'%(i)
    472508            cent = data['center']
    473             print ('for ring # %2i @ d-space %.4f: dist %.3f rotate %6.2f tilt %6.2f Xcent %.3f Ycent %.3f Npts %d'
    474                 %(i,dsp,dist,phi-90.,Tilt,cent[0],cent[1],numZ))
     509#            print ('for ring # %2i @ d-space %.4f: dist %.3f rotate %6.2f tilt %6.2f Xcent %.3f Ycent %.3f Npts %d'
     510#                %(i,dsp,dist,phi-90.,Tilt,cent[0],cent[1],numZ))
    475511            data['ellipses'].append(copy.deepcopy(ellipse+('r',)))
    476512        else:
     
    507543    p0 = [dist,cent[0],cent[1],phi,tilt]
    508544    result,newWave = FitDetector(rings,p0,wave)
    509     print 'Suggested new wavelength = ',('%.5f')%(newWave),' (not reliable if distance > 2m)'
     545#    print 'Suggested new wavelength = ',('%.5f')%(newWave),' (not reliable if distance > 2m)'
    510546    data['distance'] = result[0]
    511547    data['center'] = result[1:3]
Note: See TracChangeset for help on using the changeset viewer.