-
Notifications
You must be signed in to change notification settings - Fork 1
/
process_serialEM_montage
169 lines (169 loc) · 6.82 KB
/
process_serialEM_montage
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
#!/bin/bash
echo -e "Batch to process SerialEM mrc stacks: unstack in raw_images folder\n---------------------------------------------------------"
# check if EMAN2 is sourced
if [ -z $EMAN2DIR ]; then
echo "EMAN2 NOT found!" && exit 1
fi
# user input
echo -n "Do you want also to calculate FFTs for each image? (1 - yes, 0 - no) [default: 0] "
read fft
fft=${fft:-0}
re='^[0-1]+$'
if ! [[ $fft =~ $re ]]; then
echo "ERROR! Not a number" >&2; exit 1
fi
# adjust threshold parameters for FFT images, if you want (for display purposes only)
UpperThres=7
LowerThres=0
# check if eman2 is sourced
if [ -z $EMAN2DIR ]; then
echo "EMAN2 NOT found!" && exit 1
fi
[ $fft -eq 1 ] && [ -d FFTs ] && rm -rf FFTs
[ -d raw_images ] && echo "ERROR! Folder raw_images already exists!" && exit 1
[ -f output.log ] && rm -f output.log
[ $fft -eq 1 ] && mkdir FFTs
mkdir raw_images
docfile=`ls -l *.mdoc | head -1 | awk '{print $NF}'`
[ $? -ne 0 ] && echo "ERROR! Need at least one .mdoc file to get pixel size!" && exit 1
grep "SerialEM: Digitized by" $docfile &> /dev/null
[ $? -ne 0 ] && echo "ERROR! SerialEM mdoc file is not recognized!" && exit 1
mag=`grep -m 1 "Magnification" $docfile | sed -e 's/[^0-9]*//g'`
case "$mag" in
"19000")
PixSize="5.5"
;;
"25000")
PixSize="4.27"
;;
"29000")
PixSize="3.65"
;;
"40000")
PixSize="2.6"
;;
"50000")
PixSize="2.12"
;;
"62000")
PixSize="1.71"
;;
*)
PixSize="0"
;;
esac
echo -n "Pixel spacing [detected ${mag}x : ${PixSize}A/px]: "
read PS
: PS=${PS:=$PixSize}
echo -n "Coarse images by a factor of [2]: "
read CO
: CO=${CO:=2}
COPS=`echo "scale=3;$PS*$CO" | bc`
if [ $fft -eq 1 ]; then
echo -n "Mask power spectrums (this is useful for IMAGIC display)? [n] "
read ans
: ans=${ans:="n"}
case $ans in
Y|y)
ans="y"
echo -n "Put a mask on power spectrum up to first zero (A): [50] "
read RES
: RES=${RES:=50}
rad=`echo "scale=2;${COPS}*2/${RES}" | bc`
;;
N|n)
ans="n"
;;
*)
echo "ERROR! Wrong answer!" && exit 1
;;
esac
fi
# end of input handling
tot=`ls *.mrc | wc -l`
echo -e "Batch has started on `date`\n" > output.log
D1=`date +%s`
echo "Pixel size (original): $PS
Coarsening factor: $CO
Number of input SerialEM mrc stacks: $tot
############################################
Micrographs converted:" >> output.log
# Converting mrc SerialEM stacks to mrc files
number=1
for i in `ls *.mrc | sed -e "s/.mrc//g"`
do
if [ ! -e "${i}.mrc.mdoc" ]
then
echo -e "WARNING: No corresponding mdoc file found for ${i}.mrc!\n" | tee -a output.log
fi
echo -en "Writing mrc images...($number out of $tot)\r"
echo "Unstacking ${i}.mrc ..." >> output.log
e2proc2d.py ${i}.mrc raw_images/${i}.mrc --unstacking --threed2twod >> output.log 2>&1
((number++))
done
echo "Writing mrc images...DONE! ($tot out of $tot)" | tee -a output.log
# Renaming files
total=`ls raw_images/*.mrc | wc -l`
zeros=${#total}
count=1
for i in `ls raw_images/*.mrc`
do
countpad=`printf "%0${zeros}d" "$count"`
mv ${i} raw_images/img-${countpad}.mrc
((count++))
echo -en "Renaming mrc images...\r"
done
echo "Renaming mrc images...DONE!" | tee -a output.log
# calculate FFTs if needed
num=`ls raw_images/img-*.mrc | wc -l`
count2=1
if [ $fft -eq 1 ]; then
echo -e "############################################\nNormalizing, coarsing and calculating power spectrum: " >> output.log
for a in `ls raw_images/img-*.mrc | sed -e 's/raw_images\///g;s/.mrc//g'`
do
if [ "$CO" != "1" ]
then
echo -ne "Normalizing, coarsing and calculating power spectrum: ${count2} of ${num}\r"
echo "Processing raw_images/${a}.mrc ..." >> output.log
e2proc2d.py raw_images/${a}.mrc FFTs/pow-stack.img --process=normalize.edgemean --process=threshold.clampminmax.nsigma:nsigma=4 --process=math.realtofft --medianshrink ${CO} >> output.log 2>&1
((count2++))
else
echo -ne "Normalizing and calculating power spectrum: ${count2} of ${num}\r"
echo "Processing raw_images/${a}.mrc ..." >> output.log
e2proc2d.py raw_images/${a}.mrc FFTs/pow-stack.img --process=normalize.edgemean --process=threshold.clampminmax.nsigma:nsigma=4 --process=math.realtofft >> output.log 2>&1
((count2++))
fi
done
echo "############################################" >> output.log
if [ -n "${IMAGIC_ROOT}" ]; then
ls raw_images/img-*.mrc | sed -e 's/raw_images\/img//g;s/.mrc//g' > mics.plt
export IMAGIC_BATCH=1
echo ""
echo "IMAGIC found! Setting micrograph number in header files..." | tee -a output.log
# checked only in IMAGIC version 110308 (08-03-2011)
${IMAGIC_ROOT}/stand/headers.e <<EOF >> output.log 2>&1
WRITE
INDEX/LABEL
LABEL_OF_INDEX
MNUMBER
FILE
mics.plt
FFTs/pow-stack
EOF
fi
if [ ${ans} == "y" ];then
echo -e "\nMasking power spectrums..." | tee -a output.log
e2proc2d.py FFTs/pow-stack.img FFTs/pow-stack-t.img --process=threshold.clampminmax:maxval=${UpperThres}:minval=${LowerThres} --process=normalize.edgemean >> output.log 2>&1
size=`e2iminfo.py -H FFTs/pow-stack.img | grep nx | awk '{print $2}'`
rad2=`echo "scale=2;${rad}*${size}" | bc | cut -d'.' -f1`
e2proc2d.py FFTs/pow-stack-t.img FFTs/pow-stack-t-masked.img --process=mask.soft:inner_radius=${rad2} --process=normalize.edgemean >> output.log 2>&1
fi
echo -e "\nDONE!" | tee -a output.log
if [ ${ans} == "y" ];then
echo "There are three stacks of power spectrums (pixel size is now ${COPS}): pow-stack, pow-stack-t (thresholded) and pow-stack-t-masked (thresholded and masked up to ${RES}A)" | tee -a output.log
else
echo "There are two stacks of power spectrums (pixel size is now ${COPS}): pow-stack and pow-stack-t (thresholded)" | tee -a output.log
fi
fi
D2=`date +%s` ; D=$(( (D2-D1) )) ; printf "Job duration = %dhr %dmin %dsec \n" $(($D/3600)) $(($D/60%60)) $(($D%60)) | tee -a output.log
echo -e "Batch has finished on `date`\n" | tee -a output.log