#!/usr/bin/perl -w
use strict;
my($ra1,$dec1,$ra2,$dec2,$r1,$r2,$f,$f1,$f2,$nu,$d,$d1,$d2,$i,$j,@list,@t,@rline,@fline,@radius,@flux,@index);
my(@flux_list,@ra_list,@dec_list);

#  inputs

#printf(" * $ARGV[0] $ARGV[1] $ARGV[2] $ARGV[3] $ARGV[4] $ARGV[5] $ARGV[6]\n"); exit 0;

#$ra1="10h44m33.00s";
#$dec1="10d00'00.11\"";

#$r1="00d00'00.00\"";
#$r2="00d24'00.00\"";

#$f1=1;
#$f2=10;

#$nu=10;

#printf("Target: $ra1 $dec1 , R1: $r1 , R2: $r2 , F1: $f1, F2: $f2\n");

#$ra1=ratorad($ra1);
#$dec1=dectorad($dec1);

#$r1=dectorad($r1);
#$r2=dectorad($r2);

# inputs from php

$ra1=$ARGV[0];
$dec1=$ARGV[1];

$r1=$ARGV[2];
$r2=$ARGV[3];

$f1=$ARGV[4];
$f2=$ARGV[5];

$nu=$ARGV[6];

printf("<br>Nvss catlog search within the radius of %2.1f deg. around the target.<br>\n",$r2*57.29577951308232);

# read nvss list
chomp(@list=     `cat nvss_radian.dat`);
chomp(@flux_list=`cat nvss_radian.dat | awk ' { print \$1 } '`);
chomp(@ra_list=  `cat nvss_radian.dat | awk ' { print \$2 } '`);
chomp(@dec_list= `cat nvss_radian.dat | awk ' { print \$3 } '`);

# search in the list 
$i=0;$j=0;
foreach (@list) { 

	$f=$flux_list[$j];
	$ra2=$ra_list[$j];
	$dec2=$dec_list[$j];

	if( ($f >= $f1) && ($f <= $f2) ) { 
		# calculate the angular distance between target and caltalog source 
		$d=sin($dec1)*sin($dec2)+cos($dec1)*cos($dec2)*cos($ra1-$ra2);
		$d=acos($d); 

		if( ($d >= $r1) && ($d <= $r2) ) { 
		$radius[$i]=$d;
		$flux[$i]=$f*(($nu/1.4)**(-0.7));
		$rline[$i]=$j;
		$fline[$i]=$j;
		$i++;
		}
	}
$j++;
}

if($i==0) { printf("## No Source found in the NVSS catlog!!! ##\n<br>");exit 0;}

# sort by radius
#@index = sort{ $radius[ $a ] <=> $radius[ $b ] } 0 .. $#rline;
#@radius=@radius[@index];
#@rline=@rline[@index];


# print sorted(radius) list and angular distance
#  printf("Sorted by radial distance<br>\n");
#$i=0;
#foreach (@rline) { 
#	$i++;
#	@t=split(/ +/,$list[$_]);
#
#	$ra2=$t[1];
#	$dec2=$t[2];
#
#	# calculate distabce between ra and dec
#	$d1=radtodeg(abs($ra1-$ra2)); if($d1 > 180) { $d1=360-$d1; }
#	$d2=radtodeg(abs($dec1-$dec2));
#
#	# calculate the angular distance between target and cal 
#	$d=sin($dec1)*sin($dec2)+cos($dec1)*cos($dec2)*cos($ra1-$ra2);
#	$d=acos($d); $d=radtodeg($d);
#
#	# printf("Nvss# %d: %-2s %-2s %-5s %-3s %-2s %-5s %-6s (R: %.1fd, RA: %.1fd, Dec: %.1fd)<br>\n",$i,$t[3],$t[4],$t[5],$t[6],$t[7],$t[8],$t[9],$d,$d1,$d2);
#} 

# sorted by flux 
@index = sort{ $flux[ $a ] <=> $flux[ $b ] } 0 .. $#fline;
@flux=@flux[@index];
@fline=@fline[@index];
@fline = reverse @fline;

# print sorted(flux) list and angular distance
# printf("Sorted by flux value<br>\n");
$i=0;
foreach (@fline) { 
	$i++;
	if($i > 100) { last; } 
	@t=split(/ +/,$list[$_]);

	$ra2=$t[1];
	$dec2=$t[2];

	# calculate distance between ra and dec
	$d1=radtodeg(abs($ra1-$ra2)); if($d1 > 180) { $d1=360-$d1; }
	$d2=radtodeg(abs($dec1-$dec2));

	# calculate the angular distance between target and cal 
	$d=sin($dec1)*sin($dec2)+cos($dec1)*cos($dec2)*cos($ra1-$ra2);
	$d=acos($d); $d=radtodeg($d);

	printf("Nvss# %d: %-2s %-2s %-5s %-3s %-2s %-5s (flux: %-6s mJy, R: %2.1f deg.)<br>\n",$i,$t[3],$t[4],$t[5],$t[6],$t[7],$t[8],$t[9],$d);
} 
# end


sub ratorad{
	use strict;
	my($ra)=@_;
	my(@t)=split(/([hms])/,$ra);
	my($deg)=$t[0]*3600.0+$t[2]*60.0+$t[4];
	$deg=$deg/240.0;
	my($rad)=$deg/57.29577951308232;
	return($rad);
	}

sub dectorad{
	use strict;
        my($dec)=@_;
        my($deg,$rad,$sign,@t);
        @t=split(//,$dec);
        $sign="+";
        if($t[0] eq "-") { $sign="-";}
        @t=split(/([d'"])/,$dec);
        if($sign eq "+"){$deg=$t[0]+$t[2]/60.0+$t[4]/3600;}
        if($sign eq "-"){$deg=$t[0]-$t[2]/60.0-$t[4]/3600;}
        $rad=$deg/57.29577951308232;
        return($rad);
        }


sub radtodeg{ 
	use strict;
	my($rad) = @_;
	my($deg)=$rad*57.29577951308232;
	return $deg;
	}

sub acos { atan2( sqrt(1 - $_[0] * $_[0]), $_[0] ) }
sub asin { atan2($_[0], sqrt(1 - $_[0] * $_[0])) }