diff -Nru astap-cli-2022.03.01/astap_command_line_iOS_arm64.lpi astap-cli-2022.12.09/astap_command_line_iOS_arm64.lpi
--- astap-cli-2022.03.01/astap_command_line_iOS_arm64.lpi 1970-01-01 00:00:00.000000000 +0000
+++ astap-cli-2022.12.09/astap_command_line_iOS_arm64.lpi 2022-12-10 01:38:10.000000000 +0000
@@ -0,0 +1,76 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff -Nru astap-cli-2022.03.01/astap_command_line.lpi astap-cli-2022.12.09/astap_command_line.lpi
--- astap-cli-2022.03.01/astap_command_line.lpi 2022-03-02 22:39:48.000000000 +0000
+++ astap-cli-2022.12.09/astap_command_line.lpi 2022-12-10 01:38:10.000000000 +0000
@@ -1,16 +1,16 @@
-
+
+
-
@@ -29,7 +29,6 @@
-
@@ -59,6 +58,7 @@
+
@@ -68,7 +68,7 @@
-
+
@@ -78,6 +78,9 @@
+
+
+
diff -Nru astap-cli-2022.03.01/astap_command_line.lpr astap-cli-2022.12.09/astap_command_line.lpr
--- astap-cli-2022.03.01/astap_command_line.lpr 2022-03-02 22:39:48.000000000 +0000
+++ astap-cli-2022.12.09/astap_command_line.lpr 2022-12-10 01:38:10.000000000 +0000
@@ -54,27 +54,27 @@
if ((paramcount=0) or (hasOption('h','help'))) then
begin
writeln(
- 'ASTAP astrometric solver, '+version+#10+
+ 'ASTAP astrometric solver version CLI-'+astap_version+#10+
'(C) 2018, 2022 by Han Kleijn. License MPL 2.0, Webpage: www.hnsky.org'+#10+
'Usage:'+#10+
'-f filename {fits, tiff, png, jpg files}'+#10+
- '-f stdin {read raw image from stdin}'+#10+
'-r radius_area_to_search[degrees]'+#10+ {changed}
'-z downsample_factor[0,1,2,3,4] {Downsample prior to solving. 0 is auto}'+#10+
'-fov diameter_field[degrees]'+#10+ {changed}
'-ra center_right ascension[hours]'+#10+
'-spd center_south_pole_distance[degrees]'+#10+
- '-s max_number_of_stars {default 500}'+#10+
+ '-s max_number_of_stars {default 500, 0 is auto}'+#10+
'-t tolerance {default 0.007}'+#10+
'-m minimum_star_size["] {default 1.5}'+#10+
- '-speed mode[auto/slow] {Slow is forcing small search steps to improve detection.}'+#10+
+ '-check apply[y/n] {Apply check pattern filter prior to solving. Use for raw OSC images only when binning is 1x1}' +#10+
+ '-speed mode[auto/slow] {Slow is forcing reading a larger area from the star database (more overlap) to improve detection}'+#10+
'-o file {Name the output files with this base path & file name}'+#10+
'-d path {specify a path to the star database}'+#10+
'-analyse snr_min {Analyse only and report median HFD and number of stars used}'+#10+
'-extract snr_min {As -analyse but additionally write a .csv file with the detected stars info}'+#10+
'-log {Write the solver log to file}'+#10+
'-progress {Log all progress steps and messages}'+#10+
- '-update {update the FITS header with the found solution. Jpg, png, tiff will be written as fits}' +#10+
+ '-update {update the FITS header with the found solution. Jpeg, png, tiff will be written as fits}' +#10+
'-wcs {Write a .wcs file in similar format as Astrometry.net. Else text style.}' +#10+
'Preference will be given to the command line values.'
);
@@ -96,16 +96,7 @@
if filespecified then
begin
filename2:=GetOptionValue('f');
- stdin_mode:=filename2='stdin';
- if stdin_mode=false then {file mode}
- file_loaded:=load_image2 {load file first to give commandline parameters later priority}
- else
- begin
- memo1.clear;{prepare for some info wcs file}
- memo1.add(head1[0]);{add SIMPLE for case option -update is used}
- memo1.add(head1[27]);{add the END to memo1 for stdin and prevent runtime error since all data is inserted for END}
- file_loaded:=read_stdin_data;
- end;
+ file_loaded:=load_image2; {load file first to give commandline parameters later priority}
if file_loaded=false then errorlevel:=16;{error file loading}
end
else
@@ -131,10 +122,11 @@
{else dec from fits header}
if hasoption('z') then downsample_for_solving1:=strtoint(GetOptionValue('z'));
- if hasoption('s') then max_stars1:=GetOptionValue('s');
+ if hasoption('s') then max_stars:=strtoint(GetOptionValue('s'));
if hasoption('t') then quad_tolerance1:=GetOptionValue('t');
if hasoption('m') then min_star_size1:=GetOptionValue('m');
if hasoption('speed') then force_oversize1:=pos('slow',GetOptionValue('speed'))<>0;
+ if hasoption('check') then check_pattern_filter1:=('y'=GetOptionValue('check'));
extractspecified:=hasoption('extract');
analysespecified:=hasoption('analyse');
@@ -159,7 +151,7 @@
if hasoption('d') then
database_path:=GetOptionValue('d')+DirectorySeparator; {specify a different database path}
- if ((file_loaded) and (solve_image(img_loaded,true {get hist}) )) then {find plate solution, filename2 extension will change to .fit}
+ if ((file_loaded) and (solve_image(img_loaded ))) then {find plate solution, filename2 extension will change to .fit}
begin
if hasoption('o') then filename2:=GetOptionValue('o');{change file name for .ini file}
write_ini(true);{write solution to ini file}
@@ -168,8 +160,7 @@
if hasoption('update') then
begin
- if fits_file_name(filename2) then
- SaveFITSwithupdatedheader1
+ if fits_file_name(filename2) then SaveFITSwithupdatedheader1
else
save_fits16bit(img_loaded,ChangeFileExt(filename2,'.fits'));{save original png,tiff jpg to 16 bits fits file}
end;
@@ -244,7 +235,7 @@
var
Application: Tastap;
-{$R *.res}
+//{$R *.res}
begin
Application:=Tastap.Create(nil);
diff -Nru astap-cli-2022.03.01/astap_command_line_win10_arm32.lpi astap-cli-2022.12.09/astap_command_line_win10_arm32.lpi
--- astap-cli-2022.03.01/astap_command_line_win10_arm32.lpi 1970-01-01 00:00:00.000000000 +0000
+++ astap-cli-2022.12.09/astap_command_line_win10_arm32.lpi 2022-12-10 01:38:10.000000000 +0000
@@ -0,0 +1,86 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff -Nru astap-cli-2022.03.01/astap_command_line_Win11_arm64.lpi astap-cli-2022.12.09/astap_command_line_Win11_arm64.lpi
--- astap-cli-2022.03.01/astap_command_line_Win11_arm64.lpi 1970-01-01 00:00:00.000000000 +0000
+++ astap-cli-2022.12.09/astap_command_line_Win11_arm64.lpi 2022-12-10 01:38:10.000000000 +0000
@@ -0,0 +1,86 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff -Nru astap-cli-2022.03.01/_compile_all.sh astap-cli-2022.12.09/_compile_all.sh
--- astap-cli-2022.03.01/_compile_all.sh 2022-03-02 22:39:48.000000000 +0000
+++ astap-cli-2022.12.09/_compile_all.sh 2022-12-10 01:38:10.000000000 +0000
@@ -1,23 +1,41 @@
-/home/h/fpcupdeluxe2/lazarus/lazbuild astap_command_line_linux_arm64.lpi
+rm astap_cli
+
+/home/h/fpcupdeluxe/lazarus/lazbuild astap_command_line_linux_arm64.lpi
zip astap_command-line_version_Linux_arm64.zip astap_cli
+rm astap_cli
+
-/home/h/fpcupdeluxe2/lazarus/lazbuild astap_command_line_linux_armhf.lpi
+/home/h/fpcupdeluxe/lazarus/lazbuild astap_command_line_linux_armhf.lpi
zip astap_command-line_version_Linux_armhf.zip astap_cli
+rm astap_cli
-/home/h/fpcupdeluxe2/lazarus/lazbuild astap_command_line_darwin_M1.lpi
+
+/home/h/fpcupdeluxe/lazarus/lazbuild astap_command_line_darwin_M1.lpi
zip astap_command-line_version_macOS_M1.zip astap_cli
+rm astap_cli
+
-/home/h/fpcupdeluxe2/lazarus/lazbuild astap_command_line_darwin_X86_64.lpi
+/home/h/fpcupdeluxe/lazarus/lazbuild astap_command_line_darwin_X86_64.lpi
zip astap_command-line_version_macOS_X86_64.zip astap_cli
+rm astap_cli
-/home/h/fpcupdeluxe2/lazarus/lazbuild astap_command_line_android_armhf.lpi
-zip astap_command-line_version_android_armhf.zip astap_cli
+/home/h/fpcupdeluxe/lazarus/lazbuild astap_command_line_android_armhf.lpi
+zip astap_command-line_version_Android_armhf.zip astap_cli
+rm astap_cli
-/home/h/fpcupdeluxe2/lazarus/lazbuild astap_command_line_android_arm64.lpi
-zip astap_command-line_version_android_aarch64.zip astap_cli
+/home/h/fpcupdeluxe/lazarus/lazbuild astap_command_line_android_arm64.lpi
+zip astap_command-line_version_Android_aarch64.zip astap_cli
+rm astap_cli
-/home/h/fpcupdeluxe2/lazarus/lazbuild astap_command_line_linux.lpi
+/home/h/fpcupdeluxe/lazarus/lazbuild astap_command_line_linux.lpi
zip astap_command-line_version_Linux_amd64.zip astap_cli
+rm astap_cli
+
+
+/home/h/fpcupdeluxe/lazarus/lazbuild astap_command_line_Win11_arm64.lpi
+zip astap_command-line_version_Win11_arm64 astap_cli.exe
+rm astap_cli.exe
+
diff -Nru astap-cli-2022.03.01/debian/changelog astap-cli-2022.12.09/debian/changelog
--- astap-cli-2022.03.01/debian/changelog 2022-03-02 21:30:51.000000000 +0000
+++ astap-cli-2022.12.09/debian/changelog 2022-12-09 22:30:51.000000000 +0000
@@ -1,3 +1,15 @@
+astap-cli (2022.12.09-1) unstable; urgency=medium
+
+ [ Debian Janitor ]
+ * Use secure URI in Homepage field.
+
+ [ Thorsten Alteholz ]
+ * new upstream version
+ * debian/control: bump standard to 4.6.1 (no changes)
+ * debian/copyright: update license
+
+ -- Thorsten Alteholz Fri, 09 Dec 2022 23:30:51 +0100
+
astap-cli (2022.03.01-1) unstable; urgency=medium
* new upstream version
diff -Nru astap-cli-2022.03.01/debian/control astap-cli-2022.12.09/debian/control
--- astap-cli-2022.03.01/debian/control 2022-03-02 21:30:51.000000000 +0000
+++ astap-cli-2022.12.09/debian/control 2022-12-09 22:30:51.000000000 +0000
@@ -10,8 +10,8 @@
, lcl
, lcl-utils
, lcl-nogui
-Standards-Version: 4.6.0
-Homepage: http://www.hnsky.org/astap.htm
+Standards-Version: 4.6.1
+Homepage: https://www.hnsky.org/astap.htm
Rules-Requires-Root: binary-targets
Vcs-Browser: https://salsa.debian.org/debian-astro-team/astap-cli
Vcs-Git: https://salsa.debian.org/debian-astro-team/astap-cli.git
diff -Nru astap-cli-2022.03.01/debian/copyright astap-cli-2022.12.09/debian/copyright
--- astap-cli-2022.03.01/debian/copyright 2022-03-02 21:30:51.000000000 +0000
+++ astap-cli-2022.12.09/debian/copyright 2022-12-09 22:30:51.000000000 +0000
@@ -3,24 +3,14 @@
Source: https://sourceforge.net/p/astap-program/code/ci/default/tree/
Files: *
-Copyright: 2017, 2018, 2021 by Han Kleijn, www.hnsky.org
-License: LGPL-3+
+Copyright: 2017, 2018, 2021, 2022 by Han Kleijn, www.hnsky.org
+License: MPL-2.0
Files: debian/*
-Copyright: 2021 Thorsten Alteholz
-License: LGPL-3+
+Copyright: 2021-2022 Thorsten Alteholz
+License: MPL-2.0
-License: LGPL-3+
- This program is free software: you can redistribute it and/or modify it under
- the terms of the GNU Lesser General Public License (LGPL) as published by the
- Free Software Foundation, either version 3 of the License, or(at your option)
- any later version.
- .
- This program is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
- FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
- details.
- .
- On Debian systems, the complete text of the GNU Lesser General Public License
- can be found in `/usr/share/common-licenses/LGPL-3`.
+License: MPL-2.0
+ On Debian systems, the complete text of the MPL can be found in
+ `/usr/share/common-licenses/MPL-2.0`.
diff -Nru astap-cli-2022.03.01/debian/patches/series astap-cli-2022.12.09/debian/patches/series
--- astap-cli-2022.03.01/debian/patches/series 2022-03-02 21:30:51.000000000 +0000
+++ astap-cli-2022.12.09/debian/patches/series 2022-12-09 22:30:51.000000000 +0000
@@ -1 +1 @@
-spelling.patch
+#spelling.patch
diff -Nru astap-cli-2022.03.01/debian/patches/spelling.patch astap-cli-2022.12.09/debian/patches/spelling.patch
--- astap-cli-2022.03.01/debian/patches/spelling.patch 2022-03-02 21:30:51.000000000 +0000
+++ astap-cli-2022.12.09/debian/patches/spelling.patch 1970-01-01 00:00:00.000000000 +0000
@@ -1,14 +0,0 @@
-Description: spelling error found by lintian
-Index: astap-cli_2022.03.01.orig/unit_command_line_star_database.pas
-===================================================================
---- astap-cli_2022.03.01.orig.orig/unit_command_line_star_database.pas 2022-03-02 23:45:10.871404169 +0100
-+++ astap-cli_2022.03.01.orig/unit_command_line_star_database.pas 2022-03-02 23:45:10.867404168 +0100
-@@ -2762,7 +2762,7 @@
- begin
- if fileexists( database_path+'w08_0101.001') then begin name_database:='w08';database_type:=001; end
- else
-- memo2_message('Could not find w08 star database. Will try with an other database.');
-+ memo2_message('Could not find w08 star database. Will try with another database.');
- end
- else
- if fileexists( database_path+'h18_0101.1476') then begin name_database:='h18'; end
diff -Nru astap-cli-2022.03.01/linux_scripts/build_astap.sh astap-cli-2022.12.09/linux_scripts/build_astap.sh
--- astap-cli-2022.03.01/linux_scripts/build_astap.sh 2022-03-02 22:39:48.000000000 +0000
+++ astap-cli-2022.12.09/linux_scripts/build_astap.sh 2022-12-10 01:38:10.000000000 +0000
@@ -1,4 +1,4 @@
-/home/h/fpcupdeluxe2/lazarus/lazbuild ~/astap.fpc/astap_linux.lpi
+/home/h/fpcupdeluxe/lazarus/lazbuild ~/astap.fpc/astap_linux.lpi
cp ~/astap.fpc/astap ~/astap_install/astap_amd64/opt/astap
cd ~/astap_install
@@ -12,13 +12,16 @@
tar -czvf astap_amd64.tar.gz /opt/astap/astap /opt/astap/astap.ico /opt/astap/astap.ico /opt/astap/copyright.txt /opt/astap/deep_sky.csv /opt/astap/variable_stars.csv /usr/share/applications/ASTAP.desktop /usr/local/bin/astap /opt/astap/dcraw-astap /opt/astap/unprocessed_raw-astap
-/home/h/fpcupdeluxe2/lazarus/lazbuild ~/astap.fpc/astap_linux_qt5.lpi
+/home/h/fpcupdeluxe/lazarus/lazbuild ~/astap.fpc/astap_linux_qt5.lpi
sudo cp ~/astap.fpc/astap /opt/astap
tar -czvf astap_amd64_qt5.tar.gz /opt/astap/astap /opt/astap/astap.ico /opt/astap/astap.ico /opt/astap/copyright.txt /opt/astap/deep_sky.csv /opt/astap/variable_stars.csv /usr/share/applications/ASTAP.desktop /usr/local/bin/astap /opt/astap/dcraw-astap /opt/astap/unprocessed_raw-astap
-
-
-
-
-
+/home/h/fpcupdeluxe/lazarus/lazbuild ~/astap.fpc/astap_linux_cross_compile_to_Darwin_M1.lpi
+sudo cp ~/astap.fpc/astap
+zip astap_mac_M1.zip ~/astap.fpc/astap
+
+
+/home/h/fpcupdeluxe/lazarus/lazbuild ~/astap.fpc/astap_linux_cross_compile_to_Darwin_X86_64.lpi
+sudo cp ~/astap.fpc/astap
+zip astap_mac_X86_64.zip ~/astap.fpc/astap
diff -Nru astap-cli-2022.03.01/unit_command_line_general.pas astap-cli-2022.12.09/unit_command_line_general.pas
--- astap-cli-2022.03.01/unit_command_line_general.pas 2022-03-02 22:39:48.000000000 +0000
+++ astap-cli-2022.12.09/unit_command_line_general.pas 2022-12-10 01:38:10.000000000 +0000
@@ -7,7 +7,6 @@
uses
Classes, SysUtils,
math,
- iostream, {for using stdin for data}
FPImage,
strutils,
fpreadTIFF, {all part of fcl-image}
@@ -16,16 +15,16 @@
var {################# initialised variables #########################}
- stdin_mode : boolean=false;{file send via stdin}
- version: string=' CLI-2022-2-04';
+ astap_version: string='2022.11.11';
ra1 : string='0';
dec1 : string='0';
search_fov1 : string='0';{search FOV}
radius_search1 : string='180';
- max_stars1: string='500';
+ max_stars: integer=500;
quad_tolerance1: string='0.007';
min_star_size1: string='1.5'{arcsec};
- force_oversize1 : boolean=false;
+ force_oversize1 : boolean=false;
+ check_pattern_filter1 : boolean=false;
commandline_log : boolean=false;{file log request in command line}
solve_show_log : boolean=false;{log all progress steps}
errorlevel : integer=0;{report errors when shutdown}
@@ -42,27 +41,6 @@
var
memo1,memo2 :tstrings; {settings for save and loading}
- type
- timgbackup = record
- crpix1 : double;{could be modified by crop}
- crpix2 : double;
- crval1 : double;
- crval2 : double;
- crota1 : double;{for 90 degrees rotate}
- crota2 : double;
- cdelt1 : double;
- cdelt2 : double;
- cd1_1 : double;
- cd1_2 : double;
- cd2_1 : double;
- cd2_2 : double;
- header : string;
- img : array of array of array of single;
- end;
-
- var
- img_backup : array of timgbackup;{dynamic so memory can be freed}
-
user_path : string;{c:\users\name\appdata\local\astap or ~/home/.config/astap}
img_loaded,img_temp,img_dark,img_flat,img_bias,img_average,img_variance,img_buffer,img_final : image_array;
filename2: string;
@@ -81,7 +59,7 @@
his_mean,noise_level : array[0..2] of integer;
esc_pressed, fov_specified {, last_extension }: boolean;
star_level : double;
- exposure : double;
+ exposure,focallen,equinox : double;
cblack, cwhite,
gain :double; {from FITS}
date_obs : string;
@@ -89,23 +67,20 @@
datamin_org, datamax_org :double;
warning_str : string;{for solver}
+ xpixsz,ypixsz: double;//Pixel Width in microns (after binning)
+ ra0,dec0 : double; {plate center values}
+ cdelt1,cdelt2: double;{deg/pixel for x}
+ naxis : integer;{number of dimensions}
+ naxis3 : integer;{number of colors}
const
hist_range {range histogram 255 or 65535 or streched} : integer=255;
width2 : integer=100;
height2: integer=100;{do not use reserved word height}
- naxis : integer=2;{number of dimensions}
- naxis3 : integer=1;{number of colors}
fits_file: boolean=false;
- ra0 : double=0;
- dec0: double=0; {plate center values}
crpix1: double=0;{reference pixel}
crpix2: double=0;
- cdelt1: double=0;{deg/pixel for x}
- cdelt2: double=0;
- xpixsz: double=0;//Pixel Width in microns (after binning)
- ypixsz: double=0;//Pixel height in microns (after binning)
- focallen: double=0;
+
const bufwide=1024*120;{buffer size in bytes}
@@ -172,11 +147,9 @@
function prepare_ra(rax:double; sep:string):string; {radialen to text, format 24: 00 00.0 }
function prepare_dec(decx:double; sep:string):string; {radialen to text, format 90d 00 00}
procedure new_to_old_WCS;{convert new style FITsS to old style}
-procedure old_to_new_WCS;{ convert old WCS to new}
procedure write_astronomy_wcs;
procedure write_ini(solution:boolean);{write solution to ini file}
function load_image2 : boolean; {load fits or PNG, BMP, TIF}
-function read_stdin_data: boolean; {reads via stdin a raw image based on the INDI standard}
procedure SaveFITSwithupdatedheader1;
function save_fits16bit(img: image_array;filen2:ansistring): boolean;{save to 16 fits file}
@@ -372,31 +345,57 @@
end;
-procedure old_to_new_WCS;{ convert old WCS to new}
-var
- sign : integer;
+procedure new_to_old_WCS;{convert new style FITS to old style, revison 2022}
begin
- cd1_1:=cdelt1*cos(crota2*pi/180); {note 2013 should be crota1 if skewed}
- if cdelt1>=0 then sign:=+1 else sign:=-1;
- cd1_2:=abs(cdelt2)*sign*sin(crota2*pi/180);{note 2013 should be crota1 if skewed}
- if cdelt2>=0 then sign:=+1 else sign:=-1;
- cd2_1:=-abs(cdelt1)*sign*sin(crota2*pi/180);
- cd2_2:= cdelt2*cos(crota2*pi/180);
-end;
+ // https://www.aanda.org/articles/aa/full/2002/45/aah3860/aah3860.right.html
+ // Representations of World Coordinates in FITS paper II aah3860
+ // formula 191
+ if cd2_1>0 then crota1:=arctan2(-cd2_1,-cd1_1)*180/pi
+ else
+ if cd2_1<0 then crota1:=arctan2(+cd2_1,+cd1_1)*180/pi
+ else
+ crota1:=0;
-procedure new_to_old_WCS;{convert new style FITsS to old style}
-var
- sign : integer;
-begin
- { convert to old WCS. Based on draft 1988 , do not use conversion article Alain Klotz, give sometimes zero CROTA}
- if (cd1_1*cd2_2-cd1_2*cd2_1)>=0 then sign:=+1 else sign:=-1;
+ if cd1_2>0 then crota2:=arctan2(-cd1_2,cd2_2)*180/pi
+ else
+ if cd1_2<0 then crota2:=arctan2(cd1_2,-cd2_2)*180/pi
+ else
+ crota2:=0;
- cdelt1:=sqrt(sqr(cd1_1)+sqr(cd2_1))*sign;{if no old wcs header use cd2_2 of new WCS style for pixel size}
- cdelt2:=sqrt(sqr(cd1_2)+sqr(cd2_2));{if no old wcs header use cd2_2 of new WCS style for pixel size}
+ // https://www.aanda.org/articles/aa/full/2002/45/aah3860/aah3860.right.html
+ // Representations of World Coordinates in FITS paper II aah3860
+ // Formula 193 improved for crota close to or equal to +90 or -90 degrees
+ // Calculate cdelt1, cdelt2 values using the longest side of the triangle
+ if abs(cd1_1)>abs(cd2_1) then
+ begin
+ cdelt1:=+cd1_1/cos(crota1*pi/180);
+ cdelt2:=+cd2_2/cos(crota2*pi/180);
+ end
+ else
+ begin
+ cdelt1:=+cd2_1/sin(crota1*pi/180);
+ cdelt2:=-cd1_2/sin(crota2*pi/180);
+ end;
- crota1:= +arctan2(sign*cd1_2,cd2_2)*180/pi;
- crota2:= -arctan2(cd2_1,sign*cd1_1)*180/pi; // crota2old := (atn_2(sign*cd1_1,cd2_1)-pi/2)*180/pi;
+ //Solutions for CROTA2 come in pairs separated by 180degr. The other solution is obtained by subtracting 180 from CROTA2 and negating CDELT1 and CDELT2.
+ //While each solution is equally valid, if one makes CDELT1 < 0 and CDELT2 > 0 then it would normally be the one chosen.
+ if cdelt2<0 then //CDELT2 is always kept positive and if not the solution is flipped by negating both CDELT2, CDELT2 and shifting the angle 180 degrees. So if the image is flipped the solution is reporting "flipped horizontal" and not an equivalent "flipped vertical".
+ begin
+ if crota2<0 then
+ begin
+ crota2:=crota2+180;
+ crota1:=crota1+180;
+ end
+ else
+ begin
+ crota2:=crota2-180;
+ crota1:=crota1-180;
+ end;
+
+ cdelt2:=-cdelt2;
+ cdelt1:=-cdelt1;
+ end;//make cdelt2 always positive
end;
@@ -439,19 +438,7 @@
var
f: text;
begin
- if stdin_mode then {raw file send via stdin}
- begin
- {$IFDEF MSWINDOWS}
- if isConsole then {console available}
- AssignFile(f,'')
- else {no console available, compiler option -WH is checked}
- assignfile(f,filename2+'.ini'); {filename2 could be long due to -o option}
- {$ELSE}
- AssignFile(f,''); {write to console, unix and Darwin}
- {$ENDIF}
- end
- else
- assignfile(f,ChangeFileExt(filename2,'.ini'));
+ assignfile(f,ChangeFileExt(filename2,'.ini'));
rewrite(f);
if solution then
begin
@@ -534,7 +521,7 @@
end;
- Reader := TReader.Create (theFile3,$4000);{number of hnsky records}
+ Reader := TReader.Create (theFile3,$60000);// 393216 byte buffer
{thefile3.size-reader.position>sizeof(hnskyhdr) could also be used but slow down a factor of 2 !!!}
I:=0;
reader_position:=0;
@@ -690,11 +677,58 @@
end;
+procedure precession_jnow_to_J2000(equinox : double; var ra1,dec1 : double); {simple precession correction, new Meeus chapter precession formula 20.1}
+var
+ t,dra,ddec,m,n,n2 : double;
+begin
+ t:=(equinox-2000)/100;{time in julian centuries since j2000 }
+ m:=3.07496+0.00186*t;{seconds}
+ n:=1.33621-0.00057*t; {seconds}
+ n2:=20.0431-0.0085*t;{arcsec}
+ dra:=(m + n *sin(ra1)*tan(dec1))*pi/(3600*12);{yearly ra drift in radians}
+ ddec:=n2*cos(ra1)*pi/(3600*180); {yearly dec drift in radians}
+ ra1:=ra1-(dra*t*100);{multiply with number of years is t*100. Subtract because we go back to J2000}
+ dec1:=dec1-(ddec*t*100);
+end;
+
+
+procedure reset_fits_global_variables; {reset the global variable}
+begin
+ ra0:=0;
+ dec0:=0;
+ ra_mount:=99999;
+ dec_mount:=99999;
+ cdelt1:=0;
+ cdelt2:=0;
+ xpixsz:=0;
+ ypixsz:=0;
+ focallen:=0;
+ cd1_1:=0;{just for the case it is not available}
+ cd1_2:=0;{just for the case it is not available}
+ cd2_1:=0;{just for the case it is not available}
+ cd2_2:=0;{just for the case it is not available}
+ xbinning:=1;{normal}
+ ybinning:=1;
+ ra1:='';
+ dec1:='';
+ equinox:=2000;
+
+ naxis:=1;
+ naxis3:=1;
+ datamin_org:=0;
+ datamax_org:=$FFFF;
+
+ date_obs:='';
+ exposure:=0;
+ gain:=999;{assume no data available}
+end;
+
+
function load_fits(filen:string;out img_loaded2: image_array): boolean;{load fits file}
var
header : array[0..2880] of ansichar;
i,j,k,error3,naxis1, reader_position : integer;
- dummy,scale : double;
+ dummy : double;
col_float,bscale,measured_max,scalefactor : single;
bzero : integer;{zero shift. For example used in AMT, Tricky do not use int64, maxim DL writes BZERO value -2147483647 as +2147483648 !! }
aline : ansistring;
@@ -740,22 +774,11 @@
begin
result:='';
r:=I+11;{pos12, single quotes should for fix format should be at position 11 according FITS standard 4.0, chapter 4.2.1.1}
- // while ((header[r-1]<>#39) and (r<=I+20)) do inc(r);
- repeat
- result:=result+header[r];
- inc(r);
- until ((header[r]=#39){last quote} or (r>=I+79));{read string up to position 80 equals 79}
- end;
-
- Function get_as_string:string;{read float as string values. Universal e.g. for latitude and longitude which could be either string or float}
- var r: integer;
- begin
- result:='';
- r:=I+11;{pos12, single quotes should for fix format should be at position 11 according FITS standard 4.0, chapter 4.2.1.1}
- repeat
+ while ((header[r]<>#39){last quote} and (r=I+30));{read string up to position 30}
+ end;
end;
begin
@@ -774,38 +797,14 @@
memo1.clear;{clear memo for new header}
- Reader := TReader.Create (theFile3,500*2880);{number of records. Buffer but not speed difference between 6*2880 and 1000*2880}
- {thefile3.size-reader.position>sizeof(hnskyhdr) could also be used but slow down a factor of 2 !!!}
+ Reader := TReader.Create (theFile3,128*2880);{number of records. 128*2880 is 2% faster then 8* 2880}
{Reset variables for case they are not specified in the file}
- crota2:=99999;{just for the case it is not available, make it later zero}
- crota1:=99999;
- ra0:=0;
- dec0:=0;
- ra_mount:=99999;
- dec_mount:=99999;
- cdelt1:=0;
- cdelt2:=0;
- scale:=0; {SGP files}
- xpixsz:=0;
- ypixsz:=0;
- focallen:=0;
- cd1_1:=0;{just for the case it is not available}
- cd1_2:=0;{just for the case it is not available}
- cd2_1:=0;{just for the case it is not available}
- cd2_2:=0;{just for the case it is not available}
- xbinning:=1;{normal}
- ybinning:=1;
- ra1:='';
- dec1:='';
+ reset_fits_global_variables; {reset the global variable}
- naxis:=-1;
- naxis1:=0;
- naxis3:=1;
bzero:=0;{just for the case it is not available. 0.0 is the default according https://heasarc.gsfc.nasa.gov/docs/fcg/standard_dict.html}
bscale:=1;
- datamin_org:=0;
-
+ naxis1:=0;
measured_max:=0;
reader_position:=0;
@@ -814,7 +813,7 @@
I:=0;
try
reader.read(header[I],2880);{read file header, 2880 bytes}
- inc(reader_position,2880);
+ inc(reader_position,2880); {thefile3.size-reader.position>sizeof(hnskyhdr) could also be used but slow down a factor of 2 !!!}
if ((reader_position=2880) and (header[0]='S') and (header[1]='I') and (header[2]='M') and (header[3]='P') and (header[4]='L') and (header[5]='E') and (header[6]=' ')) then
begin
simple:=true;
@@ -858,6 +857,12 @@
height2:=naxis3;
naxis3:=1;
end;
+ if naxis3>3 then {panic, more then three colours. Program https://github.com/cbassa/stvid is storing the mean, st, max and argmax values of each pixel respectively from multiple files }
+ begin
+ naxis3:=1; {display only the first colour}
+ memo2_message('Warning more then three colours. Will use only the first one.');
+ end;
+
end;
end;
@@ -890,35 +895,29 @@
if ((header[i]='Y') and (header[i+1]='B') and (header[i+2]='I') and (header[i+3]='N') and (header[i+4]='N') and (header[i+5]='I')) then
ybinning:=round(validate_double);{binning}
- if ((header[i]='C') and (header[i+1]='R') and (header[i+2]='O') and (header[i+3]='T') and (header[i+4]='A')) then {crota2}
- begin
- if (header[i+5]='2') then crota2:=validate_double else {read double value}
- if (header[i+5]='1') then crota1:=validate_double;{read double value}
- end;
- if ((header[i]='C') and (header[i+1]='R') and (header[i+2]='P') and (header[i+3]='I') and (header[i+4]='X')) then {crpix1}
- begin
- if header[i+5]='1' then crpix1:=validate_double else{ref pixel for x}
- if header[i+5]='2' then crpix2:=validate_double; {ref pixel for y}
- end;
if ((header[i]='C') and (header[i+1]='D') and (header[i+2]='E') and (header[i+3]='L') and (header[i+4]='T')) then {cdelt1}
begin
if header[i+5]='1' then cdelt1:=validate_double else{deg/pixel for RA}
if header[i+5]='2' then cdelt2:=validate_double; {deg/pixel for DEC}
end;
if ( ((header[i]='S') and (header[i+1]='E') and (header[i+2]='C') and (header[i+3]='P') and (header[i+4]='I') and (header[i+5]='X')) or {secpix1/2}
+ ((header[i]='S') and (header[i+1]='C') and (header[i+2]='A') and (header[i+3]='L') and (header[i+4]='E') and (header[i+5]=' ')) or {SCALE value for SGP files}
((header[i]='P') and (header[i+1]='I') and (header[i+2]='X') and (header[i+3]='S') and (header[i+4]='C') and (header[i+5]='A')) ) then {pixscale}
begin
if cdelt2=0 then
begin cdelt2:=validate_double/3600; {deg/pixel for RA} cdelt1:=cdelt2; end; {no CDELT1/2 found yet, use alternative}
end;
+ if ((header[i]='E') and (header[i+1]='Q') and (header[i+2]='U') and (header[i+3]='I') and (header[i+4]='N') and (header[i+5]='O') and (header[i+6]='X')) then
+ equinox:=validate_double;
+
if ((header[i]='X') and (header[i+1]='P') and (header[i+2]='I') and (header[i+3]='X') and (header[i+4]='S') and (header[i+5]='Z')) then {xpixsz}
xpixsz:=validate_double;{Pixel Width in microns (after binning), maxim DL keyword}
if ((header[i]='Y') and (header[i+1]='P') and (header[i+2]='I') and (header[i+3]='X') and (header[i+4]='S') and (header[i+5]='Z')) then {xpixsz}
- ypixsz:=validate_double;{Pixel Width in microns (after binning), maxim DL keyword}
+ ypixsz:=validate_double;{Pixel Width in microns (after binning), maxim DL keyword}
-// if ((header[i]='F') and (header[i+1]='O') and (header[i+2]='C') and (header[i+3]='A') and (header[i+4]='L') and (header[i+5]='L')) then {focall}
-// focallen:=validate_double;{Focal length of telescope in mm, maxim DL keyword}
+ if ((header[i]='F') and (header[i+1]='O') and (header[i+2]='C') and (header[i+3]='A') and (header[i+4]='L') and (header[i+5]='L')) then {focall}
+ focallen:=validate_double;{Focal length of telescope in mm, maxim DL keyword}
if ((header[i]='C') and (header[i+1]='R') and (header[i+2]='V') and (header[i+3]='A') and (header[i+4]='L')) then {crval1/2}
begin
@@ -961,7 +960,6 @@
if ((header[i+2]='2') and (header[i+3]='_') and (header[i+4]='2')) then cd2_2:=validate_double;
end;
- if ((header[i]='S') and (header[i+1]='C') and (header[i+2]='A') and (header[i+3]='L') and (header[i+4]='E')) then scale:=validate_double; {SCALE value for SGP files}
end; {image header}
@@ -987,21 +985,14 @@
nrbits:=24; {threat RGB fits as 2 dimensional with 24 bits data}
naxis3:=3; {will be converted while reading}
end;
- if ( ((cdelt1=0) or (crota2>=999)) and (cd1_1<>0)) then
- begin
- new_to_old_WCS;{ convert old WCS to new}
- end
- else
- if ((crota2<999) and (cd1_1=0) and(cdelt1<>0)) then {valid crota2 value}
- begin
- old_to_new_WCS;{ convert old WCS to new}
- end;
-
- if crota2>999 then crota2:=0;{not defined, set at 0}
- if crota1>999 then crota1:=crota2; {for case crota1 is not used}
if ((ra0<>0) or (dec0<>0)) then
begin
+ if equinox<>2000 then //e.g. in SharpCap
+ begin
+ precession_Jnow_to_J2000(equinox,ra0,dec0); {precession, from unknown equinox to J2000}
+ if dec_mount<999 then precession_Jnow_to_J2000(equinox,ra_mount,dec_mount); {precession, from unknown equinox to J2000}
+ end;
ra1:=prepare_ra(ra0,' ');
dec1:=prepare_dec(dec0,' ');
end
@@ -1012,17 +1003,17 @@
dec_text_to_radians( dec1,dec0,error1); {convert dec text to dec0 in radians}
end;
-
- if ((cd1_1=0) and (cdelt2=0)) then {no scale, try to fix it}
+ if cdelt2=0 then {simple code for astap-cli only}
begin
- if scale<>0 then {sgp file, use scale to find image dimensions}
- cdelt2:=scale/3600 {scale is in arcsec/pixel }
- else
+ if cd1_1=0 then {no scale, try to fix it}
+ begin
if ((focallen<>0) and (xpixsz<>0)) then
- cdelt2:=180/(pi*1000)*xpixsz/focallen; {use maxim DL key word. xpixsz is including binning}
+ cdelt2:=180/(pi*1000)*xpixsz/focallen; {use maxim DL key word. xpixsz is including binning}
+ end
+ else
+ cdelt2:=sqrt(sqr(cd1_2)+sqr(cd2_2));
end;
-
{############################## read image}
i:=round(bufwide/(abs(nrbits/8)));{check if buffer is wide enough for one image line}
if width2>i then
@@ -1406,6 +1397,137 @@
end;
+procedure read_keys_memo;{for tiff, header in the describtion decoding}
+var
+ key : string;
+ count1,index : integer;
+ error1 : boolean;
+
+ function read_float: double;
+ var
+ err: integer;
+ begin
+ val(copy(Memo1[index],11,20),result,err);
+ end;
+ function read_integer: integer;
+ var
+ err: integer;
+ begin
+ val(copy(Memo1[index],11,20),result,err);
+ end;
+ function read_string: string;
+ var
+ p1,p2 :integer;
+ begin
+ result:=copy(Memo1[index],11,80-11);
+ p1:=pos(char(39),result);
+ p2:=posex(char(39),result,p1+1);
+ if p2=0 then p2:=20;
+ result:=trim(copy(result,p1+1,p2-p1-1));{remove all spaces}
+ end;
+
+begin
+ {variables are already reset}
+ count1:=Memo1.Count-1-1;
+// ccd_temperature:=999;
+
+ index:=1;
+ while index<=count1 do {read keys}
+ begin
+ key:=copy(Memo1[index],1,9);
+ if key='CD1_1 =' then cd1_1:=read_float else
+ if key='CD1_2 =' then cd1_2:=read_float else
+ if key='CD2_1 =' then cd2_1:=read_float else
+ if key='CD2_2 =' then cd2_2:=read_float else
+
+ if key='CRVAL1 =' then ra0:=read_float*pi/180 {degrees -> radians} else
+ if key='CRVAL2 =' then dec0:=read_float*pi/180 else
+ if key='RA =' then
+ begin
+ ra_mount:=read_float*pi/180;{degrees -> radians}
+ if ra0=0 then ra0:=ra_mount; {ra telescope, read double value only if crval is not available}
+ end else
+ if key='DEC =' then
+ begin
+ dec_mount:=read_float*pi/180;
+ if dec0=0 then dec0:=dec_mount; {ra telescope, read double value only if crval is not available}
+ end else
+ if ((key='OBJCTRA =') and (ra_mount>=999)) {ra_mount value is unfilled, preference for keyword RA} then
+ begin
+ ra1:=read_string;{triggers an onchange event which will convert the string to ra_radians}
+ ra_mount:=ra_radians;{preference for keyword RA}
+ end else
+ if ((key='OBJCTDEC=') and (dec_mount>=999)) {dec_mount value is unfilled, preference for keyword DEC} then
+ begin
+ dec1:=read_string;{triggers an onchange event which will convert the string to dec_radians}
+ dec_mount:=dec_radians;
+ end else
+
+ if (key='XBINNING=') then xbinning:=read_integer else
+ if (key='YBINNING=') then ybinning:=read_integer else
+
+ if (key='FOCALLEN=') then focallen:=read_float else
+ if (key='XPIXSZ =') then xpixsz:=read_float else {pixelscale in microns}
+ if (key='YPIXSZ =') then ypixsz:=read_float else
+ if (key='CDELT2 =') then cdelt2:=read_float else {deg/pixel}
+ if (key='EQUINOX =') then equinox:=read_float else
+
+
+ if ((key='SECPIX2 =') or
+ (key='PIXSCALE=') or
+ (key='SCALE =')) then begin if cdelt2=0 then cdelt2:=read_float/3600; end {no cdelt1/2 found yet, use alternative, image scale arcseconds per pixel}
+ else
+
+ if key='DATE-OBS=' then date_obs:=read_string else
+
+ if index=1 then if key<>'BITPIX =' then begin Memo1.insert(index,'BITPIX = 16 / Bits per entry '); inc(count1); end;{data will be added later}
+ if index=2 then if key<>'NAXIS =' then begin Memo1.insert(index,'NAXIS = 2 / Number of dimensions ');inc(count1); end;{data will be added later}
+ if index=3 then if key<>'NAXIS1 =' then begin Memo1.insert(index,'NAXIS1 = 100 / length of x axis ');inc(count1); end;{data will be added later}
+ if index=4 then if key<>'NAXIS2 =' then begin Memo1.insert(index,'NAXIS2 = 100 / length of y axis ');inc(count1); end;{data will be added later}
+ if ((index=5) and (naxis>1)) then if key<>'NAXIS3 =' then
+ begin Memo1.insert(index,'NAXIS3 = 3 / length of z axis (mostly colors) ');inc(count1); end;
+ index:=index+1;
+ end;
+
+ if ((ra0<>0) or (dec0<>0)) then
+ begin
+ if equinox<>2000 then //e.g. in SharpCap
+ begin
+ precession_Jnow_to_J2000(equinox,ra0,dec0); {precession, from unknown equinox to J2000}
+ if dec_mount<999 then precession_Jnow_to_J2000(equinox,ra_mount,dec_mount); {precession, from unknown equinox to J2000}
+ end;
+ ra1:=prepare_ra(ra0,' ');{this will create Ra_radians for solving}
+ dec1:=prepare_dec(dec0,' ');
+ end
+ else
+ if ra1<>'' then
+ begin
+ ra_text_to_radians ( ra1 ,ra0,error1); {convert ra text to ra0 in radians}
+ dec_text_to_radians( dec1,dec0,error1); {convert dec text to dec0 in radians}
+ end;
+
+
+ { condition keyword to
+ if ra_mount>999 then objctra--->ra1.text--------------->ra_radians--->ra_mount
+ ra--->ra_mount if ra0=0 then ra_mount--->ra0
+ crval1--->ra0
+
+ if ra0<>0 then ra0--->ra1.text------------------->ra_radians}
+
+
+
+ if cdelt2=0 then {simple code for astap-cli only}
+ begin
+ if cd1_1=0 then {no scale, try to fix it}
+ begin
+ if ((focallen<>0) and (xpixsz<>0)) then
+ cdelt2:=180/(pi*1000)*xpixsz/focallen; {use maxim DL key word. xpixsz is including binning}
+ end
+ else
+ cdelt2:=sqrt(sqr(cd1_2)+sqr(cd2_2));
+ end;
+end;
+
function load_PPM_PGM_PFM(filen:string; var img_loaded2: image_array) : boolean;{load PPM (color),PGM (gray scale)file or PFM color}
var
@@ -1442,40 +1564,14 @@
end;
memo1.clear;{clear memo for new header}
- Reader := TReader.Create (theFile3,$4000);{number of hnsky records}
+ Reader := TReader.Create (theFile3,$60000);// 393216 byte buffer
{thefile3.size-reader.position>sizeof(hnskyhdr) could also be used but slow down a factor of 2 !!!}
{Reset variables}
- crota2:=99999;{just for the case it is not available, make it later zero}
- crota1:=99999;
- ra0:=0;
- dec0:=0;
- ra_mount:=99999;
- dec_mount:=99999;
- cdelt1:=0;
- cdelt2:=0;
- xpixsz:=0;
- ypixsz:=0;
- focallen:=0;
- cd1_1:=0;
- cd1_2:=0;
- cd2_1:=0;
- cd2_2:=0;
- date_obs:='';
-
- naxis:=1;
- naxis3:=1;
-
- xbinning:=1;{default}
- ybinning:=1;
- exposure:=0;
-
- gain:=999;{assume no data available}
+ reset_fits_global_variables; {reset the global variable}
I:=0;
reader_position:=0;
-
-
aline:='';
try
for i:=0 to 2 do begin reader.read(ch,1); aline:=aline+ch; inc(reader_position,1);end;
@@ -1567,9 +1663,7 @@
exit;
end; {should contain 255 or 65535}
- datamin_org:=0;
fits_file:=true;
-
cblack:=datamin_org;{for case histogram is not called}
cwhite:=datamax_org;
@@ -1698,7 +1792,6 @@
end;
-
function load_TIFFPNGJPEG(filen:string; var img_loaded2: image_array) : boolean;{load 8 or 16 bit TIFF, PNG, JPEG, BMP image}
var
i,j : integer;
@@ -1708,7 +1801,7 @@
tiff, png,jpeg,colour,saved_header : boolean;
ext,descrip : string;
begin
- naxis:=0; {0 dimensions}
+ naxis:=0; {0 dimensions for case failure}
result:=false; {assume failure}
tiff:=false;
jpeg:=false;
@@ -1775,6 +1868,9 @@
end;
{$ENDIF}
+ {Reset variables}
+ reset_fits_global_variables; {reset the global variable}
+
if colour=false then
begin
naxis:=2;
@@ -1788,31 +1884,10 @@
memo1.clear;{clear memo for new header}
- {Reset variables}
- crota2:=99999;{just for the case it is not available, make it later zero}
- crota1:=99999;
- ra0:=0;
- dec0:=0;
- ra_mount:=99999;
- dec_mount:=99999;
- cdelt1:=0;
- cdelt2:=0;
- xpixsz:=0;
- ypixsz:=0;
- focallen:=0;
- cd1_1:=0;
- cd1_2:=0;
- cd2_1:=0;
- cd2_2:=0;
- date_obs:='';
- xbinning:=1;{default}
- ybinning:=1;
-
{set data}
fits_file:=true;
nrbits:=16;
- datamin_org:=0;
- datamax_org:=$FFFF;
+
cblack:=datamin_org;{for case histogram is not called}
cwhite:=datamax_org;
@@ -1842,12 +1917,28 @@
begin
descrip:=image.Extra['TiffImageDescription']; {restore full header in TIFF !!!}
end;
+
for j:=0 to 10 do {create an header with fixed sequence}
if ((j<>5) or (naxis3<>1)) then {skip naxis3 for mono images}
memo1.add(head1[j]); {add lines to empthy memo1}
memo1.add(head1[27]); {add end}
if descrip<>'' then add_long_comment(descrip);{add TIFF describtion}
+ if copy(descrip,1,6)='SIMPLE' then {fits header included}
+ begin
+ memo1.text:=descrip;
+ read_keys_memo;
+ saved_header:=true;
+ end
+ else {no fits header in tiff file available}
+ begin
+ for j:=0 to 10 do {create an header with fixed sequence}
+ if ((j<>5) or (naxis3<>1)) then {skip naxis3 for mono images}
+ memo1.add(head1[j]); {add lines to empthy memo1}
+ memo1.add(head1[27]); {add end}
+ if descrip<>'' then add_long_comment(descrip);{add TIFF describtion}
+ end;
+
update_integer('BITPIX =',' / Bits per entry ' ,nrbits);
update_integer('NAXIS =',' / Number of dimensions ' ,naxis);{2 for mono, 3 for colour}
update_integer('NAXIS1 =',' / length of x axis ' ,width2);
@@ -1904,86 +1995,6 @@
end;
-function read_stdin_data: boolean; {reads via stdin a raw image based on the INDI standard}
-var
- count,format_type,i,h,bufferlength: Integer;
- rgbdummy : byteX3;
- x_longword : longword;
- x_single : single absolute x_longword;{for conversion 32 bit float}
- InputStream : TIOStream;
-begin
- result:=false;{assume failure}
- cd1_1:=0;{no solution}
-
- try
-
- InputStream := TIOStream.Create(iosInput);
- {read header}
- Count := InputStream.Read(format_type, 4);
-
- naxis3:=1; {assume mono}
- if format_type=$32574152 then nrbits:=16 {RAW2, 16 bit mono identifier}
- else
- if format_type=$31574152 then nrbits:=8 {RAW1, 8 bit mono}
- else
- if format_type=$33574152 then begin nrbits:=24; {24bit RGB = RAW3 = 0x33574152, rgb,rgb,rgb} naxis3:=3;end
- else
- if format_type=$34574152 then nrbits:=32 {RAW4, 32 bit mono}
- else
- if format_type=$66574152 then nrbits:=-32 {RAWf, -32 bit float mono}
- else
- exit;
-
- Count := InputStream.Read(width2, 4);
- Count := InputStream.Read(height2, 4);
-
- {read image data}
- setlength(img_loaded,naxis3,width2,height2);
- h:=0;
- bufferlength:=width2*(nrbits div 8);
- repeat
- count:=0;
- repeat {adapt reading to one image line}
- Count :=count + InputStream.Read(fitsBuffer[count],bufferlength-count );
- until count>=bufferlength;
- begin
- if nrbits=16 then
- for i:=0 to width2-1 do {fill array} img_loaded[0,i,h]:=fitsbuffer2[i]
- else
- if nrbits=8 then
- for i:=0 to width2-1 do {fill array} img_loaded[0,i,h]:=fitsbuffer[i]
- else
- if nrbits=24 then
- begin
- for i:=0 to width2-1 do {fill array}
- begin
- rgbdummy:=fitsbufferRGB[i];{RGB fits with naxis1=3, treated as 24 bits coded pixels in 2 dimensions}
- img_loaded[0,i,h]:=rgbdummy[0];{store in memory array}
- img_loaded[1,i,h]:=rgbdummy[1];{store in memory array}
- img_loaded[2,i,h]:=rgbdummy[2];{store in memory array}
- end;
- end
- else
- if nrbits=-32 then {floats}
- for i:=0 to width2-1 do {fill array}
- begin
- x_longword:=fitsbuffer4[i];
- img_loaded[0,i,h]:= x_single;{for conversion 32 bit float}
- end
- else
- if nrbits=32 then
- for i:=0 to width2-1 do {fill array} img_loaded[0,i,h]:= fitsbuffer4[i];
- end;
- inc(h)
- until ((h>=height2-1) or (count=0));
- result:=true;
- datamin_org:=0;
- datamax_org:=65535;
- except
- end;
-end;
-
-
procedure get_hist(colour:integer; img :image_array);
var
i,j,col,his_total,count, width5, height5,offsetW,offsetH : integer;
diff -Nru astap-cli-2022.03.01/unit_command_line_solving.pas astap-cli-2022.12.09/unit_command_line_solving.pas
--- astap-cli-2022.03.01/unit_command_line_solving.pas 2022-03-02 22:39:48.000000000 +0000
+++ astap-cli-2022.12.09/unit_command_line_solving.pas 2022-12-10 01:38:10.000000000 +0000
@@ -29,7 +29,7 @@
function SMedian(list: array of double; leng: integer): double;{get median of an array of double. Taken from CCDciel code but slightly modified}
-function solve_image(img :image_array; get_hist{update hist}:boolean) : boolean;{find match between image and star database}
+function solve_image(img :image_array ) : boolean;{find match between image and star database}
procedure bin_and_find_stars(img :image_array;binning:integer;cropping,hfd_min:double;get_hist{update hist}:boolean; out starlist3:star_list; out short_warning : string);{bin, measure background, find stars}
function report_binning(height:double) : integer;{select the binning}
var
@@ -172,53 +172,95 @@
end; {lsq_fit}
-procedure find_quads(starlist :star_list; min_leng:double; out quad_smallest:double; out quad_star_distances :star_list); {build quads using closest stars, revised 2020-9-28}
+
+procedure QuickSort_starlist(var A: star_list; iLo, iHi: Integer) ;{ Fast quick sort. Sorts elements in the array list with indices between lo and hi, sort in X only}
var
- i,j,k,nrstars_min_one,j_used1,j_used2,j_used3,nrquads,buffersize : integer;
- distance,distance1,distance2,distance3{,dummy },x1,x2,x3,x4,xt,y1,y2,y3,y4,yt,
- dist1,dist2,dist3,dist4,dist5,dist6,distx, dummy :double;
+ Lo, Hi : integer;
+ Pivot, Tx,Ty: double;{ pivot, T are the same type as the elements of array }
+begin
+ Lo := iLo;
+ Hi := iHi;
+ Pivot := A[0,(Lo + Hi) div 2];
+ repeat
+ while A[0,Lo] < Pivot do Inc(Lo) ; {sort in X only}
+ while A[0,Hi] > Pivot do Dec(Hi) ;
+ if Lo <= Hi then
+ begin {swap}
+ Tx := A[0,Lo];
+ Ty := A[1,Lo];
+ A[0,Lo] := A[0,Hi];
+ A[1,Lo] := A[1,Hi];
+ A[0,Hi] := Tx;
+ A[1,Hi] := Ty;
+ Inc(Lo) ;
+ Dec(Hi) ;
+ end;
+ until Lo > Hi;
+ if Hi > iLo then QuickSort_starlist(A, iLo, Hi) ; {executes itself recursively}
+ if Lo < iHi then QuickSort_starlist(A, Lo, iHi) ; {executes itself recursively}
+end;
+
+
+procedure find_quads(starlist :star_list; min_leng:double; out quad_smallest:double; out quad_star_distances :star_list); {build quads using closest stars, revised 2022-4-10}
+var
+ i,j,k,nrstars,j_used1,j_used2,j_used3,nrquads,Sstart,Send,tolerance : integer;
+ distance,distance1,distance2,distance3,x1,x2,x3,x4,xt,y1,y2,y3,y4,yt,
+ dist1,dist2,dist3,dist4,dist5,dist6,dummy,disty : double;
identical_quad : boolean;
begin
- nrstars_min_one:=Length(starlist[0])-1;
- buffersize:=nrstars_min_one;{number of quads will be lower}
+
+ nrstars:=Length(starlist[0]);{number of quads will be equal (super rare) or lower}
quad_smallest:=9999999;
- if nrstars_min_one<3 then
+
+ if nrstars<4 then
begin {not enough stars for quads}
SetLength(quad_star_distances,8,0);
exit;
end;
+ if nrstars>=150 then
+ begin
+ quickSort_starlist(starlist,0,nrstars-1); {sort in X only}
+ tolerance:=round(0.5*sqrt(nrstars));{tolerance band is about twice the every star distance}
+ end
+ else
+ tolerance:=1;{switch pre-filtering in X off}
+
nrquads:=0;
- SetLength(quad_star_distances,8,buffersize);{will contain the six distances and the central position}
+ SetLength(quad_star_distances,8,nrstars);{will contain the six distances and the central position}
j_used1:=0;{give it a default value}
j_used2:=0;
j_used3:=0;
- for i:=0 to nrstars_min_one do
+ for i:=0 to nrstars-1 do
begin
distance1:=1E99;{distance closest star}
distance2:=1E99;{distance second closest star}
distance3:=1E99;{distance third closest star}
- for j:=0 to nrstars_min_one do {find closest stars}
+
+ Sstart:=max(0,i-(nrstars div tolerance));
+ Send:=min(nrstars-1,i+(nrstars div tolerance)); {search in a limited X band only. The stars list is sorted in X. Search speed increases with about 30%}
+
+ for j:=Sstart to Send do {find closest stars}
begin
if j<>i{not the first star} then
begin
- distx:=sqr(starlist[0,j]-starlist[0,i]);
- if distx1 then {not an identical star. Mod 2021-6-25}
begin
if distancemin_leng then {large enough for earth based telescope}
begin
quad_star_distances[6,nrquads]:=xt;{store mean x position}
@@ -414,41 +452,7 @@
end;
result:=true;{3 or more references}
end;
-// else
-// if ((nr_references>=2) and (nrquads1<10) and (nrquads2<10)) then {use 8 stars of 2 quads as reference. Solver requires 3 equations minimum}
-// begin
-// {fill equations}
-// setlength(A_XYpositions,nr_references*4,3);{use 4 stars as reference instead of quads center}
-// setlength(b_Xrefpositions,nr_references*4);
-// setlength(b_Yrefpositions,nr_references*4);
-
-// for k:=0 to (nr_references*4)-1 do {use 2 * 4 stars of 2 quads}
-// begin
-// kd4:=k div 4;
-// A_XYpositions[kd4,0]:=starlistquads2[0,matchlist1[1,kd4]]; {x position of star}
-// A_XYpositions[kd4,1]:=starlistquads2[1,matchlist1[1,kd4]]; {y position of star}
-// A_XYpositions[kd4,2]:=1;
-// A_XYpositions[kd4+1,0]:=starlistquads2[2,matchlist1[1,kd4]]; {x position of star}
-// A_XYpositions[kd4+1,1]:=starlistquads2[3,matchlist1[1,kd4]]; {y position of star}
-// A_XYpositions[kd4+1,2]:=1;
-// A_XYpositions[kd4+2,0]:=starlistquads2[4,matchlist1[1,kd4]]; {x position of star}
-// A_XYpositions[kd4+2,1]:=starlistquads2[5,matchlist1[1,kd4]]; {y position of star}
-// A_XYpositions[kd4+2,2]:=1;
-// A_XYpositions[kd4+3,0]:=starlistquads2[6,matchlist1[1,kd4]]; {x position of star}
-// A_XYpositions[kd4+3,1]:=starlistquads2[7,matchlist1[1,kd4]]; {y position of star}
-// A_XYpositions[kd4+3,2]:=1;
-
-// b_Xrefpositions[kd4]:=starlistquads1[0,matchlist1[0,kd4]]; {x position of ref star}
-// b_Yrefpositions[kd4]:=starlistquads1[1,matchlist1[0,kd4]]; {Y position of ref star}
-// b_Xrefpositions[kd4+1]:=starlistquads1[2,matchlist1[0,kd4]]; {x position of ref star}
-// b_Yrefpositions[kd4+1]:=starlistquads1[3,matchlist1[0,kd4]]; {Y position of ref star}
-// b_Xrefpositions[kd4+2]:=starlistquads1[4,matchlist1[0,kd4]]; {x position of ref star}
-// b_Yrefpositions[kd4+2]:=starlistquads1[5,matchlist1[0,kd4]]; {Y position of ref star}
-// b_Xrefpositions[kd4+3]:=starlistquads1[6,matchlist1[0,kd4]]; {x position of ref star}
-// b_Yrefpositions[kd4+3]:=starlistquads1[7,matchlist1[0,kd4]]; {Y position of ref star}
-// result:=true;{8 star references from 2 quads}
-// end;
-// end;
+ // else It is possible to use one quad and the four star positions but it is not reliable.
matchlist2:=nil;
matchlist1:=nil;
end;
@@ -502,7 +506,7 @@
procedure find_stars(img :image_array;hfd_min:double;out starlist1: star_list);{find stars and put them in a list}
var
- fitsX, fitsY,nrstars,radius,i,j,max_stars,retries,m,n,xci,yci,sqr_radius : integer;
+ fitsX, fitsY,nrstars,radius,i,j,retries,m,n,xci,yci,sqr_radius : integer;
hfd1,star_fwhm,snr,xc,yc,highest_snr,flux, detection_level : double;
img_sa : image_array;
snr_list : array of double;
@@ -530,7 +534,7 @@
// solve_show_log:=stackmenu1.solve_show_log1.Checked;{show details, global variable}
// if solve_show_log then begin memo2_message('Start finding stars'); startTick2 := gettickcount64;end;
- max_stars:=strtoint(max_stars1);{maximum star to process, if so filter out brightest stars later}
+// max_stars:=strtoint(max_stars1);{maximum star to process, if so filter out brightest stars later}
if solve_show_log then begin memo2_message('Start finding stars'); startTick2 := gettickcount64;end;
SetLength(starlist1,2,buffersize);{set array length}
setlength(snr_list,buffersize);{set array length}
@@ -671,8 +675,6 @@
end;
-
-
function floattostrF2(const x:double; width1,decimals1 :word): string;
begin
str(x:width1:decimals1,result);
@@ -863,6 +865,61 @@
end;
+procedure check_pattern_filter(var img: image_array); {normalize bayer pattern. Colour shifts due to not using a white light source for the flat frames are avoided.}
+var
+ fitsX,fitsY,col,h,w,counter1,counter2, counter3,counter4 : integer;
+ value1,value2,value3,value4,maxval : double;
+ oddx, oddy :boolean;
+begin
+ col:=length(img);{the real number of colours}
+ h:=length(img[0,0]);{height}
+ w:=length(img[0]);{width}
+
+ if col>1 then
+ begin
+ memo2_message('Skipping check pattern filter. This filter works only for raw OSC images!');
+ exit;
+ end
+ else
+ memo2_message('Applying check pattern filter.');
+
+ value1:=0; value2:=0; value3:=0; value4:=0;
+ counter1:=0; counter2:=0; counter3:=0; counter4:=0;
+
+ for fitsY:=(h div 4) to (h*3) div 4 do {use one quarter of the image to find factors. Works also a little better if no dark-flat is subtracted. It also works better if boarder is black}
+ for fitsX:=(w div 4) to (w*3) div 4 do
+ begin
+ oddX:=odd(fitsX);
+ oddY:=odd(fitsY);
+ if ((oddX=false) and (oddY=false)) then begin value1:=value1+img[0,fitsX,fitsY]; inc(counter1) end else {separate counters for case odd() dimensions are used}
+ if ((oddX=true) and (oddY=false)) then begin value2:=value2+img[0,fitsX,fitsY]; inc(counter2) end else
+ if ((oddX=false) and (oddY=true)) then begin value3:=value3+img[0,fitsX,fitsY]; inc(counter3) end else
+ if ((oddX=true) and (oddY=true)) then begin value4:=value4+img[0,fitsX,fitsY]; inc(counter4) end;
+ end;
+
+ {now normalise the bayer pattern pixels}
+ value1:=value1/counter1;
+ value2:=value2/counter2;
+ value3:=value3/counter3;
+ value4:=value4/counter4;
+ maxval:=max(max(value1,value2),max(value3,value4));
+ value1:=maxval/value1;
+ value2:=maxval/value2;
+ value3:=maxval/value3;
+ value4:=maxval/value4;
+
+ for fitsY:=0 to h-1 do
+ for fitsX:=0 to w-1 do
+ begin
+ oddX:=odd(fitsX);
+ oddY:=odd(fitsY);
+ if ((value1<>1) and (oddX=false) and (oddY=false)) then img[0,fitsX,fitsY]:=round(img[0,fitsX,fitsY]*value1) else
+ if ((value2<>1) and (oddX=true) and (oddY=false)) then img[0,fitsX,fitsY]:=round(img[0,fitsX,fitsY]*value2) else
+ if ((value3<>1) and (oddX=false) and (oddY=true)) then img[0,fitsX,fitsY]:=round(img[0,fitsX,fitsY]*value3) else
+ if ((value4<>1) and (oddX=true) and (oddY=true)) then img[0,fitsX,fitsY]:=round(img[0,fitsX,fitsY]*value4);
+ end;
+end;
+
procedure binX1_crop(crop {0..1}:double; img : image_array; var img2: image_array);{crop image, make mono, no binning}
var fitsX,fitsY,k, w,h, shiftX,shiftY: integer;
@@ -1089,32 +1146,29 @@
end;
-function solve_image(img :image_array;get_hist{update hist}:boolean) : boolean;{find match between image and star database}
+function solve_image(img :image_array) : boolean;{find match between image and star database}
var
- nrstars,nrstars_required,count,max_distance,nr_quads, minimum_quads,database_stars,distance,binning,match_nr,
+ nrstars,nrstars_required,count,max_distance,nr_quads, minimum_quads,database_stars,binning,match_nr,
spiral_x, spiral_y, spiral_dx, spiral_dy,spiral_t : integer;
search_field,step_size,telescope_ra,telescope_dec,telescope_ra_offset,radius,fov2,fov_org, max_fov,fov_min,
- oversize,sep,seperation,ra7,dec7,centerX,centerY,cropping, min_star_size_arcsec,hfd_min,delta_ra,current_dist,
- quad_tolerance,dummy, extrastars,flip,extra : double;
- solution, go_ahead ,autoFOV : boolean;
+ oversize,sep_search,seperation,ra7,dec7,centerX,centerY,cropping, min_star_size_arcsec,hfd_min,delta_ra,current_dist,
+ quad_tolerance,dummy, extrastars,flip,extra,distance : double;
+ solution, go_ahead ,autoFOV,autoMaxstars : boolean;
startTick : qword;{for timing/speed purposes}
distancestr,oversize_mess,mess,suggest_str, warning_downsample, solved_in, offset_found,ra_offset,dec_offset,mount_info,mount_offset : string;
-
begin
result:=false;
esc_pressed:=false;
warning_str:='';{for header}
startTick := GetTickCount64;
-
quad_tolerance:=strtofloat2(quad_tolerance1);
- if ((fov_specified=false) and (cdelt2<>0)) then {no fov in native command line and cdelt2 in header}
+ if ((fov_specified=false) and (cdelt2<>0)) then {no FOV in native command line and cdelt2 in header}
fov_org:=min(180,height2*abs(cdelt2)) {calculate FOV. PI can give negative CDELT2}
else
fov_org:=min(180,strtofloat2(search_fov1));{use specfied FOV in stackmenu. 180 max to prevent runtime errors later}
-
if select_star_database(star_database1,fov_org)=false then {select database prior to cropping selection}
begin
result:=false;
@@ -1136,6 +1190,10 @@
if warning_str<>'' then memo2_message(warning_str);
end;
+ if check_pattern_filter1 then {for OSC images with low dimensions only}
+ check_pattern_filter(img);
+
+
if database_type=1476 then {.1476 database}
max_fov:=5.142857143 {warning FOV should be less the database tiles dimensions, so <=5.142857143 degrees. Otherwise a tile beyond next tile could be selected}
else {.1476 database}
@@ -1149,218 +1207,228 @@
min_star_size_arcsec:=strtofloat2(min_star_size1); {arc sec};
autoFOV:=(fov_org=0);{specified auto FOV}
- repeat {autoFOV loop}
- if autoFOV then
+ if max_stars=0 then
+ begin
+ autoMaxstars:=true;{try several values of max stars}
+ max_stars:=30; {will be doubled to 60 in the beginning}
+ end
+ else autoMaxstars:=false;
+
+ repeat {auto max star loop}
+ if autoMaxstars then
begin
- if fov_org=0 then
+ max_stars:=max_stars*2;{try with 60, 120, 240, 480 stars max}
+ memo2_message('Solving with '+inttostr(max_stars)+' stars maximum.');
+ end;
+ repeat {autoFOV loop}
+ if autoFOV then
begin
- if database_type<>001 then
- begin
- fov_org:=9.5;
- fov_min:=0.38;
- end
- else
+ if fov_org=0 then
begin
- fov_org:=90;
- fov_min:=12;
+ if database_type<>001 then
+ begin
+ fov_org:=9.5;
+ fov_min:=0.38;
+ end
+ else
+ begin
+ fov_org:=90;
+ fov_min:=12;
+ end
end
+ else fov_org:=fov_org/1.5;
+ memo2_message('Trying FOV: '+floattostrF(fov_org,ffFixed,0,1));
+ end;
+ if fov_org>max_fov then
+ begin
+ cropping:=max_fov/fov_org;
+ fov2:=max_fov; {temporary cropped image, adjust FOV to adapt}
end
- else fov_org:=fov_org/1.5;
- memo2_message('Trying FOV: '+floattostrF(fov_org,ffFixed,0,1));
- end;
- if fov_org>max_fov then
- begin
- cropping:=max_fov/fov_org;
- fov2:=max_fov; {temporary cropped image, adjust FOV to adapt}
- end
- else
- begin
- cropping:=1;
- fov2:=fov_org;
- end;
+ else
+ begin
+ cropping:=1;
+ fov2:=fov_org;
+ end;
- binning:=report_binning(height2*max_fov/max(fov_org{could be zero},max_fov)); {select binning on dimensions of cropped image only}
- hfd_min:=max(0.8,min_star_size_arcsec/(binning*fov_org*3600/height2) );{to ignore hot pixels which are too small}
- bin_and_find_stars(img,binning,cropping,hfd_min,get_hist{update hist}, starlist2, warning_downsample);{bin, measure background, find stars. Do this every repeat since hfd_min is adapted}
- nrstars:=Length(starlist2[0]);
-
- {prepare popupnotifier1 text}
- if force_oversize1=false then mess:=' normal' else mess:=' slow';
- memo2_message('ASTAP solver '+version+#10+
- 'Search radius: '+ radius_search1+' degrees, '+#10+
- 'Start position: '+prepare_ra(ra0,': ')+', '+prepare_dec(dec0,'d ')+#10+
- 'Image height: '+floattostrf2(fov_org,0,2)+' degrees'+#10+
- 'Binning: '+inttostr(binning)+'x'+inttostr(binning)+#10+
- 'Image dimensions:'+inttostr(width2)+'x'+inttostr(height2)+#10+
- 'Quad tolerance:'+quad_tolerance1+#10+
- 'Minimum star size:'+min_star_size1+'"' +#10+
- 'Speed:'+mess);
-
- nrstars_required:=round(nrstars*(height2/width2));{square search field based on height.}
-
- solution:=false; {assume no match is found}
- go_ahead:=(nrstars>=6); {bare minimum for three quads}
-
- if go_ahead then {enough stars, lets find quads}
- begin
- find_quads(starlist2,0 {min length}, quad_smallest,quad_star_distances2);{find star quads for new image. Quads and quad_smallest are binning independend}
- nr_quads:=Length(quad_star_distances2[0]);
- go_ahead:=nr_quads>=3; {enough quads?}
+ binning:=report_binning(height2*max_fov/max(fov_org{could be zero},max_fov)); {select binning on dimensions of cropped image only}
+ hfd_min:=max(0.8,min_star_size_arcsec/(binning*fov_org*3600/height2) );{to ignore hot pixels which are too small}
+ bin_and_find_stars(img,binning,cropping,hfd_min,true{update hist}, starlist2, warning_downsample);{bin, measure background, find stars. Do this every repeat since hfd_min is adapted}
+ nrstars:=Length(starlist2[0]);
+
+ {prepare popupnotifier1 text}
+ if force_oversize1=false then mess:=' normal' else mess:=' slow';
+ memo2_message('ASTAP solver version CLI-'+astap_version+#10+
+ 'Search radius: '+ radius_search1+' degrees, '+#10+
+ 'Start position: '+prepare_ra(ra0,': ')+', '+prepare_dec(dec0,'d ')+#10+
+ 'Image height: '+floattostrf2(fov_org,0,2)+' degrees'+#10+
+ 'Binning: '+inttostr(binning)+'x'+inttostr(binning)+#10+
+ 'Image dimensions: '+inttostr(width2)+'x'+inttostr(height2)+#10+
+ 'Quad tolerance: '+quad_tolerance1+#10+
+ 'Minimum star size: '+min_star_size1+'"' +#10+
+ 'Speed:'+mess);
- {The step size is fixed. If a low amount of quads are detected, the search window (so the database read area) is increased up to 200% guaranteeing that all quads of the image are compared with the database quads while stepping through the sky}
- if nr_quads<25 then oversize:=2 {make dimensions of square search window twice then the image height}
- else
- if nr_quads>100 then oversize:=1 {make dimensions of square search window equal to the image height}
- else
- oversize:=2*sqrt(25/nr_quads);{calculate between 25 th=2 and 100 th=1, quads are area related so take sqrt to get oversize}
+ nrstars_required:=round(nrstars*(height2/width2));{square search field based on height.}
- if force_oversize1 then
- begin
- oversize:=2;
- oversize_mess:='Search window at 200%'
- end
- else
- oversize_mess:='Search window at '+ inttostr(round((oversize)*100)) +'% based on the number of quads. Step size at 100% of image height';
+ solution:=false; {assume no match is found}
+ go_ahead:=(nrstars>=6); {bare minimum for three quads}
- radius:=strtofloat2(radius_search1);{radius search field}
+ if go_ahead then {enough stars, lets find quads}
+ begin
+ find_quads(starlist2,0 {min length}, quad_smallest,quad_star_distances2);{find star quads for new image. Quads and quad_smallest are binning independend}
+ nr_quads:=Length(quad_star_distances2[0]);
+ go_ahead:=nr_quads>=3; {enough quads?}
- memo2_message(inttostr(nrstars)+' stars, '+inttostr(nr_quads)+' quads selected in the image. '+inttostr(nrstars_required)+' database stars, '+inttostr(round(nr_quads*nrstars_required/nrstars))+' database quads required for the square search field of '+floattostrF2(fov2,0,1)+'d. '+oversize_mess );
+ {The step size is fixed. If a low amount of quads are detected, the search window (so the database read area) is increased up to 200% guaranteeing that all quads of the image are compared with the database quads while stepping through the sky}
+ if nr_quads<25 then oversize:=2 {make dimensions of square search window twice then the image height}
+ else
+ if nr_quads>100 then oversize:=1 {make dimensions of square search window equal to the image height}
+ else
+ oversize:=2*sqrt(25/nr_quads);{calculate between 25 th=2 and 100 th=1, quads are area related so take sqrt to get oversize}
- if nr_quads>500 then minimum_quads:=10 else {prevent false detections for star rich images}
- if nr_quads>200 then minimum_quads:=6 else {prevent false detections for star rich images}
- minimum_quads:=3; {3 quads giving 3 center quad references}
+ if force_oversize1 then
+ begin
+ oversize:=2;
+ oversize_mess:='Search window at 200%'
+ end
+ else
+ oversize_mess:='Search window at '+ inttostr(round((oversize)*100)) +'% based on the number of quads. Step size at 100% of image height';
- end
- else
- begin
- memo2_message('Only '+inttostr(nrstars)+' stars found in image. Abort');
- errorlevel:=2;
- end;
+ radius:=strtofloat2(radius_search1);{radius search field}
+ memo2_message(inttostr(nrstars)+' stars, '+inttostr(nr_quads)+' quads selected in the image. '+inttostr(nrstars_required)+' database stars, '+inttostr(round(nr_quads*nrstars_required/nrstars))+' database quads required for the square search field of '+floattostrF2(fov2,0,1)+'d. '+oversize_mess );
- if go_ahead then
- begin
- search_field:=fov2*(pi/180);
- STEP_SIZE:=search_field;{fixed step size search spiral. Prior to version 0.9.211 this was reduced for small star counts}
- if database_type=1 then
- begin {make smal steps for wide field images. Much more reliable}
- step_size:=step_size*0.25;
- max_distance:=round(radius/(0.25*fov2+0.00001)); {expressed in steps}
- memo2_message('Wide field, making small steps for reliable solving.');
+ minimum_quads:=3 + nr_quads div 100; {prevent false detections for star rich images, 3 quads give the 3 center quad references and is the bare minimum. It possible to use one quad and four star positions but it in not reliable}
end
else
- max_distance:=round(radius/(fov2+0.00001));{expressed in steps}
+ begin
+ memo2_message('Only '+inttostr(nrstars)+' stars found in image. Abort');
+ errorlevel:=2;
+ end;
- match_nr:=0;
- repeat {Maximum accuracy loop. In case math is found on a corner, do a second solve. Result will be more accurate using all stars of the image}
- count:=0;{search field counter}
- distance:=0; {required for reporting no too often}
- {spiral variables}
- spiral_x :=0;
- spiral_y :=0;
- spiral_dx := 0;{first step size x}
- spiral_dy := -1;{first step size y}
-
- repeat {search in squared spiral}
- {begin spiral routine, find a new squared spiral position position}
- if count<>0 then {first do nothing, start with [0 0] then start with [1 0],[1 1],[0 1],[-1 1],[-1 0],[-1 -1],[0 -1],[1 -1],[2 -1].[2 0] ..............}
- begin {start spiral around [0 0]}
- if ( (spiral_x = spiral_y) or ((spiral_x < 0) and (spiral_x = -spiral_y)) or ((spiral_x > 0) and (spiral_x = 1-spiral_y))) then {turning point}
- begin {swap dx by negative dy and dy by negative dx}
- spiral_t:=spiral_dx;
- spiral_dx := -spiral_dy;
- spiral_dy := spiral_t;
- end;
- spiral_x :=spiral_x+ spiral_dx;{walk through square}
- spiral_y :=spiral_y+ spiral_dy;{walk through square}
- end;{end spiral around [0 0]}
-
- {adapt search field to matrix position, +0+0/+1+0,+1+1,+0+1,-1+1,-1+0,-1-1,+0-1,+1-1..}
- telescope_dec:=STEP_SIZE*spiral_y+dec_radians;
- flip:=0;
- if telescope_dec>+pi/2 then begin telescope_dec:=pi-telescope_dec; flip:=pi; end {crossed the pole}
- else
- if telescope_dec<-pi/2 then begin telescope_dec:=-pi-telescope_dec; flip:=pi; end;
+ if go_ahead then
+ begin
+ search_field:=fov2*(pi/180);
+ STEP_SIZE:=search_field;{fixed step size search spiral. Prior to version 0.9.211 this was reduced for small star counts}
+ if database_type=1 then
+ begin {make smal steps for wide field images. Much more reliable}
+ step_size:=step_size*0.25;
+ max_distance:=round(radius/(0.25*fov2+0.00001)); {expressed in steps}
+ memo2_message('Wide field, making small steps for reliable solving.');
+ end
+ else
+ max_distance:=round(radius/(fov2+0.00001));{expressed in steps}
- if telescope_dec>0 then extra:=step_size/2 else extra:=-step_size/2;{use the distance furthest away from the pole}
+ match_nr:=0;
+ repeat {Maximum accuracy loop. In case math is found on a corner, do a second solve. Result will be more accurate using all stars of the image}
+ count:=0;{search field counter}
+ distance:=0; {required for reporting no too often}
+ {spiral variables}
+ spiral_x :=0;
+ spiral_y :=0;
+ spiral_dx := 0;{first step size x}
+ spiral_dy := -1;{first step size y}
+
+ repeat {search in squared spiral}
+ {begin spiral routine, find a new squared spiral position position}
+ if count<>0 then {first do nothing, start with [0 0] then start with [1 0],[1 1],[0 1],[-1 1],[-1 0],[-1 -1],[0 -1],[1 -1],[2 -1].[2 0] ..............}
+ begin {start spiral around [0 0]}
+ if ( (spiral_x = spiral_y) or ((spiral_x < 0) and (spiral_x = -spiral_y)) or ((spiral_x > 0) and (spiral_x = 1-spiral_y))) then {turning point}
+ begin {swap dx by negative dy and dy by negative dx}
+ spiral_t:=spiral_dx;
+ spiral_dx := -spiral_dy;
+ spiral_dy := spiral_t;
+ end;
+ spiral_x :=spiral_x+ spiral_dx;{walk through square}
+ spiral_y :=spiral_y+ spiral_dy;{walk through square}
+ end;{end spiral around [0 0]}
+
+ {adapt search field to matrix position, +0+0/+1+0,+1+1,+0+1,-1+1,-1+0,-1-1,+0-1,+1-1..}
+ telescope_dec:=STEP_SIZE*spiral_y+dec_radians;
+ flip:=0;
+ if telescope_dec>+pi/2 then begin telescope_dec:=pi-telescope_dec; flip:=pi; end {crossed the pole}
+ else
+ if telescope_dec<-pi/2 then begin telescope_dec:=-pi-telescope_dec; flip:=pi; end;
- telescope_ra_offset:= (STEP_SIZE*spiral_x/cos(telescope_dec-extra));{step larger near pole. This telescope_ra is an offsett from zero}
- if ((telescope_ra_offset<=+pi/2+step_size/2) and (telescope_ra_offset>=-pi/2)) then {step_size for overlap}
- begin
- telescope_ra:=fnmodulo(flip+ra_radians+telescope_ra_offset,2*pi);{add offset to ra after the if statement! Otherwise no symmetrical search}
+ if telescope_dec>0 then extra:=step_size/2 else extra:=-step_size/2;{use the distance furthest away from the pole}
- ang_sep(telescope_ra,telescope_dec,ra_radians,dec_radians, {out}seperation);{calculates angular separation. according formula 9.1 old Meeus or 16.1 new Meeus, version 2018-5-23}
- if seperation<=radius*pi/180+step_size/2 then {Use only the circular area withing the square area}
+ telescope_ra_offset:= (STEP_SIZE*spiral_x/cos(telescope_dec-extra));{step larger near pole. This telescope_ra is an offsett from zero}
+ if ((telescope_ra_offset<=+pi/2+step_size/2) and (telescope_ra_offset>=-pi/2)) then {step_size for overlap}
begin
+ telescope_ra:=fnmodulo(flip+ra_radians+telescope_ra_offset,2*pi);{add offset to ra after the if statement! Otherwise no symmetrical search}
- {info reporting}
- //stackmenu1.field1.caption:= '['+inttostr(spiral_x)+','+inttostr(spiral_y)+']';{show on stackmenu what's happening}
- if ((spiral_x>distance) or (spiral_y>distance)) then {new distance reached. Update once in the square spiral, so not too often since it cost CPU time}
+ ang_sep(telescope_ra,telescope_dec,ra_radians,dec_radians, {out}seperation);{calculates angular separation. according formula 9.1 old Meeus or 16.1 new Meeus, version 2018-5-23}
+ if seperation<=radius*pi/180+step_size/2 then {Use only the circular area withing the square area}
begin
- distance:=max(spiral_x,spiral_y);{update status}
- distancestr:=inttostr(round(seperation*180/pi))+'d,';{show on stackmenu what's happening}
- write(distancestr);
- end; {info reporting}
-
- {If a low amount of quads are detected, the search window (so the database read area) is increased up to 200% guaranteeing that all quads of the image are compared with the database quads while stepping through the sky}
- {read nrstars_required stars from database. If search field is oversized, number of required stars increases with the power of the oversize factor. So the star density will be the same as in the image to solve}
- extrastars:=1/1.1;{star with a factor of one}
- repeat {loop to add extra stars if too many too small quads are excluding. Note the database is made by a space telescope with a resolution exceeding all earth telescopes}
- extrastars:=extrastars*1.1;
- if read_stars(telescope_ra,telescope_dec,search_field*oversize,round(nrstars_required*oversize*oversize*extrastars) ,{var}database_stars)= false then
+ {info reporting}
+ //stackmenu1.field1.caption:= '['+inttostr(spiral_x)+','+inttostr(spiral_y)+']';{show on stackmenu what's happening}
+ if seperation*180/pi>distance+fov_org then {new distance reached. Update once in the square spiral, so not too often since it cost CPU time}
begin
- memo2_message('Error, no star database found at '+database_path+' ! Download the h18 (or h17, v17) and install.');
- errorlevel:=33;{read error star database}
- exit; {no stars}
- end;
- find_quads(starlist1,quad_smallest*(fov_org*3600/height2 {pixelsize in"})*0.99 {filter value to exclude too small quads, convert pixels to arcsec as in database}, dummy,quad_star_distances1);{find quads for reference image/database. Filter out too small quads for Earth based telescopes}
- {Note quad_smallest is binning independent value. Don't use cdelt2 for pixelsize calculation since fov_specified could be true making cdelt2 unreliable or fov=auto}
- until ((nrstars_required>database_stars) {No more stars available in the database}
- or (nr_quads<1.1*Length(quad_star_distances1[0])*nrstars/nrstars_required) {Enough quads found. The amount quads could be too low because due to filtering out too small database quads (center m13, M16)in routine find_quads}
- or (extrastars>15)) {Go up this factor maximum};
-
- if ((solve_show_log) and (extrastars>1)) then memo2_message('Too many small quads excluded due to higher resolution database, increased the number of stars with '+inttostr(round((extrastars-1)*100))+'%');
-
- if solve_show_log then {global variable set in find stars}
- memo2_message('Search '+ inttostr(count)+', ['+inttostr(spiral_x)+','+inttostr(spiral_y)+'], position: '+ prepare_ra(telescope_ra,': ')+prepare_dec(telescope_dec,'d ')+#9+' Down to magn '+ floattostrF2(mag2/10,0,1) +#9+' '+inttostr(database_stars)+' database stars' +#9+' '+inttostr(length(quad_star_distances1[0]))+' database quads to compare.');
-
- // for testing purposes
- // create supplement lines for sky coverage testing and write to log using -log
- // memo2.add(floattostr(telescope_ra*12/pi)+',,,'+floattostr(telescope_dec*180/pi)+',,,,'+inttostr(count)+',,-99'); {create hnsky supplement to test sky coverage}
-
- solution:=find_offset_and_rotation(minimum_quads {>=3},quad_tolerance);{find an solution}
- end; {within search circle. Otherwise the search is within a kind of square}
- end;{ra in range}
- inc(count);{step further in spiral}
- until ((solution) or (spiral_x>max_distance));{squared spiral search}
+ distance:=seperation*180/pi;
+ distancestr:=inttostr(round(seperation*180/pi))+'d,';{show on stackmenu what's happening}
+ write(distancestr);
+ end; {info reporting}
+
+ {If a low amount of quads are detected, the search window (so the database read area) is increased up to 200% guaranteeing that all quads of the image are compared with the database quads while stepping through the sky}
+ {read nrstars_required stars from database. If search field is oversized, number of required stars increases with the power of the oversize factor. So the star density will be the same as in the image to solve}
+ extrastars:=1/1.1;{star with a factor of one}
+ repeat {loop to add extra stars if too many too small quads are excluding. Note the database is made by a space telescope with a resolution exceeding all earth telescopes}
+ extrastars:=extrastars*1.1;
+ if read_stars(telescope_ra,telescope_dec,search_field*oversize,round(nrstars_required*oversize*oversize*extrastars) ,{var}database_stars)= false then
+ begin
+ memo2_message('Error, no star database found at '+database_path+' ! Download the h18 (or h17, v17) and install.');
+ errorlevel:=33;{read error star database}
+ exit; {no stars}
+ end;
+ find_quads(starlist1,quad_smallest*(fov_org*3600/height2 {pixelsize in"})*0.99 {filter value to exclude too small quads, convert pixels to arcsec as in database}, dummy,quad_star_distances1);{find quads for reference image/database. Filter out too small quads for Earth based telescopes}
+ {Note quad_smallest is binning independent value. Don't use cdelt2 for pixelsize calculation since fov_specified could be true making cdelt2 unreliable or fov=auto}
+ until ((nrstars_required>database_stars) {No more stars available in the database}
+ or (nr_quads<1.1*Length(quad_star_distances1[0])*nrstars/nrstars_required) {Enough quads found. The amount quads could be too low because due to filtering out too small database quads (center m13, M16)in routine find_quads}
+ or (extrastars>15)) {Go up this factor maximum};
+
+ if ((solve_show_log) and (extrastars>1)) then memo2_message('Too many small quads excluded due to higher resolution database, increased the number of stars with '+inttostr(round((extrastars-1)*100))+'%');
+
+ if solve_show_log then {global variable set in find stars}
+ memo2_message('Search '+ inttostr(count)+', ['+inttostr(spiral_x)+','+inttostr(spiral_y)+'], position: '+ prepare_ra(telescope_ra,': ')+prepare_dec(telescope_dec,'d ')+#9+' Down to magn '+ floattostrF2(mag2/10,0,1) +#9+' '+inttostr(database_stars)+' database stars' +#9+' '+inttostr(length(quad_star_distances1[0]))+' database quads to compare.');
+
+ // for testing purposes
+ // create supplement lines for sky coverage testing and write to log using -log
+ // memo2.add(floattostr(telescope_ra*12/pi)+',,,'+floattostr(telescope_dec*180/pi)+',,,,'+inttostr(count)+',,-99'); {create hnsky supplement to test sky coverage}
+
+ solution:=find_offset_and_rotation(minimum_quads {>=3},quad_tolerance);{find an solution}
+ end; {within search circle. Otherwise the search is within a kind of square}
+ end;{ra in range}
+ inc(count);{step further in spiral}
+ until ((solution) or (spiral_x>max_distance));{squared spiral search}
- if solution then
- begin
- centerX:=(width2-1)/2 ;{center image in 0..width2-1 range}
- centerY:=(height2-1)/2;{center image in 0..height2-1 range}
- crpix1:=centerX+1;{center image in fits coordinate range 1..width2}
- crpix2:=centery+1;
-
- standard_equatorial( telescope_ra,telescope_dec,
- (solution_vectorX[0]*(centerX) + solution_vectorX[1]*(centerY) +solution_vectorX[2]), {x}
- (solution_vectorY[0]*(centerX) + solution_vectorY[1]*(centerY) +solution_vectorY[2]), {y}
- 1, {CCD scale}
- ra0 ,dec0{center equatorial position});
- if match_nr=0 then ang_sep(ra_radians,dec_radians,ra0,dec0, sep);{offset found}
- ra_radians:=ra0;
- dec_radians:=dec0;
- current_dist:=sqrt(sqr(solution_vectorX[0]*(centerX) + solution_vectorX[1]*(centerY) +solution_vectorX[2]) + sqr(solution_vectorY[0]*(centerX) + solution_vectorY[1]*(centerY) +solution_vectorY[2]))/3600; {current distance telescope and image center in degrees}
- inc(match_nr);
- end;
- until ((solution=false) or (current_dist=2));{Maximum accurcy loop. After match possible on a corner do a second solve using the found ra0,dec0 for maximum accuracy USING ALL STARS}
+ if solution then
+ begin
+ centerX:=(width2-1)/2 ;{center image in 0..width2-1 range}
+ centerY:=(height2-1)/2;{center image in 0..height2-1 range}
+ crpix1:=centerX+1;{center image in fits coordinate range 1..width2}
+ crpix2:=centery+1;
+
+ standard_equatorial( telescope_ra,telescope_dec,
+ (solution_vectorX[0]*(centerX) + solution_vectorX[1]*(centerY) +solution_vectorX[2]), {x}
+ (solution_vectorY[0]*(centerX) + solution_vectorY[1]*(centerY) +solution_vectorY[2]), {y}
+ 1, {CCD scale}
+ ra_radians ,dec_radians {center equatorial position});
+ current_dist:=sqrt(sqr(solution_vectorX[0]*(centerX) + solution_vectorX[1]*(centerY) +solution_vectorX[2]) + sqr(solution_vectorY[0]*(centerX) + solution_vectorY[1]*(centerY) +solution_vectorY[2]))/3600; {current distance telescope and image center in degrees}
+ inc(match_nr);
+ end;
+ until ((solution=false) or (current_dist=2));{Maximum accurcy loop. After match possible on a corner do a second solve using the found ra0,dec0 for maximum accuracy USING ALL STARS}
- end; {enough quads in image}
- until ((autoFOV=false) or (solution) or (fov2<=fov_min)); {loop for autoFOV from 9.5 to 0.37 degrees. Will lock between 9.5*1.25 downto 0.37/1.25 or 11.9 downto 0.3 degrees}
+ end; {enough quads in image}
+ until ((autoFOV=false) or (solution) or (fov2<=fov_min)); {loop for autoFOV from 9.5 to 0.37 degrees. Will lock between 9.5*1.25 downto 0.37/1.25 or 11.9 downto 0.3 degrees}
+ until ((autoMaxstars=false) or (solution) or (max_stars>=480) or (max_stars-5>nrstars){no more stars to find});{auto max star loop}
if solution then
begin
+ ang_sep(ra_radians,dec_radians,ra0,dec0, sep_search);{calculate search offset}
+ ra0:=ra_radians;//store solution
+ dec0:=dec_radians;
+
memo2_message(#10+inttostr(nr_references)+ ' of '+ inttostr(nr_references2)+' quads selected matching within '+quad_tolerance1+' tolerance.' {2 quads are required giving 8 star references or 3 quads giving 3 center quad references}
+#10+'Solution["] x:='+floattostr6(solution_vectorX[0])+'*x+ '+floattostr6(solution_vectorX[1])+'*y+ '+floattostr6(solution_vectorX[2])
+', y:='+floattostr6(solution_vectorY[0])+'*x+ '+floattostr6(solution_vectorY[1])+'*y+ '+floattostr6(solution_vectorY[2]) );
@@ -1400,12 +1468,11 @@
new_to_old_WCS;
solved_in:='Solved in '+ floattostr(round((GetTickCount64 - startTick)/100)/10)+' sec';{make string to report in FITS header.}
- offset_found:=distance_to_string(sep ,sep)+'.';
-
+ offset_found:=distance_to_string(sep_search ,sep_search)+'.';
if ra_mount<99 then {mount position known and specified}
begin
- ra_offset:=distance_to_string(sep, fnmodulo(ra_mount-ra0,pi)*cos((dec0+dec_mount)*0.5 {average dec}));
- dec_offset:=distance_to_string(sep,dec_mount-dec0);
+ ra_offset:=distance_to_string(sep_search, pi*frac((ra_mount-ra0)/pi) * cos((dec0+dec_mount)*0.5 {average dec}));
+ dec_offset:=distance_to_string(sep_search,dec_mount-dec0);
mount_offset:=' Mount offset RA='+ra_offset+', DEC='+dec_offset;{ascii}
mount_info:=' Mount Δα='+ra_offset+ ', Δδ='+dec_offset+'. ';
@@ -1439,7 +1506,7 @@
update_float ('CD1_2 =',' / CD matrix to convert (x,y) to (Ra, Dec) ' ,cd1_2);
update_float ('CD2_1 =',' / CD matrix to convert (x,y) to (Ra, Dec) ' ,cd2_1);
update_float ('CD2_2 =',' / CD matrix to convert (x,y) to (Ra, Dec) ' ,cd2_2);
- update_text ('PLTSOLVD=',' T / ASTAP internal solver ');
+ update_text ('PLTSOLVD=',' T / Astrometric solved by ASTAP_CLI v'+astap_version+'. ');
update_text ('COMMENT 7', solved_in+' Offset was '+offset_found+mount_offset);
diff -Nru astap-cli-2022.03.01/unit_command_line_star_database.pas astap-cli-2022.12.09/unit_command_line_star_database.pas
--- astap-cli-2022.03.01/unit_command_line_star_database.pas 2022-03-02 22:39:48.000000000 +0000
+++ astap-cli-2022.12.09/unit_command_line_star_database.pas 2022-12-10 01:38:10.000000000 +0000
@@ -2760,11 +2760,11 @@
if fov>20 then
begin
- if fileexists( database_path+'w08_0101.001') then begin name_database:='w08';database_type:=001; end
+ if fileexists( database_path+'w08_0101.001') then begin name_database:='w08';database_type:=001; exit; end
else
memo2_message('Could not find w08 star database. Will try with an other database.');
- end
- else
+ end;
+
if fileexists( database_path+'h18_0101.1476') then begin name_database:='h18'; end
else
if fileexists( database_path+'g18_0101.290') then begin name_database:='g18'; database_type:=290; end