Changeset 4685 for trunk

Timestamp:

Dec 29, 2020 12:37:31 PM (2 years ago)

Author:

vondreele

Message:

Add gain map calculator

Location:

trunk

Files:

• GSASIIdataGUI.py (modified) (1 diff)
• GSASIIimage.py (modified) (2 diffs)
• GSASIIimgGUI.py (modified) (6 diffs)

trunk/GSASIIdataGUI.py

@@ -3871 +3871 @@
                         self.ErrorDialog('Image sum error','No nonzero image multipliers found')
                         return
-                        
-                        
+                                               
                     newImage = np.array(newImage,dtype=np.int32)                        
                     outname = 'IMG '+result[-1]

trunk/GSASIIimage.py

@@ -1262 +1262 @@
     return useMask
 
+def MakeGainMap(image,Ix,Iy,data,masks,blkSize=128):
+    import scipy.ndimage.filters as sdif
+    Iy *= npcosd(Ix[:-1])**3       #undo parallax 3 or 4?
+    Iy *= (1000./data['distance'])**2    #undo r^2 effect
+    Iy /= np.array(G2pwd.Polarization(data['PolaVal'][0],Ix[:-1],0.)[0])    #undo polarization
+    if data['Oblique'][1]:
+        Iy *= G2pwd.Oblique(data['Oblique'][0],Ix[:-1])     #undo penetration
+    IyInt = scint.interp1d(Ix[:-1],Iy[0],bounds_error=False)
+    GainMap = np.zeros_like(image,dtype=float)
+    #do interpolation on all points - fills in the empty bins; leaves others the same
+    Nx,Ny = data['size']
+    nXBlks = (Nx-1)//blkSize+1
+    nYBlks = (Ny-1)//blkSize+1
+    for iBlk in range(nYBlks):
+        iBeg = iBlk*blkSize
+        iFin = min(iBeg+blkSize,Ny)
+        for jBlk in range(nXBlks):
+            jBeg = jBlk*blkSize
+            jFin = min(jBeg+blkSize,Nx)
+            TA = Make2ThetaAzimuthMap(data,(iBeg,iFin),(jBeg,jFin))           #2-theta & azimuth arrays & create position mask
+            Ipix = IyInt(TA[0])
+            GainMap[iBeg:iFin,jBeg:jFin] = image[iBeg:iFin,jBeg:jFin]/(Ipix*TA[3])
+    GainMap = np.nan_to_num(GainMap)
+    GainMap = sdif.gaussian_filter(GainMap,3.,order=0)
+    return 1./GainMap
+
 def ImageIntegrate(image,data,masks,blkSize=128,returnN=False,useTA=None,useMask=None):
     'Integrate an image; called from OnIntegrateAll and OnIntegrate in G2imgGUI'    #for q, log(q) bins need data['binType']
@@ -1373 +1399 @@
         H0 *= np.array(G2pwd.Polarization(data['PolaVal'][0],H2[:-1],0.)[0])
     H0 /= np.abs(npcosd(H2[:-1]-np.abs(data['det2theta'])))           #parallax correction
+    H0 *= (data['distance']/1000.)**2    #remove r^2 effect
     if 'SASD' in data['type']:
         H0 /= npcosd(H2[:-1])           #one more for small angle scattering data?

trunk/GSASIIimgGUI.py

@@ -122 +122 @@
             GMid = G2gd.GetGPXtreeItemId(G2frame, G2frame.root, gainMap)
             if GMid:
-                Npix,gainfile,imagetag = G2IO.GetCheckImageFile(G2frame,GMid)
-                GMimage = G2IO.GetImageData(G2frame,gainfile,True,ImageTag=imagetag,FormatName=formatName)
+                Npix,gainfile,imagetag = G2IO.GetCheckImageFile(G2frame,GMid)                
+                Gdata = G2frame.GPXtree.GetItemPyData(G2gd.GetGPXtreeItemId(G2frame,GMid,'Image Controls'))
+                gformat = Gdata['formatName']
+                GMimage = G2IO.GetImageData(G2frame,gainfile,True,ImageTag=imagetag,FormatName=gformat)
                 sumImg = sumImg*GMimage/1000
     sumImg -= int(data.get('Flat Bkg',0))
@@ -154 +156 @@
         dlg.Destroy()
         
