From: Dennis Shea <shea_at_nyahnyahspammersnyahnyah>

Date: Wed Oct 19 2011 - 08:15:09 MDT

Received on Wed Oct 19 08:15:19 2011

Date: Wed Oct 19 2011 - 08:15:09 MDT

You have set opt=1. From the 'pattern_cor' documentation:

http://www.ncl.ucar.edu/Document/Functions/Contributed/pattern_cor.shtml

opt

opt=0 means compute the centered pattern correlation;

opt=1 means compute the uncentered pattern correlation.

=====

Hence the NCL 'pattern_cor' function which, under the hood'

uses the 'pattern_cor2' function, is going to calculate the

*uncentered* pattern correlation because you set opt=1

Your code is calculating the *centered* pattern correlation.

There is a difference.

=====

In your current code, change the last argument to 0.

=====

The best debugging is via print statements.

I have attached a test 'pattern_cor2' function that has print

statements. You can print the equivalent numbers from your computations.

=====

load "./tst_pattern_cor2.ncl"

r_x = tst_pattern_cor2(X,O,1.0,0)

r_y = tst_pattern_cor2(Y,O,1.0,0)

=====

Finally, you are use fortran-77 'do' [or C 'for'] style.

There is nothing wrong with this BUT in any interpreted language

using multiple do loops can be slow. NCL is an array language, like

fortran-90. Array operatons are much quicker but, at times' it can take

some thinking to make sure it is done correctly.

Another side point, if your X,Y,O arrays had missing values

(_FillValue), you would have to add and if statement

if (.not.ismissing(X(i,j)) .or. .not.ismissing(Y(i,j)) .or. \

.not.ismissing(O(i,j)) ) then

:

:

end if

Good luck

On 10/18/11 10:06 PM, Ivanna Mo wrote:

*> Hi,
*

*> I use the built-in function pattern_cor to calculate the pattern
*

*> correlation coefficient of my result. I find that the r is strange and
*

*> different from the result I compute myself. Here is the script that I use:
*

*>
*

*> r_x = pattern_cor(X,O,1.0,1) ; r_x = 0.9998
*

*> r_y = pattern_cor(Y,O,1.0,1) ; r_y= 0.99996
*

*>
*

*>
*

*> Avg_X = avg(X)
*

*> Avg_Y = avg(Y)
*

*> Avg_O = avg(O)
*

*>
*

*> Sum_X = 0.0
*

*> Sum_y = 0.0
*

*> Sum_dev_X = 0.0
*

*> Sum_dev_Y = 0.0
*

*> Sum_dev_O = 0.0
*

*>
*

*> do i=0, 67
*

*> do j=0, 67
*

*>
*

*> Sum_x = Sum_x + (X(i,j) - Avg_X)*(O(i,j) - Avg_O)
*

*> Sum_y = Sum_y + (Y(i,j) - Avg_Y)*(O(i,j) - Avg_O)
*

*> Sum_dev_X = Sum_dev_X + (X(i,j) - Avg_X)*(X(i,j) - Avg_X)
*

*> Sum_dev_Y = Sum_dev_Y +(Y(i,j) - Avg_Y)*(Y(i,j) - Avg_Y)
*

*> Sum_dev_O = Sum_dev_O +(O(i,j) - Avg_O)*(O(i,j) - Avg_O)
*

*>
*

*> end do
*

*> end do
*

*>
*

*> R_x = Sum_x/(sqrt(Sum_dev_X)*sqrt(Sum_dev_O)) ; R_x = 0.70
*

*>
*

*> R_y = Sum_y/(sqrt(Sum_dev_Y)*sqrt(Sum_dev_O)) ; R_y = 0.98
*

*>
*

*>
*

*>
*

*>
*

*>
*

*>
*

*> _______________________________________________
*

*> ncl-talk mailing list
*

*> List instructions, subscriber options, unsubscribe:
*

*> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
*

_______________________________________________

ncl-talk mailing list

List instructions, subscriber options, unsubscribe:

http://mailman.ucar.edu/mailman/listinfo/ncl-talk

- text/plain attachment: tst_pattern_cor2.ncl

*
This archive was generated by hypermail 2.1.8
: Wed Oct 19 2011 - 13:36:10 MDT
*