rename compiled to cImageD11

ImageD11/__init__.py

@@ -18,6 +18,6 @@
 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 """
 
-__version__ = "1.7.0"
+__version__ = "1.8.0"
 __author__ = 'Jon Wright',
 __author_email__ = '[email protected]'

ImageD11/fft_index_refac.py

@@ -12,7 +12,7 @@
 import logging, time, sys
 import numpy
 
-from ImageD11 import labelimage, closest
+from ImageD11 import labelimage, cImageD11
 
 def get_options(parser):
     parser.add_option( '-n', '--ngrid',
@@ -100,7 +100,7 @@ def gv_to_grid_new(self, gv):
                   thkl[:,1]*grid.shape[1] + \
                   thkl[:,2]
 
-            closest.put_incr( flatgrid , ind.astype(numpy.intp), vol )
+            cImageD11.put_incr( flatgrid , ind.astype(numpy.intp), vol )
             # print thkl,vol
 
             # vals = numpy.take(grid.flat, ind)

ImageD11/indexing.py

@@ -17,7 +17,7 @@
 
 
 import numpy as np
-from ImageD11 import closest
+from ImageD11 import cImageD11
 from xfab.tools import ubi_to_u, u_to_rod, ubi_to_rod
 
 import math, time, sys, logging
@@ -477,7 +477,7 @@ def find(self):
                       % (i*100./len(i1),len(hits)),
             costheta=np.dot(n2,n1[i])
             if tol > 0: # This is the original algorithm - the closest angle
-                best,diff = closest.closest(costheta,cs)
+                best,diff = cImageD11.closest(costheta,cs)
                 if diff < tol:
                     hits.append( [ diff, i1[i], i2[best] ])
             else:
@@ -562,10 +562,10 @@ def scorethem(self):
                 print self.gv[i]
                 print self.gv[j]
                 raise
-            npk = closest.score(self.unitcell.UBI,gv,tol)
+            npk = cImageD11.score(self.unitcell.UBI,gv,tol)
             if npk > self.minpks:
                 self.unitcell.orient(self.ring_1, self.gv[i,:], self.ring_2, self.gv[j,:],verbose=0,all=True)
-                npks=[closest.score(UBItest,gv,tol) for UBItest in self.unitcell.UBIlist]
+                npks=[cImageD11.score(UBItest,gv,tol) for UBItest in self.unitcell.UBIlist]
                 choice = np.argmax(npks)
                 UBI = self.unitcell.UBIlist[choice]
                 npk = npks[choice]
@@ -614,7 +614,7 @@ def fight_over_peaks(self):
         i = -1
         for ubi in self.ubis:
             i += 1
-            npk = closest.score_and_assign( ubi, self.gv, self.hkl_tol,
+            npk = cImageD11.score_and_assign( ubi, self.gv, self.hkl_tol,
                                             self.drlv2, labels, i)
 
         self.ga = labels
@@ -653,12 +653,11 @@ def saveindexing(self, filename, tol = None):
             i += 1
             # Each ubi has peaks in self.ga
             uinverses.append( np.linalg.inv(ubi) )
-            npk , mdrlv = closest.refine_assigned(
+            npk , mdrlv = cImageD11.refine_assigned(
                 ubi.copy(),
                 self.gv,
                 self.ga,
-                i,
-                -1)
+                i)
             assert npk == self.gas[i]
             f.write("Grain: %d   Npeaks=%d   <drlv>=%f\n"%(
                 i, self.gas[i], np.sqrt(mdrlv) ))
@@ -779,9 +778,9 @@ def score(self,UBI,tol=None):
         """
 #      t0=time.time()
         if tol==None:
-            return closest.score(UBI,self.gvflat,self.hkl_tol)
+            return cImageD11.score(UBI,self.gvflat,self.hkl_tol)
         else:
-            return closest.score(UBI,self.gvflat,tol)
+            return cImageD11.score(UBI,self.gvflat,tol)
         # NONE OF THIS FOLLOWING CODE IS EXECUTED, EVER!
 #        t1=time.time()
 #        h=matrixmultiply(UBI,transpose(self.gv))

ImageD11/refinegrains.py

@@ -23,7 +23,7 @@
 import numpy
 
 from ImageD11 import transform, indexing, parameters
-from ImageD11 import grain, columnfile, closest
+from ImageD11 import grain, columnfile, cImageD11
 
 import simplex, xfab.tools
 
@@ -340,7 +340,7 @@ def compute_gv(self, thisgrain , update_columns = False  ):
             self.gv = gv.T
         if self.OMEGA_FLOAT:
             mat = thisgrain.ubi.copy()
