#!/usr/bin/env perl
# Calculate the self-energy and the spectral function via the self-energy trick
# General approach, for QS symmetry type
# Rok Zitko, rok.zitko@ijs.si, Oct 2008, June 2012

use warnings;
use strict;
use Math::Complex;

my @omega = readcol("c-imG-d.dat", 1);

my @imG = readcol("c-imG-d.dat", 2);
my @reG = readcol("c-reG-d.dat", 2);
my @imF = readcol("c-imF-d.dat", 2);
my @reF = readcol("c-reF-d.dat", 2);

my @imdelta = readcol("Delta-d.dat", 2); # not prefixed!
my @redelta = readcol("Delta-d-re.dat", 2);

my $selffn = "c-self-d.dat";
open(SELF, ">$selffn") or die "Can't open $selffn for writing, stopped";

my $imsigma = "imsigma-d.dat";
open(IM, ">$imsigma") or die "Can't open $imsigma for writing, stopped";

my $resigma = "resigma-d.dat";
open(RE, ">$resigma") or die "Can't open $resigma for writing, stopped";

my $len = @omega;
for (my $i = 0; $i < $len; $i++) {
  my $o = $omega[$i];
  my $G = $reG[$i] + i * $imG[$i];
  my $F = $reF[$i] + i * $imF[$i];
  my $sigma = $F/$G;
  my $delta = i * $imdelta[$i] + $redelta[$i];
 
  my $gf = 1.0/($o + $delta - $sigma);
  my $aw = -1.0/pi * Im($gf);

  print SELF "$o $aw\n";

  my $re = Re($sigma);
  my $im = Im($sigma);
  
  print IM "$o $im\n";
  print RE "$o $re\n";
}

# Read a column from a file
sub readcol {
    my ($filename, $column) = @_;
    my @l;
    open(F, "<$filename") or die "Can't open $filename for reading, stopped";
    while (<F>) {
	chomp;
	my $x = (split)[$column-1];
	push(@l, $x);
    }
    close(F);
    @l
};