Poscar2xcrysden.pl

go back to Main Page, Computational Resources, Scripts, Scripts for VASP

#!/usr/bin/perl

#----------------------------------------------------------------
# PROGRAM FOR  GENERATE XCRYSDEN INPUTS    
# OF VASP POSCAR/CONTCAR FILES   
#
# UTILIZATION:
# perl Poscar2xcrysden.pl file1 file2 .... fileN
#
# By Javier Carrasco (carrasco@qf.ub.es)
# Modified By Monica Garcia  (mgarcia@iciq.es)
#----------------------------------------------------------------

#READ LABELS
open(IN,$ARGV[0]);
while (<IN>) {
if ($.==6){
 @at=split(' ');
          }; 
};
close(IN);
for ($p=0;$p<@at;$p++){
  $k=$p+1;
  print "Label $k, $at[$p] atom/s of:\n";
  $label[$p]=<STDIN>;chomp($label[$p]);
};

#READ poscar/contcar FILES 
for ($r=0; $r<@ARGV;$r++) {
#for ($r=0;$r<$#ARGV+1;$r++) { #$numArgs = $#ARGV + 1; 
#print "$ARGV[$r] : $ENV{$ARGV[$r]}\n";
open (IN,$ARGV[$r]);
        $i=0;
        @at=();
        while (<IN>) {
 		if ($. ==1){ $name=$_;chomp($name);};
 		if ($. ==2) { $escala=$_;chomp($escala);};
 		if ($. ==3) {
     			@fields=split(' ');
     			$ax=$fields[0];
     			$ay=$fields[1];
     			$az=$fields[2];
 		};
  		if ($.==4){
     			@fields=split(' ');
     			$bx=$fields[0];
      			$by=$fields[1];
     			$bz=$fields[2];
 		};
  		if ($.==5){
     			@fields=split(' ');
     			$cx=$fields[0];
      			$cy=$fields[1];
     			$cz=$fields[2];
		 };
 		if ($.==6){
     			@at=split(' ');
                       #print "@at";
                       $NAT=0;
      			$suma=0;
      			$sum=0;
      			for ($j=0;$j<@at;$j++){
        			$at[$j]=$at[$j];
         			$sum=$sum+1;
         			$NAT=$NAT+$at[$j];
      			};
                       #print "$NAT";
 		};
 		if (($.>8)&&($.<$NAT+9)){
     			@fields=split (' ');
     			$x[$i]=$fields[0];$y[$i]=$fields[1];$z[$i]=$fields[2];
     			$optx[$i]=$fields[3];$opty[$i]=$fields[4];$optz[$i]=$fields[5];
     			$i++;
 		};
};
close (IN);

# DATA PROCESSING

# Convert fractional coordinates to positive values 
for ($p=0;$p<$NAT;$p++){
  if ($x[$p]<0.0) {
     $x[$p]=$x[$p]+1;
  };
  if ($y[$p]<0.0) {
     $y[$p]=$y[$p]+1;
  };
  if ($z[$p]<0.0) {
     $z[$p]=$z[$p]+1;
  };
};

# 1.Transformation to cartesian coordinates  
for ($i=0;$i<$NAT;$i++){
  $xc[$i]=($escala)*($x[$i]*$ax+$y[$i]*$bx+$z[$i]*$cx);
  #print "$xc[$i]";
  $yc[$i]=($escala)*($x[$i]*$ay+$y[$i]*$by+$z[$i]*$cy);
  #print "$yc[$i]";
  $zc[$i]=($escala)*($x[$i]*$az+$y[$i]*$bz+$z[$i]*$cz);
  #print "$zc[$i]";
};

# 2.Write results in suitable format 
@output=<$ARGV[$r].xcrysden>; chomp(@output);
open (OUT1, ">@output"); 
print OUT1 "set xsfStructure {\n";
print OUT1 " DIM-GROUP\n";
print OUT1 " 3 1\n";
print OUT1 " PRIMVEC\n";
$eax=$escala*$ax;$eay=$escala*$ay;$eaz=$escala*$az;
printf OUT1 "%3.10f %3.10f %3.10f\n",$eax,$eay,$eaz;
$ebx=$escala*$bx;$eby=$escala*$by;$ebz=$escala*$bz;
printf OUT1 "%3.10f %3.10f %3.10f\n",$ebx,$eby,$ebz;
$ecx=$escala*$cx;$ecy=$escala*$cy;$ecz=$escala*$cz;
printf OUT1 "%3.10f %3.10f %3.10f\n",$ecx,$ecy,$ecz;
print OUT1 " PRIMCOORD\n";
printf OUT1 "%5.0f 1\n",$NAT;

#Initialization  $index
 $index=0;
 for ($p=0;$p<@at;$p++){
  for ($i=$index;$i<$index+$at[$p];$i++){
  printf OUT1 "$label[$p] %2.8f %2.8f %2.8f\n",$xc[$i],$yc[$i],$zc[$i];
  };
  $index=$index+$at[$p];
 };

print OUT1 "}\n";
print OUT1 "WriteFile   \$system(SCRDIR)/xc_tmp.\$system(PID) \$xsfStructure w\n";
print OUT1 "catch {exec \$system(BINDIR)/xsf2xsf \$system(SCRDIR)/xc_tmp.\$system(PID) \$system(SCRDIR)/xc_struc.\$system(PID)}\n";
print OUT1 "xc_openstr xcr \$system(SCRDIR)/xc_struc.\$system(PID) .mesa PL\n";
print OUT1 "UpdateWMTitle  tempest.xcrysden\n";
print OUT1 "set working_XSF_file \$system(SCRDIR)/xc_struc.\$system(PID)\n";
print OUT1 "Get_sInfoArray\n";
print OUT1 "DisplayDefaultMode\n";
print OUT1 " xcAppendState render\n";
print OUT1 " xcUpdateState\n";
print OUT1 "SetWatchCursor\n";
print OUT1 " set xcCursor(dont_update) 1\n";
print OUT1 "set nxdir 1\n";
print OUT1 "set nydir 1\n";
print OUT1 "set nzdir 1\n";
print OUT1 "array set radio {space {SpaceFill based on covalent radii} .mesa,bg #000000 cellmode prim frames rods unitrep cell hexamode parapipedal ball      {Balls based on covalent radii}}\n";
print OUT1 "CellMode 1\n";
print OUT1 "DispFramesAs\n";
print OUT1 "xc_newvalue .mesa 2\n";
print OUT1 "xc_newvalue .mesa 0\n";
print OUT1 "set xcCursor(dont_update) 0\n";
print OUT1 "ResetCursor\n";
close (OUT1);
print "\n";
print "Result in @output\n";
print "\n";
print "                                      J.C.R\n";
};