function [PVal,D] = k_count(loc,sims,d,b,M,poly)

% K_COUNT computes the raw K-function for a point population in an areal
% lattice (with boundary file and reference measure given). The P-values
% under the 'Randomness Hypothesis' with respect to the reference measure
% are both saved and plotted.[The P-Value here is the chance of getting a
% count at least as large as the observed count under the Randomness
% Hypothesis.]

%Written by: TONY E. SMITH, 5/11/99 {Modified 2/04/08]

%Programs Called: dist_vec.m, count.m, polyform.m, rand_loc.m
%
% INPUTS:
%     (i)   loc  = population location file [loc(i)=(Xi, Yi),i=1..n]
%     (ii)  sims = number of random patterns generated
%     (iii) d    = maximum distance for k function
%     (iv)  b    = number of ticks (bins) to use on the distance axis 
%     (v)   M    = k-vector of measure values for each of k polygons
%     (vi)  s    = seed to be used in random number generator  
%     (vii) poly = (n:2) matrix describing boundaries of k polygons in the form:
%
%          [1   n1
%           x11 y11 
%           x12 y12
%           .
%           .
%           x11 y11
%           2   n2
%           x21 x22
%           .
%           .
%           .
%           k   nk
%           xk1 yk1
%           .
%           .
%           xk1 yk1]

% DATA OUTPUTS: D = Array of distance values for bins.
%                   [Here each bin, k = 1,.., b, is associated 
%                   with distance value, dist(k) = (k/b)*d.]

%            PVal = Clustering P-Value under null hypothesis at each distance in D


% SCREEN OUTPUT: Plots Clustering PVal against D


n = length(loc(:,1)) ; % number of locations

k = length(M) ; % number of polygons

disp('Construct Count for given points in loc');

Z = polyform(poly) ;

C0 = count(loc,b,d);

disp('Construct Polygon Sampling Distribution');

u = ones(1,k);

m = u*M ;

P = M/m ;

% Form cdf

F = zeros(k,1);

F(1) = P(1) ;

i = 2 ;

while i <= k
   
   F(i) = P(i) + F(i-1) ;
   
   i = i + 1;
   
end


%Compute Sample Counts plus P-Values

Freq = zeros(b,1);

ctr = 1 ;

disp('Begin Simulations');

index = 10;


while ctr <= sims
   
	% Construct Polygon Sample Numbers

	Num = zeros(k,1) ;

	i = 1 ;

	while i <= n
   
        x = rand(1,1) ;
        
   	j = 1 ;
   
   	% Locate Polygon Interval containing x
   
   	while j <= k
      
      	if x <= F(j)
         
         	Num(j) = Num(j) + 1 ;
         
         	break
         
      	end
      
      	j = j + 1 ;
      
   	end
   
   	i = i + 1 ;
   
	end

   
   %Form random point sample
   
  LL = rand_loc(poly,Z,Num) ;        
      
   %Construct Count
   
   C = count(LL,b,d) ;
        
   j = 1 ;
   
   while j <= b
      
      if C(j) >= C0(j)
         
         Freq(j) = Freq(j) + 1;         
               
      end
      
      j = j + 1 ;
      
   end
   
     
   if 100*(ctr/sims)>= index
      
      disp(['percent_done = ',num2str(index)]);
      
      index = index + 10;
      
   end
        
   ctr = ctr + 1;
   
end

PVal = (1 + Freq) ./(sims + 1) ;

D = [1:b].*(d/b);