-            junk = closest.score_and_refine(mat , gv.T,
+            junk = cImageD11.score_and_refine(mat , gv.T,
                                             self.tolerance)
             hklf = numpy.dot( mat, gv )
             hkli = numpy.floor( hklf + 0.5 )
@@ -397,7 +397,7 @@ def refine(self, ubi, quiet=True):
         mat=ubi.copy()
         # print "In refine",self.tolerance, self.gv.shape
         # First time fits the mat
-        self.npks, self.avg_drlv2 = closest.score_and_refine(mat, self.gv,
+        self.npks, self.avg_drlv2 = cImageD11.score_and_refine(mat, self.gv,
                                                              self.tolerance)
         # apply symmetry to mat:
         if self.latticesymmetry is not triclinic:
@@ -406,7 +406,7 @@ def refine(self, ubi, quiet=True):
             mat = xfab.tools.u_to_ubi( U, self.latticesymmetry( cp ) )
 
         # Second time updates the score with the new mat
-        self.npks, self.avg_drlv2 = closest.score_and_refine(mat, self.gv,
+        self.npks, self.avg_drlv2 = cImageD11.score_and_refine(mat, self.gv,
                                                              self.tolerance)
         # apply symmetry to mat:
         if self.latticesymmetry is not triclinic:
@@ -423,7 +423,7 @@ def refine(self, ubi, quiet=True):
 #                raise
 
         #print self.tolerance
-#        self.npks, self.avg_drlv2 = closest.score_and_refine(mat, self.gv,
+#        self.npks, self.avg_drlv2 = cImageD11.score_and_refine(mat, self.gv,
 #                                                             self.tolerance)
         #tm = indexing.refine(ubi,self.gv,self.tolerance,quiet=quiet)
         #print ubi, tm,ubi-tm,mat-tm
@@ -666,7 +666,7 @@ def assignlabels(self, quiet=False):
                     gr.translation,
                     gv)
                 fImageD11.assign( gr.ubi, gv, self.tolerance, drlv2, int_tmp, int(g)) 
-#                closest.score_and_assign( gr.ubi,
+#                cImageD11.score_and_assign( gr.ubi,
 #                                          gv.T,
 #                                          self.tolerance,
 #                                          drlv2_2,

ImageD11/rsv_mapper.py

@@ -6,7 +6,7 @@
 
 import numpy, logging
 from ImageD11 import parameters, transform, indexing, \
-    closest, blobcorrector, rsv
+    cImageD11, blobcorrector, rsv
 
 
 class rsv_mapper(object):
@@ -224,17 +224,17 @@ def add_image( self, om, data ):
 
         if msk is not None:
             numpy.multiply(dat, msk, dat)
-            closest.put_incr( self.rsv.SIG,
+            cImageD11.put_incr( self.rsv.SIG,
                               ind,
                               dat )
-            closest.put_incr( self.rsv.MON,
+            cImageD11.put_incr( self.rsv.MON,
                               ind,
                               self.lorfac * msk)
         else:
-            closest.put_incr( self.rsv.SIG,
+            cImageD11.put_incr( self.rsv.SIG,
                               ind,
                               dat )
-            closest.put_incr( self.rsv.MON,
+            cImageD11.put_incr( self.rsv.MON,
                               ind,
                               self.lorfac)
         return

sandbox/art.py

@@ -4,7 +4,7 @@
 Some simplistic ART algorithms
 """
 import numpy as np
-from ImageD11 import closest
+from ImageD11 import cImageD11
 import pylab as pl
 
 
@@ -45,15 +45,15 @@ def update_wtd( recon, proj, angle, msk, dbg=True ):
     r = recon.ravel()
     np.subtract( idx_lo, mn, idx_lo) # Move these to offsets in calc_proj
     np.subtract( idx_hi, mn, idx_hi)
-    closest.put_incr( calc_proj, idx_lo, r*wt_lo )
-    closest.put_incr( calc_proj, idx_hi, r*wt_hi )
+    cImageD11.put_incr( calc_proj, idx_lo, r*wt_lo )
+    cImageD11.put_incr( calc_proj, idx_hi, r*wt_hi )
     error = np.zeros( calc_proj.shape, np.float32 )
     start = (len(calc_proj)- len(proj))/2
     error[ start:start+len(proj) ] = proj
     error = error - calc_proj
     npcalc_proj = np.zeros( (mx-mn), np.float32 )
-    closest.put_incr( npcalc_proj, idx_hi, wt_hi )
-    closest.put_incr( npcalc_proj, idx_lo, wt_lo )
+    cImageD11.put_incr( npcalc_proj, idx_hi, wt_hi )
+    cImageD11.put_incr( npcalc_proj, idx_lo, wt_lo )
     npcalc_proj = np.where(npcalc_proj >0 ,npcalc_proj, 1)
     update =  np.zeros( r.shape, np.float32)
     update =  wt_hi*np.take( error/npcalc_proj, idx_hi ) 
@@ -66,8 +66,8 @@ def update_wtd( recon, proj, angle, msk, dbg=True ):
     # Now take this error off of recon according to the weights
 
     #tot    =  np.zeros( (mx-mn), np.float32 )
-    #closest.put_incr( tot, idx_hi, wt_hi)
-    #closest.put_incr( tot, idx_lo, wt_lo)
+    #cImageD11.put_incr( tot, idx_hi, wt_hi)
+    #cImageD11.put_incr( tot, idx_lo, wt_lo)
 
 
 

sandbox/ftomo.py

@@ -1,5 +1,5 @@
 import math, numpy, time
-from ImageD11 import closest
+from ImageD11 import cImageD11
 from fabio.openimage import openimage
 print "Using class version"
 
@@ -53,7 +53,7 @@ def make_indices(self):
         nim  = numpy.zeros( ( nv* nh), numpy.float32 )
         wons = numpy.ones( (len(inds)), dtype=numpy.float32 )
         # This is now more dense - bincount?
-        closest.put_incr( nim , inds, wons )
+        cImageD11.put_incr( nim , inds, wons )
         nim = nim.astype(numpy.int)
         self.nim_div = nim + (nim==0)
 
@@ -76,8 +76,8 @@ def process_sinogram( self,
         numpy.multiply( faprojc, self.conjer, faprojc)
         fimr = numpy.zeros( self.ftimlen , numpy.float32 )
         fimc = numpy.zeros( self.ftimlen , numpy.float32 )
-        closest.put_incr( fimr, self.inds, faprojr.ravel())
-        closest.put_incr( fimc, self.inds, faprojc.ravel())
+        cImageD11.put_incr( fimr, self.inds, faprojr.ravel())
+        cImageD11.put_incr( fimc, self.inds, faprojc.ravel())
         fim = fimr + fimc*1j
         fim = numpy.divide( fimr + fimc*1j, self.nim_div)
         fim.shape = self.ftimshape 

sandbox/indexer.py

@@ -17,7 +17,7 @@
 
 
 import numpy as np
-from ImageD11 import closest, grain, transform
+from ImageD11 import cImageD11, grain, transform
 import unitcell
 
 import math, time, sys, logging
@@ -192,7 +192,7 @@ def pairs(self, hkl1, hkl2, cos_tol = 0.02, hkl_tol = 0.05):
                T_g[1] = unit( np.cross( this_n1, T_g[2]) )
                U = np.dot( T_g.T, T_c)
                ubi =  np.dot( U, self.unitcell.B) 
-#               npks = closest.score(ubi,gvf,0.1)
+#               npks = cImageD11.score(ubi,gvf,0.1)
                pairs.append( (ind1[i], ind2[k], U, ubi ) )
 #               print npks, ubi
         self.pairs=pairs

sandbox/no_pixel_split_fazit.py

@@ -1,7 +1,7 @@
 
 
 import sys, numpy, time
-from ImageD11 import closest, blobcorrector, ImageD11_file_series
+from ImageD11 import cImageD11, blobcorrector, ImageD11_file_series
 from fabio.openimage import openimage
 from fabio.edfimage import edfimage
 
@@ -103,7 +103,7 @@ def main():
                     on = on * fit2dmask
 
                 # Number of pixels and mask are constant
-                closest.put_incr( outnp,
+                cImageD11.put_incr( outnp,
                                   indices,
                                   on )
 
@@ -135,7 +135,7 @@ def main():
 
             dm = (dataim.ravel()*fit2dmask).astype(numpy.float32)
 
-            closest.put_incr( outsum,
+            cImageD11.put_incr( outsum,
                               indices,
                               dm )
 

sandbox/rad.py

@@ -1,6 +1,6 @@
 import math, numpy, time
 # , smooth
-from ImageD11 import closest
+from ImageD11 import cImageD11
 from fabio.openimage import openimage
 print "Using class version"
 
@@ -60,7 +60,7 @@ def make_indices(self):
         nim  = numpy.zeros( ( nv* nh), numpy.float32 )
         wons = numpy.ones( (len(inds)), dtype=numpy.float32 )
         # This is now more dense - bincount?
-        closest.put_incr( nim , inds, wons )
+        cImageD11.put_incr( nim , inds, wons )
         nim = nim.astype(numpy.int)
         self.nim_div = nim + (nim==0)
 
@@ -83,8 +83,8 @@ def process_sinogram( self,
         numpy.multiply( faprojc, self.conjer, faprojc)
         fimr = numpy.zeros( self.ftimlen , numpy.float32 )
         fimc = numpy.zeros( self.ftimlen , numpy.float32 )
-        closest.put_incr( fimr, self.inds, faprojr)
-        closest.put_incr( fimc, self.inds, faprojc)
+        cImageD11.put_incr( fimr, self.inds, faprojr)
+        cImageD11.put_incr( fimc, self.inds, faprojc)
         fim = fimr + fimc*1j
         fim = numpy.divide( fimr + fimc*1j, self.nim_div)
         fim.shape = self.ftimshape 

src/connectedpixels.c

@@ -7,6 +7,16 @@
 #include "blobs.h"  /* INTEGER */
 
 
+void boundscheck( jpk, n2, ipk, n1){
+    if( (jpk < 0) || jpk>=n2 ) {
+        printf("Bounds check error, jpk, n2\n");
+        exit(0);
+    }
+    if( (ipk < 0) || ipk>=n1 ) {
+        printf("Bounds check error, jpk, n1\n");
+        exit(0);
+    }
+}
 
 
 /*
@@ -277,7 +287,7 @@ int bloboverlaps( INTEGER* b1, INTEGER n1, double* res1,
                    INTEGER* b2, INTEGER n2, double* res2, 
                    int verbose, int ns, int nf){
 
-    int i, j, safelyneed, ipx ;
+    int i, j, safelyneed, ipx, k ;
     INTEGER *link, p1, p2, ipk, jpk, npk, *T;
 
   /* Initialise a disjoint set in link 
@@ -321,21 +331,23 @@ int bloboverlaps( INTEGER* b1, INTEGER n1, double* res1,
         if((i > n2+1) && (j < n2+1)){ /* linking between images */
             jpk = j-1;
             ipk = i-n2-2;
-            assert( (n2 > jpk) && (jpk >= 0 ) && ( n1 > ipk) && ( ipk >= 0 ) );
+            /* Bounds checking */
+            boundscheck( jpk, n2, ipk, n1);
             merge( &res2[NPROPERTY*jpk], &res1[NPROPERTY*ipk]);
             continue;
         }
         if((i > n2+1) && (j > n2+1) ){ /* linking on the same image (1) */
             jpk = j-n2-2;
             ipk = i-n2-2;
-            assert( (n2 > jpk) && (jpk >= 0 ) && ( n1 > ipk) && ( ipk >= 0 ) );
+            boundscheck( jpk, n1, ipk, n1);
+            assert( (n1 > jpk) && (jpk >= 0 ) && ( n1 > ipk) && ( ipk >= 0 ) );
             merge( &res1[NPROPERTY*jpk], &res1[NPROPERTY*ipk]);
             continue;
         }
         if( i < n2+1 && j < n2+1 ){ /* linking on the same image (2) */
             jpk = j-1;
             ipk = i-1;
-            assert( (n2 > jpk) && (jpk >= 0 ) && ( n1 > ipk) && ( ipk >= 0 ) );
+            boundscheck( jpk, n2, ipk, n2);
             merge( &res2[NPROPERTY*jpk], &res2[NPROPERTY*ipk]);
             continue;
         }

test/jpa/uniqify_ubi.py

@@ -2,7 +2,7 @@
 
 import numpy as np
 
-from ImageD11 import closest, indexing, transformer
+from ImageD11 import cImageD11, indexing, transformer
 
 
 # open up a set of ubi files
@@ -57,7 +57,7 @@
         seen = 0
         for j in range(len(uniq_ubis)):
             u = uniq_ubis[j][0]
-            npk=closest.score(u,gv.T,tol)
+            npk=cImageD11.score(u,gv.T,tol)
             if npk == len(pks):
                 # seen you before
                 uniq_ubis[j][1].append((name, i))
@@ -114,7 +114,7 @@ def get_z(name):
         avg = 4
         npix = 3
         intensity = np.sum((t.finalpeaks[avg,:]*t.finalpeaks[avg,:])[ind])
-        print closest.score(entry[0],t.gv.T,tol), intensity
+        print cImageD11.score(entry[0],t.gv.T,tol), intensity
 
 
 # now re