#!/usr/bin/perl
#
use POSIX;
sub print_help {
    print STDERR <<EOF;
  usage:  multifit filename [options]
      This script splits a two-column x-sorted datafile containing x y
      into blocks of data of a given x-length, and
      returns the average y and the slope of the linear fits in those
      intervals.  It allows practically to extract a smoothed function
      <x>, <y>, <dy/dx>
      This script requires that the program fit is somewhere in the path
  options:
      -l x-width of interval for data splitting (==1 by default)
  output: 4 columns of data:
      average-x   average-y   fit-slope   standard deviation on slope
                                    this version by Nick Manini, 03-02-2021
EOF
return;
}

$tmpfile="/tmp/multi".$$;

$length=1;

while($#ARGV+1){
    $p=shift(@ARGV);
    if($p eq "-h"){# print help
        print_help;
        exit;
    }
    if($p eq "-l"){ # quiet
        $length=shift(@ARGV)+0;
        next;
    }
    if($#filename < 0 ){
	if ($p eq "-"){
	    push(@filename,$p);
	    next;
	}else{
	    if(-e $p){
		push(@filename,$p);
		next;
	    }else{
		die "non existing file: $p\n";
	    }
	}
    }
    if($tin eq "unset"){
	$tin=$p;
        next;
    }
    if($tfin eq "unset"){
	$tfin=$p;
        next;
    }
    if($step eq "unset"){
	$step=$p;
	$step>1 or die "better have at lest 2 steps!!\n";
        next;
    }
    die "error: \"$p\" is an unknown option\n";
}



$file=$filename[0];

$nint=9999999999;
$length2=0.5*$length;

open(OUT, ">", $tmpfile) || die "Cannot open $tmpfile\n";

open(INPUT, $file) || die "Cannot open $file\n";
print "#<x>	<y>	slope dy/dx   standard deviation on slope\n";

while(<INPUT>){
    s/^\s*//;       # remove leading spaces
    if(!/^\#/ && !/^\s*$/){
	@li = split('\s+');
	$x=$li[0]+0;
	if($nint==9999999999){
	    $nint=POSIX::floor($x/$length);
	}
	$xshift=$x-$length*$nint;
	#	print "DEBUG ",$nint," ",$x," ",$xshift,"\n";
	if($xshift>=0.){
	    if($count>1){
		close(OUT);
		$result=`fit -l -w $tmpfile`;
		@boh=split('\s+',$result);
		print $boh[18+23],"\t",$boh[18+37],"\t",$boh[18+11],"\t",$boh[18+13],"\n";
		open(OUT, ">", $tmpfile) || die "Cannot open $tmpfile\n";
	    }
	    $nint++;
	    #    print "DEEEE ",$nint," ",$x," ",$xshift,"\n";
	    $count=0;
	}
	$errorbar=1+2*(($xshift+$length2)/$length2)**8;# weight 1 to 1/9 
	print OUT $x," ",$li[1]," ",$errorbar,"\n"; 
	$count++;
    }
}
close(INPUT);
close(OUT);
# final bunch of data is lost: to correct this here one could copy line from $result=`fit... until print... not really worth!

if (-e $tmpfile){
    `rm $tmpfile`;
}