Révision 4 bin/image2geometry/orient_sepal.py

orient_sepal.py (revision 4)
19 19 
import numpy as np
20 20 
from scipy.ndimage.morphology import binary_closing, binary_opening
21 21 

  		





  22
  
  
    
from titk_tools.io import imread, imsave, SpatialImage


  		





  23
  
  
    
from titk_tools.rw.registration import isotropic_resampling_seg


  		





  24
  
  
    
from titk_tools.rw.registration import apply_transfo_on_seg1


  		





  
  22
  
    
from titk_tools.io import imread, imsave, tfread, SpatialImage


  		





  
  23
  
    
from titk_tools.registration import isotropic_resampling_seg


  		





  
  24
  
    
from titk_tools.registration import apply_transfo_on_seg1


  		





  25
  25
  
    

  		





  26
  26
  
    
from bib_sepals import *


  		





  27
  27
  
    

  		





  28
  
  
    

  		





  29
  
  
    

  		





  30
  28
  
    
print('============== orient_sepal.py =================')


  		





  31
  29
  
    

  		





  32
  30
  
    
filename = argv[1]


  		





  
  31
  
    
outputdir=argv[2]+'/'


  		





  33
  32
  
    

  		





  34
  33
  
    
print(argv)


  		





  35
  34
  
    

  		





  36
  35
  
    
imname = filename.split('/')[-1]


  		





  37
  36
  
    
imnameroot = imname.split('.')[0]


  		





  38
  37
  
    

  		





  39
  
  
    
outputdir = imnameroot+'_oriented/'


  		





  40
  
  
    
os.system("mkdir -p "+outputdir)


  		





  41
  
  
    

  		





  42
  
  
    
isodir = imnameroot+'_iso/'


  		





  43
  
  
    
os.system("mkdir -p "+isodir)


  		





  44
  
  
    

  		





  45
  
  
    
#outname = outputfile.split('/')[-1]


  		





  46
  
  
    
#outroot = outname.split('.')[0]


  		





  47
  
  
    

  		





  48
  
  
    
inrname = outputdir+imnameroot+'.inr.gz'


  		





  49
  38
  
    
outfilename = outputdir+imnameroot+'_oriented.inr.gz'


  		





  50
  
  
    
isoname = isodir+imnameroot+'_iso.inr.gz'


  		





  51
  
  
    
sname = isodir+imnameroot+'_iso_closed.inr.gz'


  		





  52
  39
  
    

  		





  
  40
  
    
# read the image


  		





  
  41
  
    
sepale=imread(filename)


  		





  53
  42
  
    
# resample isotropically


  		





  54
  43
  
    
print("----------------------")


  		





  55
  44
  
    
print("isotropic resampling of : filename=", filename)


  		





  56
  45
  
    
print("----------------------")


  		





  57
  
  
    
sepale=imread(filename)


  		





  58
  46
  
    
voxelsize = sepale.voxelsize


  		





  59
  47
  
    
print("voxelsize=",voxelsize)


  		





  60
  48
  
    
s = sepale.shape


  		





  ......




  66
  54
  
    
print("new voxelsize=",isovs)


  		





  67
  55
  
    

  		





  68
  56
  
    
sepale=SpatialImage(sepale>0, dtype='uint8', voxelsize=voxelsize, origin=[0, 0, 0])


  		





  69
  
  
    
imsave(inrname, sepale)


  		





  70
  
  
    
isotropic_resampling_seg(inrname, isoname, isovs)


  		





  
  57
  
    
sepale=isotropic_resampling_seg(sepale,isovs)


  		





  71
  58
  
    

  		





  72
  
  
    
sepale=imread(isoname)


  		





  73
  59
  
    
s = sepale.shape


  		





  74
  60
  
    
voxelsize = sepale.voxelsize


  		





  75
  
  
    
print('new voxelsize read=', voxelsize)


  		





  76
  61
  
    
print('new shape', s)


  		





  77
  62
  
    

  		





  
  63
  
    

  		





  78
  64
  
    
# alongate the target image so that a sepal originally along the diagonal fits into the target image


  		





  79
  65
  
    
print("----------------------")


  		





  80
  66
  
    
print("elongating field : filename=", filename)


  		





  81
  67
  
    
print("----------------------")


  		





  82
  68
  
    
s_addy=int((np.sqrt(s[0]*s[1]*2)-np.sqrt(s[0]*s[1]))/2)


  		





  83
  
  
    
