OpenStructure
spoke_pattern.py

This script can be used to create a "spoke pattern" image in the style of the ones featured in the following paper:

Downing K.H., Glaeser R.M., Restoration of weak phase-contrast images recorded with a high degree of defocus: The "twin image" problem associated with CTF correction, Ultramicroscopy 108 (2008), 921-928

All tweakable parameters are defined in the first lines of the script, which can be edited before the script is executed. The script generates the image, opens a viewer to show it and saves it in the current folder with the name specified by the user.

Usage:

dng spoke_pattern.py <output image>



1 import math
2 import sys
3 import ost.img.alg
4 
5 # Tweakable parameters
6 wedge_start = 20 # Width of wedge at the narrower point
7 wedge_end = 500 # Width of wedge at the broader point
8 number_of_bands = 11 # Number of alternating bands (must be odd)
9 size_of_image_x = 500 # Size of the image (X)
10 size_of_image_y = 500 # Size of the image (Y)
11 pixel_sampling = 1.0 * Units.A # Pixel width
12 threshold = 2.0 * Units.A # Threshold for low pass filtering
13 amplitude = 255 # amplitude of the obscillations
14 
15 # Image is created
16 image=img.CreateImage(img.Size(size_of_image_x,size_of_image_y))
17 image.CenterSpatialOrigin()
18 image.SetSpatialSampling(pixel_sampling)
19 
20 # For a cleaner code, 4 shorthands are defined
21 image_extent=image.GetExtent()
22 start_y=image_extent.GetStart()[1]
23 end_y=image_extent.GetEnd()[1]
24 start_x=image_extent.GetStart()[0]
25 end_x=image_extent.GetEnd()[0]
26 
27 # Wedge is written in the image
28 for y in range (start_y,end_y+1):
29  wedge_width=wedge_start+(y-start_y)*(wedge_end-wedge_start)/(end_y-start_y)
30  for x in range (start_x,end_x+1):
31  half_wedge=wedge_width/2.0
32  outer_bands=(number_of_bands-1)/2
33  factor=(0.5+outer_bands)*math.pi/half_wedge
34  if float(x)>-half_wedge and float(x)<half_wedge:
35  value=255*math.cos(float(x)*factor)
36  image.SetReal(img.Point(x,y),value)
37 
38 # Image is low-pass filtered
39 filter=ost.img.alg.GaussianLowPassFilter(threshold)
40 image.ApplyIP(filter)
41 
42 # Viewer is launched to show the result
43 v=gui.CreateDataViewer(image,"Modulated Image")