+    def OnMakeGainMap(event):
+        sumImg = GetImageZ(G2frame,data)
+        masks = copy.deepcopy(G2frame.GPXtree.GetItemPyData(
+            G2gd.GetGPXtreeItemId(G2frame,G2frame.Image,'Masks')))
+        Data = copy.deepcopy(data)
+        #force defaults for GainMap calc
+        Data['IOtth'] = [0.1,60.0]
+        Data['outAzimuths'] = 1
+        Data['LRazimuth'] = [0.,360.]
+        Data['outChannels'] = 2500
+        Data['SampleAbs'] = [0.0,False]
+        Data['binType'] = '2-theta'
+        G2frame.Integrate = G2img.ImageIntegrate(sumImg,Data,masks,blkSize)            
+        Iy,azms,Ix = G2frame.Integrate[:3]
+        GainMap = G2img.MakeGainMap(sumImg,Ix,Iy,Data,masks,blkSize)*1000.
+        Npix,imagefile,imagetag = G2IO.GetCheckImageFile(G2frame,G2frame.Image)
+        pth = os.path.split(os.path.abspath(imagefile))[0]
+        outname = 'GainMap'
+        dlg = wx.FileDialog(G2frame, 'Choose gain map filename', pth,outname, 
+            'G2img files (*.G2img)|*.G2img',wx.FD_SAVE|wx.FD_OVERWRITE_PROMPT)
+        if dlg.ShowModal() == wx.ID_OK:
+            newimagefile = dlg.GetPath()
+            newimagefile = G2IO.FileDlgFixExt(dlg,newimagefile)
+            Data['formatName'] = 'GSAS-II image'
+            Data['range'] = [(800,1200),[800,1200]]
+            G2IO.PutG2Image(newimagefile,[],data,Npix,GainMap)
+            Id = G2frame.GPXtree.AppendItem(parent=G2frame.root,text='IMG '+os.path.split(newimagefile)[1])
+            G2frame.GPXtree.SetItemPyData(Id,[Npix,newimagefile])
+            G2frame.GPXtree.SetItemPyData(G2frame.GPXtree.AppendItem(Id,text='Comments'),[])
+            G2frame.GPXtree.SetItemPyData(G2frame.GPXtree.AppendItem(Id,text='Image Controls'),Data)
+            G2frame.GPXtree.SetItemPyData(G2frame.GPXtree.AppendItem(Id,text='Masks'),masks)
+            G2frame.GPXtree.Expand(Id)
+            G2frame.GPXtree.SelectItem(Id)      #to show the gain map & put it in the list 
+
     G2frame.dataWindow.ClearData()
     G2frame.ImageZ = GetImageZ(G2frame,data)
@@ -194 +230 @@
     tthSizer.Add(wx.StaticText(G2frame.dataWindow,label=' Sample changer position %.2f mm '%data['samplechangerpos']),0,WACV)
     mainSizer.Add(tthSizer,0)
+    if not data['Gain map']:
+        makeGain = wx.Button(G2frame.dataWindow,label='Make gain map from this image? NB: use only an amorphous pattern for this')
+        makeGain.Bind(wx.EVT_BUTTON,OnMakeGainMap)
+        mainSizer.Add(makeGain)
     G2frame.dataWindow.SetDataSize()
 
@@ -473 +513 @@
         '''Integrate image in response to a menu event or from the AutoIntegrate
         dialog. In the latter case, event=None. 
-        TODO: think of not repeating x,y-->2th,azm calc for unchanging calibrations?
         '''
         CleanupMasks(masks)
@@ -1311 +1350 @@
             ResetThresholds()
             wx.CallAfter(G2plt.PlotExposedImage,G2frame,event=event)
-        
+            
         backSizer = wx.FlexGridSizer(0,6,5,5)
         oldFlat = data.get('Flat Bkg',0.)
@@ -1344 +1383 @@
         backSizer.Add((5,5),0)
         backSizer.Add(wx.StaticText(G2frame.dataWindow,-1,' Gain map'),0,WACV)
-        gainMap = wx.ComboBox(parent=G2frame.dataWindow,value=data['Gain map'],choices=Choices,
+        gainMap = wx.ComboBox(G2frame.dataWindow,value=data['Gain map'],choices=Choices,
             style=wx.CB_READONLY|wx.CB_DROPDOWN)
         gainMap.Bind(wx.EVT_COMBOBOX,OnGainMap)