(x31,x36,x37,x56,x59,x89)

; PURPOSE: Computes power sectrum.  Called from: power.spi
;
; SOURCE:  power_p1.spi  

FR
?micrograph image?[mic]
                                          
FR
?pw_avg?[spectrum]
                                          
FR
?roo?[roo]
                                          
FR
?tmparr?[tmparr]


FI x23,x24  ; Get the dimensions of the full image
[mic]
12,2

; x31 = large window size

;x87=((x23-2*x56)*(x24-2*x59))number of pieces in this micrograph

x38=100/(100-x36); Normalization of % of the overlap in x 
x39=100/(100-x37); Normalization of % of the overlap in y

x26=INT(x38*((x23-2*x56)/x31-1)+1)  ; Number of pieces horizontal dim.(X)
x29=INT(x39*((x24-2*x59)/x31-1)+1)  ; Number of pieces vertical dim.(Y)

x94 = 0
; -------------------- Loops over X and Y -------------- 
DO LB2 x12 = 1, x29
   x21 = (x31/x39)*(x12-1) + x59     ; Y-direction
   DO LB3 x13 = 1, x26
      x22 = (x31/x38)*(x13-1) + x56  ; X-direction

      WI         ; Window image
      [mic]      ; Get small pieces of the input image 
      _5
      x31,x31    ; Size of small pieces
      x22,x21    ; Coordinates of the upper left corner

      RA         ; Correct ramp effects 
      _5
      _6

      PW         ; Calc power spectrum of each small piece
      _6
      _7

      SQ         ; Square image             
      _7  
      _9

      x94=x94+1

      IF (x94.ge.2) THEN
         AD      ; Add images
         _3
         _9
         _3
         *
      ELSE
         CP      ; Copy image
         _9
         _3
      ENDIF

   LB3
LB2

; --------- End Loops. Resulting spectrum in _3 ------------

AR                ; Average over power spectra of small pieces 
_3
_8
P1/x94

WU                ; Calculate the square root 
_8
_2   

IF (x89.EQ.0) GOTO LB4  ; Mask out center 
   FI x41,x42
   _2
   (12,2)          ; Get image dimensions

   x43 = x41/2 + 1
   x44 = x42/2 + 1

   GP x49          ; Get pixel intensity
   _2
   (x43+10,5)

   MA              ; Mask image
   _2
   _9
   (0.0,0.05*x41)  ; Radius of the central dot. e.g., .05 x 500 = 25 pixels
   D
   E
   x49
   (x43,x44)

   CP               ; Copy image
   _9
   _2

LB4                 ; End maskout if clause. Spectrum in _2 

CP                  ; Write the output file
_2
[spectrum]

; Create the document file, with column headings & spatial freq. column
RO 
_2
_3                  ; 1D rotational average

LI D
_3
[tmparr]            ; Doc file            (output)
R
1

; Labels for columns in output file
SD /     Amplitude     X-axis    Spat.Freq.
[roo]                 ; Doc file            (output)

x71 = x31/2  ; 1/2 large window size

DO LB5 x70 = 1,x71
   UD IC x70, x72,x73
   [tmparr]                 ; Doc file      (input)

   x74 = x70 / (2.0 * x71)  ; Compute spatial freq.

   SD x70,x72,x73,x74
   [roo]                    ; Doc file         (output)
LB5

RE