#s_addy=int(np.sqrt(s[0]*s[1])(1.-np.sqrt(2.))/2. * 1.05)


  		





  84
  69
  
    
s_addz=int(s[2]*0.20)


  		





  85
  
  
    
sepale_long=add_side_slices(sepale, 0, 0, s_addy, s_addy, s_addz, s_addz)


  		





  86
  
  
    
sepale=SpatialImage(sepale_long, dtype='uint8', voxelsize=voxelsize, origin=[0, 0, 0])


  		





  87
  
  
    
imsave(isoname, sepale)


  		





  88
  
  
    
print('--------------------')


  		





  
  70
  
    
sepale=add_side_slices(sepale, 0, 0, s_addy, s_addy, s_addz, s_addz)


  		





  89
  71
  
    

  		





  
  72
  
    
'''


  		





  90
  73
  
    
# smooth it with opening


  		





  91
  
  
    
sepale=imread(isoname)


  		





  
  74
  
    
#sepale=imread(isoname)


  		





  92
  75
  
    
s = sepale.shape


  		





  93
  76
  
    
voxelsize = sepale.voxelsize


  		





  94
  77
  
    
struct_el=struct_element_sphere(radius)


  		





  ......




  96
  79
  
    
struct_el=struct_element_sphere(radius/2)


  		





  97
  80
  
    
sepale = binary_closing(sepale, structure=struct_el, iterations=int(radius/2))


  		





  98
  81
  
    
sepale=SpatialImage(sepale*255, dtype='uint8', voxelsize=voxelsize, origin=[0, 0, 0])


  		





  99
  
  
    
imsave(sname, sepale)


  		





  
  82
  
    
#imsave(sname, sepale)


  		





  
  83
  
    
'''


  		





  100
  84
  
    

  		





  101
  
  
    

  		





  102
  85
  
    
# put the center of mass in the center of the image


  		





  103
  86
  
    
sepale_centered=center_on_cm(sepale)


  		





  104
  
  
    
voxelsize = sepale_centered.resolution


  		





  105
  
  
    
sepale_centered_filename = outputdir+imnameroot+'_centered.inr.gz'


  		





  106
  
  
    
imsave(sepale_centered_filename, sepale_centered)


  		





  107
  87
  
    

  		





  108
  
  
    

  		





  109
  88
  
    
# compute principal directions


  		





  110
  
  
    
pv_flo=principal_directions(sepale_centered)


  		





  
  89
  
    
pv_flo=principal_directions(sepale_centered*255)


  		





  111
  90
  
    
print(pv_flo)


  		





  112
  91
  
    

  		





  113
  92
  
    
# sepals originally are oriented with top in more-or less "up" direction on the images


  		





  ......




  137
  116
  
    
trsf_file=outputdir+imnameroot+'_rigid.txt'


  		





  138
  117
  
    
write_superposing_transfo(pv_flo, pv_ref, sepale_centered, trsf_file)


  		





  139
  118
  
    

  		





  
  119
  
    
trsf=tfread(trsf_file)


  		





  
  120
  
    
os.system("rm "+trsf_file)


  		





  
  121
  
    

  		





  140
  122
  
    
# apply transformation


  		





  141
  
  
    
sepale_filename = outputdir+imnameroot+'_oriented.inr.gz'


  		





  142
  
  
    
apply_transfo_on_seg1(trsf_file, sepale_centered_filename, sepale_filename)


  		





  
  123
  
    
oriented=apply_transfo_on_seg1(trsf, sepale_centered)


  		





  143
  124
  
    

  		





  144
  125
  
    
# smooth the transformed image


  		





  145
  
  
    
sepale=imread(sepale_filename)


  		





  146
  
  
    
sepale = binary_closing(sepale, structure=struct_el, iterations=radius)


  		





  
  126
  
    
#sepale = binary_closing(oriented, structure=struct_el, iterations=radius)


  		





  147
  127
  
    

  		





  148
  128
  
    
# and save the oriented sepal


  		





  149
  
  
    
sepale=SpatialImage(sepale*255, dtype='uint8', voxelsize=voxelsize, origin=[0, 0, 0])


  		





  
  129
  
    
sepale_filename = outputdir+imnameroot+'_oriented.inr.gz'


  		





  
  130
  
    
sepale=SpatialImage(oriented*255, dtype='uint8', voxelsize=voxelsize, origin=[0, 0, 0])


  		





  150
  131
  
    
imsave(sepale_filename, sepale)


  		





  151
  132
  
    

  		












Formats disponibles :
  Unified diff