user3032425
user3032425

Reputation: 53

Calculating Distance Between Atoms by atomic Coordinates

I got this Script below from this Website.

use strict;
use warnings;

open(OUTPUT,">/home/guest/Desktop/test.txt") or die ("couldn't write the file");
my @line;
for (<>) {
    push @line, $_;               
    next if (@line)<2;
    print OUTPUT distance(@line), "\n";
}
sub distance {
    my @a = split ' ',shift;   
    my @b = split ' ',;   
    my $x = ($a[6]-$b[6]);      
    my $y = ($a[7]-$b[7]);
    my $z = ($a[8]-$b[8]);
    my $dist = sprintf "%.1f",sqrt($x**2+$y**2+$z**2);   
   return "$a[1] to $b[1] |$a[11]-$b[11]| Distance-\t$dist";
}

I m using this script to Identify the distance between atoms. I can achieve only with one atom.

INPUT

HETATM    1  N   MSE A   1      16.272  16.789  53.256  1.00 31.95           N  
HETATM    2  CA  MSE A   1      15.408  15.939  52.384  1.00 30.64           C  
HETATM    3  C   MSE A   1      13.919  16.266  52.544  1.00 27.80           C  
HETATM    4  O   MSE A   1      13.393  16.321  53.659  1.00 26.16           O

OUTPUT

1 to 2 |N-C| Distance-  1.5
1 to 3 |N-C| Distance-  2.5
1 to 4 |N-O| Distance-  2.9

DESIRED OUTPUT

1 to 2 |N-C| Distance- 1.5
1 to 3 |N-C| Distance- 2.5
1 to 4 |N-O| Distance- 2.9
2 to 1 |N-C| Distance- 1.5
2 to 3 |N-C| Distance- 2.5
2 to 4 |N-O| Distance- 2.9
3 to 1 |N-C| Distance- 1.5
3 to 2 |N-C| Distance- 2.5
3 to 4 |N-O| Distance- 2.9
4 to 1 |N-C| Distance- 1.5
4 to 2 |N-C| Distance- 2.5
4 to 3 |N-O| Distance- 2.9

Please help me out

Upvotes: 1

Views: 522

Answers (2)

Borodin
Borodin

Reputation: 126762

The line

my @b = split ' ',;

splits the $_ variable into array @b, when what you want is the top of the @_ array.

It also looks like you want the distance between every pair of atoms calculated, which you need a double-loop to achieve.

If you write a subroutine called distance, then it is best if it returns a distance, rather than some arbitrary string that the calling code happens to want.

This program seems to do what you want, and uses printf for the output.

use strict;
use warnings;

open my $out, '>', '/home/guest/Desktop/test.txt' or die "Couldn't open output file: $!";
my @data = map [ split ], grep /\S/, <>;

for my $i (0 .. $#data) {
  for my $j (0 .. $#data) {

    next if $i == $j;

    my @coords_i = @{$data[$i]}[6,7,8];
    my @coords_j = @{$data[$j]}[6,7,8];

    printf $out "%s to %s |%s-%s| Distance-%.3f\n",
        $data[$i][1], $data[$j][1],
        $data[$i][11], $data[$j][11],
        distance(\@coords_i, \@coords_j);
  }
}

sub distance {
  my ($aa, $bb) = @_;
  my ($x, $y, $z) = map { $aa->[$_] - $bb->[$_] } 0 .. $#$aa;
  return sqrt($x * $x + $y * $y + $z * $z);
}

output

1 to 2 |N-C| Distance-1.493
1 to 3 |N-C| Distance-2.513
1 to 4 |N-O| Distance-2.944
2 to 1 |C-N| Distance-1.493
2 to 3 |C-C| Distance-1.533
2 to 4 |C-O| Distance-2.415
3 to 1 |C-N| Distance-2.513
3 to 2 |C-C| Distance-1.533
3 to 4 |C-O| Distance-1.234
4 to 1 |O-N| Distance-2.944
4 to 2 |O-C| Distance-2.415
4 to 3 |O-C| Distance-1.234

Upvotes: 3

Leeor
Leeor

Reputation: 19736

You have 4 lines in the file, so your loop would have 4 iterations (and would skip the first one since the list is still of size 1).

If you want to match each 2 elements, you need a nested loop that would run n^2 iterations given n lines of input. Just switch your main loop into this

for (<>) {
    push @line, $_;               
}
for $a (@line) {
    for $b (@line) {
          next if ($a eq $b);
          print OUTPUT distance(@line), "\n"; 
    }
}

Upvotes: 1

Related Questions