; Estimates defocus values from power spectra
; SOURCE: ctffind4.spi
;
; PURPOSE: Get defocus using MRC program CTFFIND4
;
; USAGE: spider spi/$DATEXT @ctffind4 $RESULTS ncpus=$NUM_CPUS
;
; REQUIRES:
; load-mic.spi SPIDER procedure
; xmipp_image_convert XMIPP's conversion program (if MRC format)
; readmrc.py a python script
; ctffind4 in your PATH
; v4-ctffind.sh OR a shell script to pass parameters to CTFFIND4
; v4-ctffind-mp.sh a shell script which allocates multiple threads
;
; INPUTS:
; Micrographs (in SPIDER format)
; Start, end micrograph numbers
; Parameters to ctffind, including
; - spherical aberration
; - electron beam voltage (in kV)
; - amplitude contrast ratio
; - magnification
; - pixel size on scanner (in microns)
; (see ctffind script for more)
; OUTPUTS:
; Doc file of defocus and astigmatism information.
; Power spectrum generated by ctffind4 (in SPIDER format).
; Text report generated by ctffind4.
; --------------------- Parameters -----------------------------------
[numProcs] = 1 ; Number of CPUs (g.t. 0, can be overridden from command line)
[keepspi] = 1 ; Save converted SPIDER-format micrograph? (0 = No)
; Exclude edges of the micrograph
[xb] = 100 ; X distance from border
[yb] = 100 ; Y distance from border
; CTFFIND4 parameters
[boxsize] = 512 ; box size
[minres] = 50 ; minimum resolution
[maxres] = 3 ; maximum resolution
[minDf] = 8000 ; minimum defocus (Angstroms)
[maxDf] = 70000 ; maximum defocus (Angstroms)
[dfStep] = 1000 ; defocus step size
[astigExp] = 100 ; expected astigmatism
; For CTF ED and TF ED
[maskRadius] = 100 ; Inner mask radius for power spectrum, Angstroms
[tilesiz] = 500 ; Size of tiles (square)
[xover] = 50 ; X tile overlap %
[yover] = 50 ; Y tile overlap %
[xd] = 500 ; X tile dist. from the edge
[yd] = 500 ; Y tile dist. from the edge
; -------------- Input files -----------------
[params] = '../params' ; Parameters file
[sel_mic] = '../sel_micrograph' ; Micrograph File numbers
[raw] = '../Micrographs/raw{****[mic]}' ; Template for micrographs (w/extension if not SPIDER)
[xic] = 'xmipp_image_convert' ; Path for XMIPP's image conversion program (if MRC format)
; -------------- Output files -----------------
[spider_mic] = '../Micrographs/mic{****[mic]}' ; (Optional) SPIDER-format micrograph (one/micrograph)
[outdir] = 'power' ; Power spectra directory (one)
[mrc_out] = '[outdir]/powmrc****' ; Power spectrum template (one/micrograph)
[ctffind_stdout] = '[outdir]/out****' ; CTFFIND4 stdout (one/micrograph)
[report] = '[outdir]/report****' ; Report file, created by CTFFIND (one/micrograph)
[pow] = '[outdir]/pw_avg_****' ; Power spectrum images (one/micrograph)
[ctf] = '[outdir]/ctf_****' ; CTF noise doc files (one/micrograph)
[rod] = '[outdir]/roo_doc_****' ; Rotational average doc file (one/micrograph)
[defocus_doc] = 'defocus' ; Doc file (one)
; MIC_NUM CTFFIND4_DF ANG_ASTIG ASTIGMATISM CUTOFF SPIDER_DF
; -------------------- END BATCH HEADER -----------------------
; Temporary files
[micout] = 'tmpmic2spi' ; SPIDER temp file from conversion
[temp_winmic] = 'tmpwin{****[mic]}' ; Windowed SPIDER image
[mrc_mic] = 'tmpmrcmic.mrc' ; Temporary MRC file
[ctffind_outs] = 'tmpout_ctffind4' ; Template for CTFFIND4 outputs (.mrc, .txt, & _avrot.txt)
[ctffind_report] = 'tmpreport'
[temp_readmrc_doc] = 'tmp_ctffind_defocus'
[ctfed_out] = 'tmp_ctfed_defocus' ; Defocus values doc file (one)
[rosd_1davg] = '_1'
[temp_ctfed_doc] = 'tmpctfdoc' ; 1D power spectrum from CTF ED
[temp_masked_pws] = '_20' ; Masked power spectrum
MD ; Skip unnecessary output
VB OFF
MD ; Skip unnecessary output
TR OFF
MD
SET MP
1
; Parameters for CTFFIND4. Use symbolic vars to send to shell script.
UD 1,[sp_zflag] ; Get zip flag
[params] ; Params file (input)
UD 2,[sp_fflag] ; Get tif flag
[params] ; Params file (input)
UD 3,[sp_nx] ; Get X dimension
[params] ; Params file (input)
UD 4,[sp_ny] ; Get Y dimension
[params] ; Params file (input)
UD 5,[sp_pixsiz] ; Pixel size (ANGSTROMS)
[params]
UD 6,[sp_kev] ; Electron beam voltage in kV
[params]
UD 7,[sp_sph_abb] ; Spherical aberration
[params]
UD 12,[sp_acr] ; Amplitude contrast ratio
[params]
UD 14,[lambda] ; Wavelength
[params]
UD E
DE
[temp_readmrc_doc]
DE
[ctfed_out]
DE
[defocus_doc]
; Column labels for the doc. output file
SD / MIC_NUM CTFFIND4_DF ANG_ASTIG ASTIGMATISM CUTOFF SPIDER_DF RESOLUTION
[defocus_doc]
SD E
[defocus_doc]
SYS ; Make sure output dir. present
mkdir -pv [outdir]
UD N [numMics] ; Get micrograph filenumbers
[sel_mic]
[decimate] = 1
IF ([keepspi] .GE. 1) THEN
SYS
echo " Will keep micrographs in SPIDER format" ; echo
ENDIF
DO [key] = 1,[numMics] ; Loop over all files -------------------------
UD IC [key], [mic]
[sel_mic] ; [mic] is now the micrograph file number
@load-mic([mic],[sp_zflag],[sp_fflag],[deci],[sp_nx],[sp_ny],[keepspi])
[raw] ; Micrograph (input)
[micout] ; SPIDER file (output)
_2 ; Scratch file (output)
; @convert_p([decimate])
; [params] ; Parameter file
; [raw] ; Micrograph (input)
; [micout] ; SPIDER file (output)
FI [v21],[v22]
[micout]
(12,2)
[v31] = [v21] - (2*[xb]) ; Dimensions without border
[v32] = [v22] - (2*[yb])
WI ; Window the micrograph, removing borders
[micout] ; Micrograph (input)
[temp_winmic] ; Windowed micrograph (output)
([v31],[v32]) ; Dimensions
([xb],[yb]) ; Start coords
IF ([keepspi] .GE. 1) THEN
SYS
mv -vn [micout].$DATEXT [spider_mic].$DATEXT
; (Will not overwrite.)
ENDIF
; CP TO MRC ; Convert to MRC format
; [temp_winmic] ; Micrograph (input)
; [mrc_mic] ; Windowed mrc micrograph (output)
; 32 ; 4 byte real data
;
; ([sp_pixsiz],[sp_pixsiz],[sp_pixsiz]) ; pixel size ; Leave blanked above and below! (call default values)
;
SYS
[xic] -i [temp_winmic].$DATEXT -o [mrc_mic]
SYS
echo " Calling ctffind for micrograph: {****[mic]}"
[sp_sph_abb]
[sp_kev]
[sp_acr]
[sp_pixsiz]
IF ([ncpus] .GE. 2) THEN [numProcs]=[ncpus]
; Can pass optional #CPUs from command-line call
IF ([numProcs] .GE. 2) THEN
IF ([key] .EQ. 1) THEN
SYS
echo " Using {%i0%[numProcs]} processors"
ENDIF
SYS M
./v4-ctffind-mp.sh [mrc_mic] [ctffind_outs].mrc {%f5.3%[sp_pixsiz]} {%f6.0%[sp_kev]}
{%f5.2%[sp_sph_abb]} {%f5.2%[sp_acr]} {%i4%[boxsize]} {%i3%[minres]} {%i2%[maxres]} {%i5%[minDf]}
{%i5%[maxDf]} {%i4%[dfStep]} {%i4%[astigExp]} no {%i0%[numProcs]} > [ctffind_report].txt
.
ELSE
IF ([key] .EQ. 1) THEN
SYS
echo " Using single processor"
ENDIF
SYS M
./v4-ctffind.sh [mrc_mic] [ctffind_outs].mrc {%f5.3%[sp_pixsiz]} {%f6.0%[sp_kev]}
{%f5.2%[sp_sph_abb]} {%f5.2%[sp_acr]} {%i4%[boxsize]} {%i3%[minres]} {%i2%[maxres]} {%i5%[minDf]}
{%i5%[maxDf]} {%i4%[dfStep]} {%i4%[astigExp]} no > [ctffind_report].txt
.
ENDIF
; Variables passed to CTFFIND4:
; $1 Input image file name
; $2 Output diagnostic filename
; $3 Pixel size
; $4 Acceleration voltage
; $5 Spherical aberration
; $6 Amplitude contrast
; $7 Size of power spectrum to compute
; $8 Minimum resolution
; $9 Maximum resolution
; $10 Minimum defocus
; $11 Maximum defocus
; $12 Defocus search step
; $13 Expected (tolerated) astigmatism
; $14 Find additional phase shift?
; $15 (Optional) Number of threads
SYS
mv -f [ctffind_report].txt [report][mic].txt
CP FROM MRC ; Comment out these lines if you don't want spectra
[ctffind_outs].mrc ; MRC power spectra (input)
[mrc_out][mic] ; Power spectrum generated by CTFFIND3. (output)
N ; Do not fold data
SYS
\rm [ctffind_outs].mrc
SYS
./readmrc.py [ctffind_outs].txt {****[mic]} >> [temp_readmrc_doc].$DATEXT
SYS
rm -f [ctffind_outs]_avrot.txt
; The following can be helpful as a diagnostic.
SYS
mv [ctffind_outs].txt [ctffind_stdout][mic].txt
; Estimate CTF defocus parameters & save in defocus summary doc file
CTF ED [defocusCTFED],[cutoff],[ntiles]
[temp_winmic] ; Micrograph image (input)
[tilesiz],[xover],[yover] ; Tile size, x & y tile % overlap:
[xd],[yd] ; X & Y tiling border
[sp_pixsiz],[sp_sph_abb] ; Pixel size, Spherical aberration
[sp_kev] ; Electron voltage [kev]
[sp_acr] ; Ampl. contrast ratio
[ctf][mic] ; CTF noise doc file (output)
[ctfed_out] ; Defocus doc file, appended (output)
[mic] ; Key/image number for doc file
[pow][mic] ; 2D power spectrum (output)
DE
[rod][mic] ; PS rotational average doc file (output)
RO SD ; Rotational average
[pow][mic] ; 2D CTF file (input)
[rosd_1davg] ; 1D rotational average (output, not saved)
[rod][mic] ; PS rotational average doc file (output)
; Read CTFFIND's values
UD [mic], [micMRC],[avgDf],[majAxis],[minAxis],[angCtffind],[astigmatism],[resolution]
[temp_readmrc_doc]
UD E
IF ([key] .EQ. 1) THEN
; Calculate radius
[maskRadius] = 2*[sp_pixsiz]*[tilesiz]/[maskRadius]
[maskRadius] ; diagnostic
ENDIF
; Mask power spectrum
MA
[pow][mic]
[temp_masked_pws]
0,[maskRadius] ; inner, outer radii
D ; _D_isk
C ; set background to _C_ircumference average
; mask center (CR = image center)
DE
[temp_ctfed_doc]
; Write to combined doc
SD [mic], [micMRC],[avgDf],[angCtffind],[astigmatism],[cutoff],[defocusCTFED],[resolution]
[defocus_doc]
SYS
echo
DE
[temp_winmic] ; Windowed micrograph (deleted)
VM
\rm -f [temp_winmic].$DATEXT
VM
\rm [mrc_mic]
ENDDO
; Clean up
UD ICE
[sel_mic] ; [mic] is now the micrograph file number
DE
[micout] ; SPIDER file (deleted)
DE
[ctfed_out]
DE
[temp_readmrc_doc]
SD E
[defocus_doc]
EN
; Modified 2016-10-18
; 2016-10-18 (trs) -- reads resolution estimate from CTFFIND4 output
; 2016-05-23 (trs) -- replaced TF ED 1D profile with CTF ED 1D profile
; 2016-04-25 (trs) -- reordered columns to be compatible with ctffind.spi
; 2016-03-11 (trs) -- option to allocate subset of processors for CTFFIND4 (default uses all CPUs)
; 2016-01-27 (trs) -- added option to save SPIDER-format micrograph
; 2016-01-22 (trs) -- masking power spectrum and then running TF ED
; 2016-01-14 (trs) -- using XMIPP to convert to MRC format
;