;
;
; PURPOSE: Calculates average power spectra for a set of micrographs.
; First converts a scanned file to SPIDER format.
; Places power spectra in: power/roo****
; Estimates defocus from 2D power spectra and places in: defocus
; Uses SPIDER operation 'TF ED'
; Inputs: Window size, percentage of the overlap
; distance of the window from the micrograph border
; PS mask radius
;
; SOURCE: spider/docs/techs/recon/newprogs/powdefocus.spi
;
; Edit following parameters and filenames as needed.
;
; [deci] = Decimation factor for 'DC S' operation
; 0 = Get value from params file (key=16)
; 1 = Full sized image
; 2 = 1/2 size
; 4 = 1/4 size
[deci] = 0 ; Decimation factor (0 = Get value from params file)
[winsiz] = 500 ; Window size of small pieces (Sx=Sy)
[xover] = 50 ; X overlap %
[yover] = 50 ; Y overlap %
[xd] = 500 ; X Dist. from the edge
[yd] = 500 ; Y Dist. from the edge
[radius] = 225 ; Power spectrum mask radius (A) (0 = don't mask)
[dodefocus] = 1 ; Find defocus values
[keepspi] = 1 ; Keep the on-disk temp spider file (0 = discard)
[maskr] = 0 ; Mask output power spectrum this far from border (0 = nomask)
; ----------- Input files --------------
[params] = '../params' ; Parameter doc file
[sel_mic] = '../sel_micrograph' ; Micrograph selection doc file
[micgr] = '../Micrographs/raw{****[mic]}' ; Micrograph images
; ----------- Output files --------------
[outdir] = 'power' ; Power spectra directory
[pow] = '[outdir]/pw_avg{****[mic]}' ; Power spectra images
[roo] = '[outdir]/roo{****[mic]}' ; Rotational averagedoc file
[ctf] = '[outdir]/ctf{****[mic]}' ; Spectrum, envelope, noise doc file
[out] = 'defocus' ; Defocus values doc file
; ----------- Temp file --------------
[spi] = '_1' ; Temp Spider file
; -------------- END BATCH HEADER --------------------------
MD
TR OFF ; Loop info turned off
MD
VB OFF ; File info turned off
MD
() OFF ; No need for () in DO loops
MD
SET MP ; Use all available OMP processors
0
VM ; Make sure output dir. present
mkdir -p [outdir]
; ----- Get zip & format flags (can params vary??)
UD 1,[zflag] ; Get zip flag
[params] ; Params file (input)
UD 2,[fflag] ; Get tif flag
[params] ; Params file (input)
UD 3,[nx] ; HiScan X parameters
[params] ; Params file (input)
UD 4,[ny] ; HiScan Y dimension
[params] ; Params file (input)
UD 5,[sp_pixsiz] ; Get pixel size
[params] ; Params file (input)
IF ([deci].LT.1) THEN
UD 16,[deci] ; Get decimation factor
[params] ; Params file (input)
IF ([deci].LT.1) [deci] = 1.0 ; Should not be zero
ENDIF
;;[sp_pixsiz]=[sp_pixsiz]*[deci] ; Adjust pixel size for decimation??
;;(I think we should NOT be multiplying by the pixel size.
;; Every other usage of key #16 in PARAMS ignores it. --trs)
IF ([radius] .GT. 0) THEN
[circf] = 2*[sp_pixsiz]/[radius] ; Mask radius (now in 1/px)
ENDIF
IF ([dodefocus] .GT. 0) THEN
; Want to determine defocus parameters for each group
UD 7,[sp_sph_abb] ; Spherical aberration (mm)
[params] ; Params file (input)
UD 12,[sp_acr] ; Amplitude contrast ratio
[params] ; Params file (input)
UD 14,[sp_lambda] ; Lambda
[params] ; Params file (input)
; Put defocus parameters title in doc file
DE
[out] ; Defocus file (removed)
SD / Micrograph Defocus Astig.Ang Astig.Mag Cutoff.Freq
[out] ; Title for defocus file (output)
ENDIF
[tmparr] = 'tmparray' ; Name for temporary incore doc file
[winsizd2]=[winsiz]/2
SD IC NEW ; Create an in-core doc file array
[tmparr]
4,[winsizd2] ; Creates a 4 column doc file
UD N [nummics] ; Get number of micrographs
[sel_mic]
DO [key]=1,[nummics] ; Loop over all micrographs -------------------
UD [key], [mic] ; Get micrograph number
[sel_mic] ; Micrograph select doc file (input)
VM
echo ' 'Loading Micrograph: {*****[mic]}
; Convert micrographs and load into incore file
@loadmic([mic],[zflag],[fflag],[deci],[nx],[ny],[keepspi])
[micgr] ; Micrograph template (input)
[spi] ; Spider file (output)
_4 ; Hiscan & Nikon scratch (output)
VM
echo " Creating spectra: [pow] and [roo]"
DE ; If output doc file exists, delete
[roo] ; (removed)
FI H [x],[y] ; Get dimensions of converted image
[spi] ; Micrograph image (input)
NSAM,NROW
; [winsiz] = large window size
;[v87]=(([x]-2*[xd])*([y]-2*[yd]))number of pieces in this micrograph
[nox]=100/(100-[xover]); Normalization of % of the overlap in x
[noy]=100/(100-[yover]); Normalization of % of the overlap in y
[npx]=INT([nox]*(([x]-2*[xd])/[winsiz]-1)+1) ; No. of pieces horizontal dim.(X)
[npy]=INT([noy]*(([y]-2*[yd])/[winsiz]-1)+1) ; No. of pieces vertical dim.(Y)
[v94] = 0
; -------------------- Loops over X and Y --------------
DO [v12] = 1, [npy]
[uly] = ([winsiz]/[noy])*([v12]-1) + [yd] ; Y-direction
DO [v13] = 1, [npx]
[ulx] = ([winsiz]/[nox])*([v13]-1) + [xd] ; X-direction
WI ; Window image
[spi] ; Get small pieces of input image
_5
[winsiz],[winsiz] ; Size of small pieces
[ulx],[uly] ; Coord. of upper left corner
RA ; Correct ramp effects
_5 ; (input)
_6 ; (output)
PW ; Calc power spectrum of each small piece
_6 ; (input)
_7 ; PS (output)
SQ ; Square image
_7 ; (input)
_9 ; Squared (output)
[v94] = [v94]+1
IF ([v94].ge.2) THEN
AD ; Add images
_3 ; (input)
_9 ; (input)
_3 ; Sum (output)
*
ELSE
CP ; Copy image
_9 ; (input)
_3 ; (output)
ENDIF
ENDDO
ENDDO
; --------- End Loops. Resulting spectrum in: _3 ------------
AR ; Average over power spectra of small pieces
_3 ; Inline file (input)
_8 ; Inline file (output)
P1/[v94]
WU ; Calculate the square root
_8 ; Inline file (input)
_2 ; Inline file (output)
IF ([circf].GT.0.0) THEN ; Mask out center
FI H [nsam],[nrow]
_2 ; Inline file (input)
NSAM,NROW ; Get image dimensions
[cx] = [nsam]/2 + 1 ; x center
[cy] = [nrow]/2 + 1 ; Y center
GP [pv] ; Get pixel intensity
_2 ; Inline file (input)
([cx]+10,5) ; Location for background
[circ] =[circf]*[nsam] ; Mask radius
[circo]=[maskr] ; Outer Mask radius
if ([maskr].GT. 0) THEN
[circo]=[cx]-[maskr] ; Outer Mask radius
ENDIF
MA ; Mask image
_2 ; Inline file (input)
_9 ; Inline file (output)
[circo],[circ] ; Radii of central dot. e.g. (2*5A/px/200A)/x 500 = 25 pixels
D ; Disk
E ; External background
[pv] ; Background
[cx],[cy] ; Center coordinates
VM
echo " Masked spectra: {*****[circo]} and {***[circ]}"
CP ; Copy image
_9 ; Inline file (input)
[pow] ; Disk power spectrum (output)
ELSE
CP ; Copy to disk
_2 ; Inline spectrum file (input)
[pow] ; Disk power spectrum (output)
ENDIF
RO ; Rotational average
_2 ; File (input)
_3 ; 1D rotational average (output)
LI D ; List in doc file
_3 ; 1D rotational average (input)
[tmparr] ; Doc file (output)
R ; Row wanted
(1) ; Which row
; Labels for columns in output doc file
SD / Amplitude X-axis Spat.Freq.
[roo] ; Doc file (output)
[winsizd2t] = [winsiz]/2 ; 1/2 large window size
DO [v70] = 1,[winsizd2]
UD IC [v70], [v72],[v73]
[tmparr] ; Doc file (input)
; Compute spatial freq.
[sfreq] = [v70] / (2.0 * [winsizd2t])
SD [v70],[v72],[v73],[sfreq]
[roo] ; PS Doc file (output)
ENDDO
; ----------- Get defocus value -------------------
IF ([dodefocus] .GT. 0) THEN
; Save defocus parameters for each group in doc file
DE ; Delete group doc file
[ctf] ; (removed)
; Transfer Function -- Estimation of CTF parameters
TF ED [v12],[v13],[v14],[v15],[v16]
[pow] ; 2D spectrum (input)
([sp_pixsiz], [sp_sph_abb]) ; Pixel size, Spherical aberration
([sp_lambda]) ; Lambda
([sp_acr]) ; Ampl. contrast ratio
[ctf] ; Doc file (output)
; Save parameters for this group in summary doc file
SD [key], [mic],[v14],[v12],[v13],[v16]
[out] ; Doc file (output)
VM
echo ' 'Micrograph: {*****[mic]}' 'Defocus: {%f8.2%[v14]};echo ' '
ENDIF
ENDDO
SD E
[out] ; Doc file (output)
VM
echo ' ' ; cat [out].$DATEXT ; echo ' '
EN
;