title: Weighted correlation function in python
body: <p>I tried to implement a weighted correlation function based on the article in the link, formula number 2:</p>
  
  <p><a href="http://staff.ustc.edu.cn/~lshao/papers/paper07.pdf" rel="nofollow noreferrer">http://staff.ustc.edu.cn/~lshao/papers/paper07.pdf</a></p>
  
  <p>let suppose to have 3 vectors s, r and w each of n elements. 
  The vector w is obtained from the following formulas:</p>
  
  <pre><code>w = |r|/(1+D)
  D = |s - k*r|
  k = (r_Transpose * s)/(r_Transpose*r)
  </code></pre>
  
  <p>I would like to implement the formula for the weighted correlation function described on the article. It is correct my implementation? 
  I start from a matrix of dimension [224,640] which means that i have 640 vectors of 224 elements. I would like to calculated the weighted correlation coefficient between those 640 vector respect one other vector - r. Every one of those 640 vector is the vector s.</p>
  
  <pre><code>ref = reference  
  ref_mean = np.mean(ref)  # calcolo il valore medio dello spettro di riferimento
  sens = 190  
  
  
  frame_correlation = np.zeros((1,640))
  img_correlation = np.zeros((nf,npixels))
  
  for i in range(nf):
      frame_test = dati_new[:,:,i]   Selection of one matrix from a cell of matrices
      for j in range(npixels):
  
          spettro_test = frame_test[:,j]     # is my vector s
          spettro_test = np.reshape(spettro_test,(224,1))
          spettro_test_mean = np.mean(spettro_test) 
  
          k = np.dot(np.transpose(ref),spettro_test)/np.dot(np.transpose(ref),ref)
          k = k[0][0]
          D = np.abs(spettro_test - k*ref)
          W = np.abs(ref)/(1+D)
  
          # NUMERATOR OF FORMULA IN THE ARTICLE
          numeratore = np.sum(W*(spettro_test - spettro_test_mean)*(ref - ref_mean))
  
          # DENOMINATOR
          den1_ex = np.sqrt(np.sum(W*np.power(spettro_test - spettro_test_mean,2)))
          den2_ex = np.sqrt(np.sum(W*np.power(ref  - ref_mean,2)))
          denominatore = den1_ex * den2_ex
          rho = numeratore/denominatore
  
          if rho &lt; 0:
              rho  = 0
          if rho  &gt; 1: # for safety reason
              rho = 1
  
          if rho &gt;=0.99:
              rho = (sens*rho)/100
          frame_correlation[:,j]= rho
      img_correlation[i,:] = frame_correlation
  </code></pre>