Hi everyone!

I'm currently working on writing a little script to be able to calculate the average energy and the standard deviation of that energy of these molecular simulations I ran. There's about 600 snapshots for each molecule and I have to find the average and standard deviation in the energy of the snap shots for each molecule.

I am not only a super Perl newb, i'm new to coding as well. My PI helped me out for the first part (using it to find averages) and it worked. The second part was adding to it to be able to calculate the standard deviation.

So I was able to obtain some kind of standard deviation for 3 snapshots but they did not match the standard deviation i used in excel (i used STDEV.P since we have the entire population, N, not N-1).

I can't seem to get it right and I want to see if anyone can help me out here before I go and ask my PI. I want to understand the logic behind this as well, rather than have him just type something out for me.

Thank you! 
#!/usr/bin/env perl

$conf = $ARGV[0];
$n1 = $ARGV[1];
$n2 = $ARGV[2];

$esum0 = 0;
$esum1 = 0;
for($i=$n1;$i<=$n2;$i++){

    $i = sprintf("%0.3d",$i);


    system("analyze $conf.$i E | grep \"ACE\"  -A 1 > log ");


    open(IN,"log");
    @temp = split(" ",<IN>);
    $e0 = $temp[2];
    @temp = split(" ",<IN>);
    $e1 = $temp[1];
    close IN;



   print "$i  $e0 $e1 
";


#---average-calculation---

   $esum0 = $esum0 + $e0;
   $esum1 = $esum1 + $e1;

}

{
   $esum0=$esum0/$n2;
   $esum1=$esum1/$n2;
    $esum0= sprintf("%0.4f",$esum0);
   $esum1= sprintf("%0.4f",$esum1);

}

#averages
my@coolbeans = ('averages',  " $esum0", "$esum1");

print "
@coolbeans\n";




#---standard-deviation---


for($i=$n1;$i<=$n2;$i++){

        $estd0=($esum0-$e0) ** 2;
        $estd0=($estd0/$n2) **0.5;
        $estd0= sprintf ("%0.4f", $estd0);

        $estd1=($esum1-$e1) ** 2;
        $estd1=($estd1/$n2) **0.5;
        $estd1= sprintf ("%0.4f", $estd1);
}

#stdev
my@arroba = ('stdev', "$estd0", "$estd1");
print "
@arroba\n";

[download]

In reply to standard deviation help by vka725

Title:
Use:  <p> text here (a paragraph) </p>
and:  <code> code here </code>
to format your post, it's "PerlMonks-approved HTML":


  • Posts are HTML formatted. Put <p> </p> tags around your paragraphs. Put <code> </code> tags around your code and data!
  • Titles consisting of a single word are discouraged, and in most cases are disallowed outright.
  • Read Where should I post X? if you're not absolutely sure you're posting in the right place.
  • Please read these before you post! —
  • Posts may use any of the Perl Monks Approved HTML tags:
    a, abbr, b, big, blockquote, br, caption, center, col, colgroup, dd, del, details, div, dl, dt, em, font, h1, h2, h3, h4, h5, h6, hr, i, ins, li, ol, p, pre, readmore, small, span, spoiler, strike, strong, sub, summary, sup, table, tbody, td, tfoot, th, thead, tr, tt, u, ul, wbr
  • You may need to use entities for some characters, as follows. (Exception: Within code tags, you can put the characters literally.)
            For:     Use:
    & &amp;
    < &lt;
    > &gt;
    [ &#91;
    ] &#93;
  • Link using PerlMonks shortcuts! What shortcuts can I use for linking?
  • See Writeup Formatting Tips and other pages linked from there for more info.