Open Build Service

atmos.spec Changed

atmos-2.97.1.tar.gz/configure -> atmos-2.97.2.tar.gz/configure Changed

@@ -1,6 +1,6 @@
 #! /bin/sh
 # Guess values for system-dependent variables and create Makefiles.
-# Generated by GNU Autoconf 2.63 for atmos 2.97.1.
+# Generated by GNU Autoconf 2.63 for atmos 2.97.2.
 #
 # Report bugs to <mass-general@curl.sai.msu.ru>.
 #
@@ -596,8 +596,8 @@
 # Identity of this package.
 PACKAGE_NAME='atmos'
 PACKAGE_TARNAME='atmos'
-PACKAGE_VERSION='2.97.1'
-PACKAGE_STRING='atmos 2.97.1'
+PACKAGE_VERSION='2.97.2'
+PACKAGE_STRING='atmos 2.97.2'
 PACKAGE_BUGREPORT='mass-general@curl.sai.msu.ru'
 
 ac_unique_file="src/atmos.h"
@@ -1293,7 +1293,7 @@
   # Omit some internal or obsolete options to make the list less imposing.
   # This message is too long to be a string in the A/UX 3.1 sh.
   cat <<_ACEOF
-\`configure' configures atmos 2.97.1 to adapt to many kinds of systems.
+\`configure' configures atmos 2.97.2 to adapt to many kinds of systems.
 
 Usage: $0 OPTION... VAR=VALUE...
 
@@ -1359,7 +1359,7 @@
 
 if test -n "$ac_init_help"; then
   case $ac_init_help in
-     short | recursive ) echo "Configuration of atmos 2.97.1:";;
+     short | recursive ) echo "Configuration of atmos 2.97.2:";;
    esac
   cat <<\_ACEOF
 
@@ -1449,7 +1449,7 @@
 test -n "$ac_init_help" && exit $ac_status
 if $ac_init_version; then
   cat <<\_ACEOF
-atmos configure 2.97.1
+atmos configure 2.97.2
 generated by GNU Autoconf 2.63
 
 Copyright (C) 1992, 1993, 1994, 1995, 1996, 1998, 1999, 2000, 2001,
@@ -1463,7 +1463,7 @@
 This file contains any messages produced by compilers while
 running configure, to aid debugging if configure makes a mistake.
 
-It was created by atmos $as_me 2.97.1, which was
+It was created by atmos $as_me 2.97.2, which was
 generated by GNU Autoconf 2.63.  Invocation command line was
 
   $ $0 $@
@@ -2312,7 +2312,7 @@
 
 # Define the identity of the package.
  PACKAGE='atmos'
- VERSION='2.97.1'
+ VERSION='2.97.2'
 
 
 cat >>confdefs.h <<_ACEOF
@@ -4159,7 +4159,6 @@
 
 
 
-
 for ac_header in boost/algorithm/string.hpp \
 	boost/bind.hpp \
 	boost/cstdint.hpp \
@@ -4171,7 +4170,6 @@
 	boost/io/ios_state.hpp \
 	boost/iterator/transform_iterator.hpp \
 	boost/lexical_cast.hpp \
-	boost/math/tools/roots.hpp \
 	boost/mem_fn.hpp \
 	boost/mpl/if.hpp \
 	boost/program_options.hpp \
@@ -4333,6 +4331,160 @@
 
 
 
+for ac_header in boost/math/tools/roots.hpp \
+	         tr1/tuple
+do
+as_ac_Header=`$as_echo "ac_cv_header_$ac_header" | $as_tr_sh`
+if { as_var=$as_ac_Header; eval "test \"\${$as_var+set}\" = set"; }; then
+  { $as_echo "$as_me:$LINENO: checking for $ac_header" >&5
+$as_echo_n "checking for $ac_header... " >&6; }
+if { as_var=$as_ac_Header; eval "test \"\${$as_var+set}\" = set"; }; then
+  $as_echo_n "(cached) " >&6
+fi
+ac_res=`eval 'as_val=${'$as_ac_Header'}
+		 $as_echo "$as_val"'`
+	       { $as_echo "$as_me:$LINENO: result: $ac_res" >&5
+$as_echo "$ac_res" >&6; }
+else
+  # Is the header compilable?
+{ $as_echo "$as_me:$LINENO: checking $ac_header usability" >&5
+$as_echo_n "checking $ac_header usability... " >&6; }
+cat >conftest.$ac_ext <<_ACEOF
+/* confdefs.h.  */
+_ACEOF
+cat confdefs.h >>conftest.$ac_ext
+cat >>conftest.$ac_ext <<_ACEOF
+/* end confdefs.h.  */
+$ac_includes_default
+#include <$ac_header>
+_ACEOF
+rm -f conftest.$ac_objext
+if { (ac_try="$ac_compile"
+case "(($ac_try" in
+  *\"* | *\`* | *\\*) ac_try_echo=\$ac_try;;
+  *) ac_try_echo=$ac_try;;
+esac
+eval ac_try_echo="\"\$as_me:$LINENO: $ac_try_echo\""
+$as_echo "$ac_try_echo") >&5
+  (eval "$ac_compile") 2>conftest.er1
+  ac_status=$?
+  grep -v '^ *+' conftest.er1 >conftest.err
+  rm -f conftest.er1
+  cat conftest.err >&5
+  $as_echo "$as_me:$LINENO: \$? = $ac_status" >&5
+  (exit $ac_status); } && {
+	 test -z "$ac_cxx_werror_flag" ||
+	 test ! -s conftest.err
+       } && test -s conftest.$ac_objext; then
+  ac_header_compiler=yes
+else
+  $as_echo "$as_me: failed program was:" >&5
+sed 's/^/| /' conftest.$ac_ext >&5
+
+	ac_header_compiler=no
+fi
+
+rm -f core conftest.err conftest.$ac_objext conftest.$ac_ext
+{ $as_echo "$as_me:$LINENO: result: $ac_header_compiler" >&5
+$as_echo "$ac_header_compiler" >&6; }
+
+# Is the header present?
+{ $as_echo "$as_me:$LINENO: checking $ac_header presence" >&5
+$as_echo_n "checking $ac_header presence... " >&6; }
+cat >conftest.$ac_ext <<_ACEOF
+/* confdefs.h.  */
+_ACEOF
+cat confdefs.h >>conftest.$ac_ext
+cat >>conftest.$ac_ext <<_ACEOF
+/* end confdefs.h.  */
+#include <$ac_header>
+_ACEOF
+if { (ac_try="$ac_cpp conftest.$ac_ext"
+case "(($ac_try" in
+  *\"* | *\`* | *\\*) ac_try_echo=\$ac_try;;
+  *) ac_try_echo=$ac_try;;
+esac
+eval ac_try_echo="\"\$as_me:$LINENO: $ac_try_echo\""
+$as_echo "$ac_try_echo") >&5
+  (eval "$ac_cpp conftest.$ac_ext") 2>conftest.er1
+  ac_status=$?
+  grep -v '^ *+' conftest.er1 >conftest.err
+  rm -f conftest.er1
+  cat conftest.err >&5
+  $as_echo "$as_me:$LINENO: \$? = $ac_status" >&5
+  (exit $ac_status); } >/dev/null && {
+	 test -z "$ac_cxx_preproc_warn_flag$ac_cxx_werror_flag" ||
+	 test ! -s conftest.err
+       }; then
+  ac_header_preproc=yes
+else
+  $as_echo "$as_me: failed program was:" >&5
+sed 's/^/| /' conftest.$ac_ext >&5
+
+  ac_header_preproc=no
+fi
+
+rm -f conftest.err conftest.$ac_ext
+{ $as_echo "$as_me:$LINENO: result: $ac_header_preproc" >&5
+$as_echo "$ac_header_preproc" >&6; }
+
+# So?  What about this header?
+case $ac_header_compiler:$ac_header_preproc:$ac_cxx_preproc_warn_flag in
+  yes:no: )
+    { $as_echo "$as_me:$LINENO: WARNING: $ac_header: accepted by the compiler, rejected by the preprocessor!" >&5
+$as_echo "$as_me: WARNING: $ac_header: accepted by the compiler, rejected by the preprocessor!" >&2;}
+    { $as_echo "$as_me:$LINENO: WARNING: $ac_header: proceeding with the compiler's result" >&5
+$as_echo "$as_me: WARNING: $ac_header: proceeding with the compiler's result" >&2;}
+    ac_header_preproc=yes
+    ;;
+  no:yes:* )
+    { $as_echo "$as_me:$LINENO: WARNING: $ac_header: present but cannot be compiled" >&5
+$as_echo "$as_me: WARNING: $ac_header: present but cannot be compiled" >&2;}
+    { $as_echo "$as_me:$LINENO: WARNING: $ac_header:     check for missing prerequisite headers?" >&5
+$as_echo "$as_me: WARNING: $ac_header:     check for missing prerequisite headers?" >&2;}
+    { $as_echo "$as_me:$LINENO: WARNING: $ac_header: see the Autoconf documentation" >&5
+$as_echo "$as_me: WARNING: $ac_header: see the Autoconf documentation" >&2;}
+    { $as_echo "$as_me:$LINENO: WARNING: $ac_header:     section \"Present But Cannot Be Compiled\"" >&5
+$as_echo "$as_me: WARNING: $ac_header:     section \"Present But Cannot Be Compiled\"" >&2;}
+    { $as_echo "$as_me:$LINENO: WARNING: $ac_header: proceeding with the preprocessor's result" >&5
+$as_echo "$as_me: WARNING: $ac_header: proceeding with the preprocessor's result" >&2;}
+    { $as_echo "$as_me:$LINENO: WARNING: $ac_header: in the future, the compiler will take precedence" >&5
+$as_echo "$as_me: WARNING: $ac_header: in the future, the compiler will take precedence" >&2;}
+    ( cat <<\_ASBOX
+## ------------------------------------------- ##
+## Report this to mass-general@curl.sai.msu.ru ##
+## ------------------------------------------- ##
+_ASBOX
+     ) | sed "s/^/$as_me: WARNING:     /" >&2
+    ;;
+esac
+{ $as_echo "$as_me:$LINENO: checking for $ac_header" >&5
+$as_echo_n "checking for $ac_header... " >&6; }
+if { as_var=$as_ac_Header; eval "test \"\${$as_var+set}\" = set"; }; then
+  $as_echo_n "(cached) " >&6
+else
+  eval "$as_ac_Header=\$ac_header_preproc"
+fi
+ac_res=`eval 'as_val=${'$as_ac_Header'}
+		 $as_echo "$as_val"'`
+	       { $as_echo "$as_me:$LINENO: result: $ac_res" >&5
+$as_echo "$ac_res" >&6; }
+
+fi
+as_val=`eval 'as_val=${'$as_ac_Header'}
+		 $as_echo "$as_val"'`
+   if test "x$as_val" = x""yes; then
+  cat >>confdefs.h <<_ACEOF
+#define `$as_echo "HAVE_$ac_header" | $as_tr_cpp` 1
+_ACEOF
+
+fi
+
+done
+
+
+
+
 for ac_header in gsl/gsl_errno.h \
 	gsl/gsl_integration.h
 do
@@ -5303,7 +5455,7 @@
 # report actual input values of CONFIG_FILES etc. instead of their
 # values after options handling.
 ac_log="
-This file was extended by atmos $as_me 2.97.1, which was
+This file was extended by atmos $as_me 2.97.2, which was
 generated by GNU Autoconf 2.63.  Invocation command line was
 
   CONFIG_FILES    = $CONFIG_FILES
@@ -5366,7 +5518,7 @@
 _ACEOF
 cat >>$CONFIG_STATUS <<_ACEOF || ac_write_fail=1
 ac_cs_version="\\
-atmos config.status 2.97.1
+atmos config.status 2.97.2
 configured by $0, generated by GNU Autoconf 2.63,
   with options \\"`$as_echo "$ac_configure_args" | sed 's/^ //; s/\\""\`\$/\\\\&/g'`\\"

atmos-2.97.1.tar.gz/configure.in -> atmos-2.97.2.tar.gz/configure.in Changed

atmos-2.97.1.tar.gz/src/atmos.h -> atmos-2.97.2.tar.gz/src/atmos.h Changed

atmos-2.97.1.tar.gz/src/background_factory.cpp -> atmos-2.97.2.tar.gz/src/background_factory.cpp Changed

atmos-2.97.1.tar.gz/src/config.h.in -> atmos-2.97.2.tar.gz/src/config.h.in Changed

atmos-2.97.1.tar.gz/src/fft.h -> atmos-2.97.2.tar.gz/src/fft.h Changed

@@ -1,5 +1,5 @@
 /*
-    $Id: fft.h 180 2011-01-30 15:56:33Z matwey $
+    $Id: fft.h 199 2011-02-03 09:38:16Z matwey $
     Copyright (C) 2009  Sternberg Astronomical Institute, MSU
 
     This program is free software: you can redistribute it and/or modify
@@ -54,9 +54,7 @@
 
 		if( size < 1 ) return;
 
-		size_type n = nearest_log2( size );
-
-		assert( std::pow( size_type(2), n ) == size );
+		assert( std::pow( size_type(2), nearest_log2( size ) ) == size );
 
 		size_type max_step = size / 2;
 		size_type cur_step;
@@ -83,9 +81,7 @@
 
 		if( size < 1 ) return;
 
-		size_type n = nearest_log2( size );
-
-		assert( std::pow( size_type(2), n ) == size );
+		assert( std::pow( size_type(2), nearest_log2( size ) ) == size );
 
 		size_type max_step = size / 2;
 		size_type cur_step;
@@ -111,7 +107,6 @@
 			std::complex< typename expression_type::value_type > >::type complex_type;
 		typedef typename complex_type::value_type real_type;
 		typedef complex_type value_type;
-		typedef void const_iterator;
 	private:
 		template<class V> void copy_helper( V& v, const vector_expression<E>& e ) const {
 			copy_bitreversal( v, e );
@@ -163,7 +158,7 @@
 		e().apply( v );
 	}
 
-};
+}
 
 template<class E> detail::fft_radix2_result<E> fft_radix2( const vector_expression<E>& e, bool sign = true ) {
 	return detail::fft_radix2_result<E>(e,sign);

atmos-2.97.1.tar.gz/src/file_output.cpp -> atmos-2.97.2.tar.gz/src/file_output.cpp Changed

@@ -1,5 +1,5 @@
 /*
-    $Id: file_output.cpp 180 2011-01-30 15:56:33Z matwey $
+    $Id: file_output.cpp 200 2011-02-03 14:45:13Z matwey $
     Copyright (C) 2009  Sternberg Astronomical Institute, MSU
 
     This program is free software: you can redistribute it and/or modify
@@ -45,7 +45,7 @@
 	boost::io::ios_all_saver asaver( m_stm );
 	m_stm << "M " << time << " Background measurements" << std::endl;
 }
-void file::atmosphere( const time_type& time, double free_seeing, double free_seeing_err, double dimm_seeing, double dimm_seeing_err, double dimm_turbulence, double dimm_turbulence_err, double isoplanatic_angle, double isoplanatic_angle_err, double altitude_eff, double altitude_eff_err, double time_constant, double time_constant_err, double moment_0, double moment_0_err, double moment_2, double moment_2_err ) {
+void file::atmosphere( const time_type& time, double free_seeing, double free_seeing_err, double dimm_seeing, double dimm_seeing_err, double dimm_turbulence, double dimm_turbulence_err, double isoplanatic_angle, double isoplanatic_angle_err, double altitude_eff, double altitude_eff_err, double desi_time_constant, double desi_time_constant_err, double time_constant, double time_constant_err, double moment_0, double moment_0_err, double moment_2, double moment_2_err ) {
 	static datum skipped;
 	boost::io::ios_all_saver asaver( m_stm );
 	m_stm << "A " << time << " "
@@ -54,6 +54,7 @@
 		<< std::fixed << datum( altitude_eff, altitude_eff_err/altitude_eff ) << " " << skipped << " "
 		<< datum( isoplanatic_angle, isoplanatic_angle_err/isoplanatic_angle ) << " "
 		<< std::scientific << datum( moment_2, moment_2_err/moment_2 ) << " "
+		<< std::fixed << datum( desi_time_constant, desi_time_constant_err/time_constant ) << " "
 		<< std::fixed << datum( time_constant, time_constant_err/time_constant ) << std::endl;
 }
 void file::profile( const time_type& time, const std::string& method, double cn2seeing, double chi2, const vector_type& cn2, const vector_type& ecn2, const vector_type& grid ){
@@ -115,6 +116,6 @@
 void file::objectinfo( const time_type& time, const object_type& obj ){
 	boost::io::ios_all_saver asaver( m_stm );
 	m_stm << "O " << time << " " << obj << std::endl;
-};
+}
 
 }}

atmos-2.97.1.tar.gz/src/file_output.h -> atmos-2.97.2.tar.gz/src/file_output.h Changed

@@ -1,5 +1,5 @@
 /*
-    $Id: file_output.h 180 2011-01-30 15:56:33Z matwey $
+    $Id: file_output.h 200 2011-02-03 14:45:13Z matwey $
     Copyright (C) 2009  Sternberg Astronomical Institute, MSU
 
     This program is free software: you can redistribute it and/or modify
@@ -50,7 +50,7 @@
 	void comment( const std::string& comment );
 	void normal( const time_type& time, double mz );
 	void background( const time_type& time );
-	void atmosphere( const time_type& time, double free_seeing, double free_seeing_err, double dimm_seeing, double dimm_seeing_err, double dimm_turbulence, double dimm_turbulence_err, double isoplanatic_angle, double isoplanatic_angle_err, double altitude_eff, double altitude_eff_err, double time_constant, double time_constant_err, double moment_0, double moment_0_err, double moment_2, double moment_2_err );
+	void atmosphere( const time_type& time, double free_seeing, double free_seeing_err, double dimm_seeing, double dimm_seeing_err, double dimm_turbulence, double dimm_turbulence_err, double isoplanatic_angle, double isoplanatic_angle_err, double altitude_eff, double altitude_eff_err, double desi_time_constant, double desi_time_constant_err, double time_constant, double time_constant_err, double moment_0, double moment_0_err, double moment_2, double moment_2_err );
 	void profile( const time_type& time, const std::string& method, double cn2seeing, double chi2, const vector_type& cn2, const vector_type& ecn2, const vector_type& grid );
 	void residual( const time_type& time, const std::string& method, const vector_type& o_csi );
 	void fluxes( const time_type& time, const vector_type& result_flux );

atmos-2.97.1.tar.gz/src/indices.cpp -> atmos-2.97.2.tar.gz/src/indices.cpp Changed

@@ -1,5 +1,5 @@
 /*
-    $Id: indices.cpp 181 2011-01-31 10:07:58Z matwey $
+    $Id: indices.cpp 200 2011-02-03 14:45:13Z matwey $
     Copyright (C) 2009  Sternberg Astronomical Institute, MSU
 
     This program is free software: you can redistribute it and/or modify
@@ -36,7 +36,7 @@
 flux_type flux_background_correction( const flux_type& flux, const vector_type& scatter, const vector_type& background ) {
 	return photometry::flux_scatter_correction( flux, scatter ) - background;
 }
-boost::tuple< vars_type, rho1_type, rho2_type > vars_correction( const flux_type& flux, const vars_type& vars, const rho1_type& rho1, const rho2_type& rho2, const vector_type& deadtime, const vector_type& npoisson ) {
+boost::tuple< vars_type, rho1_type, rho2_type > vars_nonlinearity_correction( const flux_type& flux, const vars_type& vars, const rho1_type& rho1, const rho2_type& rho2, const vector_type& deadtime, const vector_type& npoisson, const value_type& exposition ) {
 	using namespace utils::ublas;
 
 	typedef matrix_vector_slice< vars_type > diagonal_type;
@@ -50,10 +50,8 @@
 	
 	result_type ret(vars,rho1,rho2);
 
-	diagonal_type variances( ret.get<0>(),
-		slice(0, 1, std::min( ret.get<0>().size1(), ret.get<0>().size2() )),
-		slice(0, 1, std::min( ret.get<0>().size1(), ret.get<0>().size2() )) );
-	variances -= photometry::variance_nonpoisson_correction( flux, deadtime, npoisson );
+	diagonal_type variances = diag( ret.get<0>() );
+	variances -= photometry::variance_nonpoisson_correction( flux, deadtime, npoisson ) / exposition;
 	
 	flux_type corr = photometry::variance_nonlinearity_correction( flux, deadtime );
 
@@ -63,34 +61,57 @@
 
 	return ret;
 }
-
-boost::tuple< indices_type, indesis_type > calculate( const flux_type& flux, const vars_type& vars, const rho1_type& rho1, const rho2_type& rho2 ) {
+vars_type vars_nonlinearity_correction( const flux_type& flux, const vars_type& vars, const vector_type& deadtime, const vector_type& npoisson, const value_type& exposition ) {
 	using namespace utils::ublas;
 
-	typedef matrix_vector_slice< const indices_type > indices_diagonal_type;
-	typedef matrix_vector_slice< const rho1_type > rho1_diagonal_type;
-	typedef boost::tuple< indices_type, indesis_type> result_type;
+	typedef matrix_vector_slice< vars_type > diagonal_type;
 
-	assert( vars.size1() == rho1.size1() );
-	assert( vars.size2() == rho1.size2() );
 	assert( vars.size1() == flux.size() );
 	assert( vars.size2() == flux.size() );
 
-	result_type ret;
-	
-	rho1_diagonal_type rho1_diag = diag( rho1 );
+	vars_type ret( vars );
+
+	diagonal_type variances = diag( ret );
+	variances -= photometry::variance_nonpoisson_correction( flux, deadtime, npoisson ) / exposition;
+
+	flux_type corr = photometry::variance_nonlinearity_correction( flux, deadtime );
+
+	return element_div( ret, outer_prod( corr, corr ) );
+}
+
+vars_type vars_exposure_correction( const flux_type& flux, const vars_type& vars, const rho1_type& rho1 ) {
+	return vars*1.17 - rho1*0.17;
+}
 
-	ret.get<0>() = utils::ublas::log( xp1( element_div( vars*1.17 - rho1*0.17, outer_prod( flux, flux ) ) ) );
-	ret.get<1>() = ( diag( element_div( vars, outer_prod( flux, flux ) ) )*3.0 + element_div( rho2 - rho1_diag * 4.0, element_prod(flux,flux) ) ) / 4.5;
+indices_type calculate_indices( const flux_type& flux, const vars_type& vars ) {
+	using namespace utils::ublas;
+
+	typedef matrix_vector_slice< indices_type >    indices_diagonal_type;
+
+	assert( vars.size1() == flux.size() );
+	assert( vars.size2() == flux.size() );
+
+	indices_type ret = utils::ublas::log( xp1( element_div( vars, outer_prod( flux, flux ) ) ) );
 
 	return ret;
 }
+indesis_type calculate_indesis( const flux_type& flux, const vars_type& vars, const rho1_type& rho1, const rho2_type& rho2 ) {
+	using namespace utils::ublas;
+
+	typedef matrix_vector_slice< const vars_type > vars_diagonal_type;
+	typedef matrix_vector_slice< const rho1_type > rho1_diagonal_type;
 
+	rho1_diagonal_type rho1_diag = diag( rho1 );
+	vars_diagonal_type variances = diag( vars );
+	
+	return ( element_div( variances*3.0 + rho2 - rho1_diag * 4.0, element_prod(flux,flux) ) ) / 4.5;
+}
 storage::storage( size_type reserve ):
 	table_storage< detail::sample >( reserve ),
 	m_flux_column( *this, boost::mem_fn(&detail::sample::flux) ),
 	m_indices_column( *this, boost::mem_fn(&detail::sample::indices) ),
-	m_indesis_column( *this, boost::mem_fn(&detail::sample::indesis) ) {
+	m_indesis_column( *this, boost::mem_fn(&detail::sample::indesis) ),
+	m_wind_indices_column( *this, boost::mem_fn(&detail::sample::wind_indices) ) {
 }
 
 namespace filter {

atmos-2.97.1.tar.gz/src/indices.h -> atmos-2.97.2.tar.gz/src/indices.h Changed

@@ -1,5 +1,5 @@
 /*
-    $Id: indices.h 180 2011-01-30 15:56:33Z matwey $
+    $Id: indices.h 200 2011-02-03 14:45:13Z matwey $
     Copyright (C) 2009  Sternberg Astronomical Institute, MSU
 
     This program is free software: you can redistribute it and/or modify
@@ -45,9 +45,11 @@
 
 	flux_type flux_nonlinearity_correction( const flux_type& flux, const vector_type& deadtime );
 	flux_type flux_background_correction( const flux_type& flux, const vector_type& scatter, const vector_type& background );
-	boost::tuple< vars_type, rho1_type, rho2_type >
-		vars_correction( const flux_type& flux, const vars_type& vars, const rho1_type& rho1, const rho2_type& rho2, const vector_type& deadtime, const vector_type& npoisson );
-	boost::tuple< indices_type, indesis_type > calculate( const flux_type& flux, const vars_type& vars, const rho1_type& rho1, const rho2_type& rho2 );
+	boost::tuple< vars_type, rho1_type, rho2_type > vars_nonlinearity_correction( const flux_type& flux, const vars_type& vars, const rho1_type& rho1, const rho2_type& rho2, const vector_type& deadtime, const vector_type& npoisson, const value_type& exposition );
+	vars_type vars_nonlinearity_correction( const flux_type& flux, const vars_type& vars, const vector_type& deadtime, const vector_type& npoisson, const value_type& exposition );
+	vars_type vars_exposure_correction( const flux_type& flux, const vars_type& vars, const rho1_type& rho1 );
+	indesis_type calculate_indesis( const flux_type& flux, const vars_type& vars, const rho1_type& rho1, const rho2_type& rho2 );
+	indices_type calculate_indices( const flux_type& flux, const vars_type& vars );
 
 namespace detail {
 	struct sample {
@@ -59,13 +61,15 @@
 		timestamp_type timestamp;
 		flux_type      flux;
 		indices_type   indices;
+		indices_type   wind_indices;
 		indesis_type   indesis;
 
 		sample( const timestamp_type& time,
 			const flux_type& flux,
 			const indices_type& indices,
-			const indesis_type& indesis ):
-			timestamp(time), flux(flux), indices(indices), indesis(indesis) {}
+			const indesis_type& indesis,
+			const indices_type& wind_indices ):
+			timestamp(time), flux(flux), indices(indices), indesis(indesis), wind_indices(wind_indices) {}
 	};
 }
 
@@ -81,6 +85,7 @@
 	flux_column_type    m_flux_column;
 	indices_column_type m_indices_column;
 	indesis_column_type m_indesis_column;
+	indices_column_type m_wind_indices_column;
 
 public:
 	storage( size_type reserve = 60 );
@@ -88,13 +93,15 @@
 	const flux_column_type&    flux()    const { return m_flux_column; }
 	const indices_column_type& indices() const { return m_indices_column; }
 	const indesis_column_type& indesis() const { return m_indesis_column; }
+	const indices_column_type& wind_indesis() const { return m_wind_indices_column; }
 
 	void push_back( const time_type& time,
 			const flux_type& flux,
 			const indices_type& indices,
-			const indesis_type& indesis ) {
+			const indesis_type& indesis,
+			const indices_type& wind_indices ) {
 		
-		table_storage< detail::sample >::push_back( detail::sample( time, flux, indices, indesis ) );
+		table_storage< detail::sample >::push_back( detail::sample( time, flux, indices, indesis, wind_indices ) );
 	}
 };

atmos-2.97.1.tar.gz/src/integrals.cpp -> atmos-2.97.2.tar.gz/src/integrals.cpp Changed

@@ -1,5 +1,5 @@
 /*
-    $Id: integrals.cpp 180 2011-01-30 15:56:33Z matwey $
+    $Id: integrals.cpp 200 2011-02-03 14:45:13Z matwey $
     Copyright (C) 2009  Sternberg Astronomical Institute, MSU
 
     This program is free software: you can redistribute it and/or modify
@@ -36,6 +36,21 @@
 namespace atmos {
 namespace integrals {
 
+namespace {
+
+// FIXME:
+vector_type reorder_flat( const atmos::indices::indices_type& indices ) {
+	assert( indices.size1() == atmos::channels );
+	assert( indices.size2() == atmos::channels );
+
+	vector_type inds( atmos::indices_size );
+	inds(0) = indices(0,0); inds(1) = indices(1,1); inds(2) = indices(2,2); inds(3) = indices(3,3);
+	inds(4) = indices(0,1); inds(5) = indices(0,2); inds(6) = indices(0,3);
+	inds(7) = indices(1,2); inds(8) = indices(1,3); inds(9) = indices(2,3);
+	return inds;
+}
+
+}
 const double KM2M = 1000;
 const double POWSEE = 0 ;       // Power of Cn2 moment for seeing
 const double POWISP = 5.0/3.0;  // Power of Cn2 moment for isoplanatic angle
@@ -45,9 +60,11 @@
 const double MINRELW = 6e-4 ;   // Minimal ratio of the singular value to the maximal one
 
 static const double powers = { POWSEE, POWISP, POWEFF, POWISK };
-
 const unsigned int npower = sizeof(powers)/sizeof(powers0);
 
+static const double wind_powers = { 0.0 };
+const unsigned int wind_npower = sizeof(wind_powers)/sizeof(wind_powers0);
+
 static const double HBOUND = 0.38;// Height of the boundary layer kilometers to fit the Cn2-integral
 static const double HTOP = 25.0; // Height of the hight layer kilometers
 
@@ -67,15 +84,16 @@
                                 // \f$ KP = (2.905*(2\pi)^2*(MKM2M)^{-2})^{-3/5}*206265 \f$*/
 const double KT = 0.175;        // Atmospheric time constant empirical calibration constant for DESI for 1,3ms
                                 // integrations in diameter=2cm apertures with -1km shift
+const double KL = 0.058 * 2.7464013582652968e-08;
 
 const double SPOW_TC  = -0.6;   // Power of DESI for atm. time constant
 
-functor_factory::dcf_matrix_type functor_factory::calculate( const std::string& spectral_type, const weight_type& awf ) {
+functor_factory::dcf_matrix_type functor_factory::calculate( const weight_type& wf, unsigned int npower, const double* powers ) {
 	using namespace boost::numeric::ublas;
 	using namespace utils::algo;
 	using namespace utils::ublas;
 
-	const vector<double> grid( vector_generator< lg_generator > ( lg_generator(0.0, 0.04, 0.25), 50 ) );
+	const vector<double> grid( utils::ublas::vector_generator< lg_generator > ( lg_generator(0.0, 0.04, 0.25), 50 ) );
 
 	typedef dcf_matrix_type            result_type;
 	typedef dcf_matrix_type::size_type size_type;
@@ -91,7 +109,7 @@
 	result_type res( npower, atmos::indices_size );
 
 	vector_range< const vector<double> > h( grid, range( lower.index(), higher.index()+1 ) );
-	matrix< double > wfv = vector_of_vector_to_matrix( project( grid, awf) );
+	matrix< double > wfv = vector_of_vector_to_matrix( project( grid, wf) );
 
 	for( unsigned int ipow = 0; ipow < npower; ipow++ ) {
 		double moment_power = powersipow;
@@ -143,17 +161,17 @@
 }
 functor_factory::~functor_factory() {
 }
-functor_factory::functor_type functor_factory::operator() ( const std::string& spectral_type, const weight_type& awf ) {
+const functor_factory::functor_type& functor_factory::operator() ( const std::string& spectral_type, const weight_type& null_wf, const weight_type& second_wf ) {
 	storage_type::const_iterator it = m_storage.find( spectral_type );
 	if( it != m_storage.end() )
-		return functor_type( it->second );
+		return it->second;
 
-	std::pair< storage_type::iterator, bool > ret = m_storage.insert( std::make_pair( spectral_type, calculate( spectral_type, awf ) ) );
+	std::pair< storage_type::iterator, bool > ret = m_storage.insert( std::make_pair( spectral_type, functor_type( calculate( null_wf, npower, powers ), calculate( second_wf, wind_npower, wind_powers ) ) ) );
 
-	return functor_type( ret.first->second );
+	return ret.first->second;
 }
 
-functor::result_type functor::operator() ( const indices_type& indices, const indesis_type& indesis, double mz ) const {
+functor::result_type functor::operator() ( const indices_type& indices, const indesis_type& indesis, const indices_type& wind_indices, double mz ) const {
 	typedef vector< double >       vector_type;
 	typedef vector_type::size_type size_type;
 	
@@ -161,22 +179,27 @@
 	assert( indices.size2() == atmos::channels );
 	assert( indesis.size()  == atmos::channels );
 
-	vector_type aints( atmos::indices_size );
-	aints(0) = indices(0,0); aints(1) = indices(1,1); aints(2) = indices(2,2); aints(3) = indices(3,3);
-	aints(4) = indices(0,1); aints(5) = indices(0,2); aints(6) = indices(0,3);
-	aints(7) = indices(1,2); aints(8) = indices(1,3); aints(9) = indices(2,3);
-
-	aints = prod( m_dcf, aints );
+	vector_type aints = reorder_flat(indices);
+	aints = prod( m_indices_dcf, aints );
 
 	assert( npower == aints.size() );
 	for( size_type i = 0; i < npower; ++i ){
 		aints(i) = aints(i) / std::pow( mz, 1.0 + powersi );
 	}
 
+	vector_type wind_aints = reorder_flat(wind_indices);
+	wind_aints = prod( m_wind_indices_dcf, wind_aints );
+		
+	assert( wind_npower == wind_aints.size() );
+	for( size_type i = 0; i < wind_npower; ++i ){
+		wind_aints(i) = wind_aints(i) / std::pow( mz, 1.0 + wind_powersi );
+	}
+
 	return result_type( KS * LEFF_SEE * std::pow( aints(0), JPOW_SEE ),
 			    KP * LEFF_ISP * std::pow( aints(1), JPOW_ISP ),
 			    0.001 * std::pow( aints(1) / aints(2), 1.0/( POWISP-POWEFF ) ),
 			    KT * std::pow( indesis(0), SPOW_TC ),
+			    1e3/*ms*/ * KL * std::pow( aints(0), -0.1 ) / std::sqrt( 6e6 / 6 * wind_aints(0) ),
 			    aints(0),
 			    aints(3) );
 }
@@ -186,6 +209,7 @@
 	m_free_seeing_column( *this, boost::mem_fn(&detail::sample::free_seeing) ),
 	m_isoplanatic_angle_column( *this, boost::mem_fn(&detail::sample::isoplanatic_angle) ),
 	m_altitude_eff_column( *this, boost::mem_fn(&detail::sample::altitude_eff) ),
+	m_desi_time_constant_column( *this, boost::mem_fn(&detail::sample::desi_time_constant) ),
 	m_time_constant_column( *this, boost::mem_fn(&detail::sample::time_constant) ),
 	m_moment_0_column( *this, boost::mem_fn(&detail::sample::moment_0) ),
 	m_moment_2_column( *this, boost::mem_fn(&detail::sample::moment_2) ) {
@@ -197,6 +221,7 @@
 	return ! ( is_regular( x.free_seeing ) &&
 		is_regular( x.isoplanatic_angle ) &&
 		is_regular( x.altitude_eff ) &&
+		is_regular( x.desi_time_constant ) &&
 		is_regular( x.time_constant ) &&
 		is_regular( x.moment_0 ) &&
 		is_regular( x.moment_2 ) );

atmos-2.97.1.tar.gz/src/integrals.h -> atmos-2.97.2.tar.gz/src/integrals.h Changed

@@ -1,5 +1,5 @@
 /*
-    $Id: integrals.h 180 2011-01-30 15:56:33Z matwey $
+    $Id: integrals.h 200 2011-02-03 14:45:13Z matwey $
     Copyright (C) 2009  Sternberg Astronomical Institute, MSU
 
     This program is free software: you can redistribute it and/or modify
@@ -33,6 +33,7 @@
 	using namespace boost::numeric::ublas;
 	
 	typedef double value_type;
+	typedef vector< value_type > vector_type;
 	typedef matrix< value_type > matrix_type;
 
 namespace detail {
@@ -42,6 +43,7 @@
 		value_type free_seeing;
 		value_type isoplanatic_angle;
 		value_type altitude_eff;
+		value_type desi_time_constant;
 		value_type time_constant;
 		value_type moment_0;
 		value_type moment_2;
@@ -49,10 +51,11 @@
 		sample( const value_type& free_seeing,
 			const value_type& isoplanatic_angle,
 			const value_type& altitude_eff,
+			const value_type& desi_time_constant,
 			const value_type& time_constant,
 			const value_type& moment_0,
 			const value_type& moment_2 ):
-			free_seeing(free_seeing), isoplanatic_angle(isoplanatic_angle), altitude_eff(altitude_eff), time_constant(time_constant), moment_0(moment_0), moment_2(moment_2) {
+			free_seeing(free_seeing), isoplanatic_angle(isoplanatic_angle), altitude_eff(altitude_eff), desi_time_constant(desi_time_constant), time_constant(time_constant), moment_0(moment_0), moment_2(moment_2) {
 		}
 	};
 }
@@ -67,6 +70,7 @@
 	typedef common_column_type free_seeing_column_type;
 	typedef common_column_type isoplanatic_angle_column_type;
 	typedef common_column_type altitude_eff_column_type;
+	typedef common_column_type desi_time_constant_column_type;
 	typedef common_column_type time_constant_column_type;
 	typedef common_column_type moment_0_column_type;
 	typedef common_column_type moment_2_column_type;
@@ -74,6 +78,7 @@
 	free_seeing_column_type       m_free_seeing_column;
 	isoplanatic_angle_column_type m_isoplanatic_angle_column;
 	altitude_eff_column_type      m_altitude_eff_column;
+	desi_time_constant_column_type m_desi_time_constant_column;
 	time_constant_column_type     m_time_constant_column;
 	moment_0_column_type          m_moment_0_column;
 	moment_2_column_type          m_moment_2_column;
@@ -83,6 +88,7 @@
 	const free_seeing_column_type&       free_seeing()       const { return m_free_seeing_column; }
 	const isoplanatic_angle_column_type& isoplanatic_angle() const { return m_isoplanatic_angle_column; }
 	const altitude_eff_column_type&      altitude_eff()      const { return m_altitude_eff_column; }
+	const desi_time_constant_column_type& desi_time_constant() const { return m_desi_time_constant_column; }
 	const time_constant_column_type&     time_constant()     const { return m_time_constant_column; }
 	const moment_0_column_type&          moment_0()          const { return m_moment_0_column; }
 	const moment_2_column_type&          moment_2()          const { return m_moment_2_column; }
@@ -110,12 +116,12 @@
 	public:
 		typedef matrix_type dcf_matrix_type;
 		typedef boost::function1< vector<double>, double > weight_type;
-		typedef std::map< std::string, dcf_matrix_type, utils::algo::nocase_less< std::string > > storage_type;
 		typedef functor functor_type;
+		typedef std::map< std::string, functor_type, utils::algo::nocase_less< std::string > > storage_type;
 	private:
 		storage_type m_storage;
 
-		dcf_matrix_type calculate( const std::string& spectral_type, const weight_type& awf );
+		dcf_matrix_type calculate( const weight_type& wf, unsigned int npower, const double* powers );
 	public:
 		functor_factory();
 		functor_factory( const functor_factory& factory );
@@ -123,7 +129,7 @@
 		functor_factory& operator=( const functor_factory& factory );
 		~functor_factory();
 
-		functor_type operator() ( const std::string& spectral_type, const weight_type& awf );
+		const functor_type& operator() ( const std::string& spectral_type, const weight_type& null_wf, const weight_type& second_wf );
 };
 
 class functor {
@@ -134,12 +140,13 @@
 		typedef atmos::indices::indesis_type indesis_type;
 		typedef atmos::integrals::detail::sample result_type;
 	private:
-		const dcf_matrix_type& m_dcf;
+		dcf_matrix_type m_indices_dcf;
+		dcf_matrix_type m_wind_indices_dcf;
 	public:
-		functor( const dcf_matrix_type& dcf ):
-			m_dcf( dcf ) {}
+		functor( const dcf_matrix_type& indices_dcf, const dcf_matrix_type& wind_indices_dcf ):
+			m_indices_dcf( indices_dcf ), m_wind_indices_dcf( wind_indices_dcf ) {}
 
-		result_type operator() ( const indices_type& indices, const indesis_type& indesis, double mz ) const;
+		result_type operator() ( const indices_type& indices, const indesis_type& indesis, const indices_type& wind_indices, double mz ) const;
 };
 
 }}

atmos-2.97.1.tar.gz/src/lambert.h -> atmos-2.97.2.tar.gz/src/lambert.h Changed

atmos-2.97.1.tar.gz/src/photometry.h -> atmos-2.97.2.tar.gz/src/photometry.h Changed

@@ -1,5 +1,5 @@
 /*
-    $Id: photometry.h 180 2011-01-30 15:56:33Z matwey $
+    $Id: photometry.h 203 2011-02-03 17:40:12Z matwey $
     Copyright (C) 2010  Sternberg Astronomical Institute, MSU
 
     This program is free software: you can redistribute it and/or modify
@@ -22,6 +22,7 @@
 #include <boost/numeric/ublas/vector.hpp>
 
 #include "apply_to_all.h"
+#include "lambert.h"
 
 namespace atmos {
 namespace photometry {
@@ -34,13 +35,21 @@
 	template<class V1,class V2>
 	vector< typename V1::value_type >
 	flux_nonlinearity_correction( const vector_expression<V1>& flux, const vector_expression<V2>& deadtime ) {
-		using namespace utils::ublas;
-		
+		using namespace utils::math;
+
 		typedef typename V1::value_type value_type;
+		typedef vector< value_type > vector_type;
+		typedef typename vector_type::size_type size_type;
 		
 		assert( flux().size() == deadtime().size() );
 
-		return -element_div( lambert_W0( -element_prod( flux, deadtime ) ), deadtime );
+		vector_type result(flux);
+
+		for( size_type i = 0; i < result.size(); ++i )
+			if( deadtime()(i) != value_type(0) )
+				result(i) = -lambert_W0( -result(i) * deadtime()(i) ) / deadtime()(i);
+
+		return result;
 	}
 
 	/**

atmos-2.97.1.tar.gz/src/profile.cpp -> atmos-2.97.2.tar.gz/src/profile.cpp Changed

atmos-2.97.1.tar.gz/src/spectral_factory.cpp -> atmos-2.97.2.tar.gz/src/spectral_factory.cpp Changed

@@ -1,5 +1,5 @@
 /*
-    $Id: spectral_factory.cpp 180 2011-01-30 15:56:33Z matwey $
+    $Id: spectral_factory.cpp 200 2011-02-03 14:45:13Z matwey $
     Copyright (C) 2009  Sternberg Astronomical Institute, MSU
 
     This program is free software: you can redistribute it and/or modify
@@ -164,9 +164,9 @@
 	boost::filesystem::path spect( m_wdir / (spectral_type + ".sp") );
 	return spect;
 }
-const factory::weight_type& factory::weight( const std::string& spectype, const vector<double>& inner_diam, const vector<double>& outer_diam, double dimmbase, double dimmsize ) {
-	storage_type::const_iterator it = m_storage.find( spectype );
-	if( it != m_storage.end() )
+const factory::weight_type& factory::nullexposure_weight( const std::string& spectype, const vector<double>& inner_diam, const vector<double>& outer_diam, double dimmbase, double dimmsize ) {
+	storage_type::const_iterator it = m_nullexposure.find( spectype );
+	if( it != m_nullexposure.end() )
 		return it->second;
 
 	std::pair< vector_type, vector_type > spectral_response = load_reaction( get_spectrum_filename( spectype ) );
@@ -174,7 +174,21 @@
 	std::pair< vector_type, vector_type > mass_response = photon_distribution( spectral_response.first, spectral_response.second, m_mass_reaction.first, m_mass_reaction.second );
 	std::pair< vector_type, vector_type > dimm_response = photon_distribution( spectral_response.first, spectral_response.second, m_dimm_reaction.first, m_dimm_reaction.second );
 
-	std::pair< storage_type::iterator, bool > ret = m_storage.insert( std::make_pair( spectype, weight_type( m_basis_grid, utils::ublas::project( m_basis_grid , atmos::weight::massdimm_weight( mass_response.first, mass_response.second, dimm_response.first, dimm_response.second, inner_diam, outer_diam, dimmsize, dimmbase ) ) ) ) );
+	std::pair< storage_type::iterator, bool > ret = m_nullexposure.insert( std::make_pair( spectype, weight_type( m_basis_grid, utils::ublas::project( m_basis_grid , atmos::weight::massdimm_weight( mass_response.first, mass_response.second, dimm_response.first, dimm_response.second, inner_diam, outer_diam, dimmsize, dimmbase ) ) ) ) );
+
+	return ret.first->second;
+}
+
+const factory::weight_type& factory::shortexposure_weight( const std::string& spectype, const vector<double>& inner_diam, const vector<double>& outer_diam ) {
+	storage_type::const_iterator it = m_shortexposure.find( spectype );
+	if( it != m_shortexposure.end() )
+		return it->second;
+
+	std::pair< vector_type, vector_type > spectral_response = load_reaction( get_spectrum_filename( spectype ) );
+
+	std::pair< vector_type, vector_type > mass_response = photon_distribution( spectral_response.first, spectral_response.second, m_mass_reaction.first, m_mass_reaction.second );
+
+	std::pair< storage_type::iterator, bool > ret = m_shortexposure.insert( std::make_pair( spectype, weight_type( m_basis_grid, utils::ublas::project( m_basis_grid , atmos::weight::mass_shortexposure_weight( mass_response.first, mass_response.second, inner_diam, outer_diam) ) ) ) );
 
 	return ret.first->second;
 }

atmos-2.97.1.tar.gz/src/spectral_factory.h -> atmos-2.97.2.tar.gz/src/spectral_factory.h Changed

@@ -1,5 +1,5 @@
 /*
-    $Id: spectral_factory.h 180 2011-01-30 15:56:33Z matwey $
+    $Id: spectral_factory.h 200 2011-02-03 14:45:13Z matwey $
     Copyright (C) 2009  Sternberg Astronomical Institute, MSU
 
     This program is free software: you can redistribute it and/or modify
@@ -60,7 +60,8 @@
 
 	vector_generator< lg_generator > m_basis_grid;
 
-	storage_type m_storage;
+	storage_type m_nullexposure;
+	storage_type m_shortexposure;
 
 	std::pair< vector_type, vector_type > m_mass_reaction;
 	std::pair< vector_type, vector_type > m_dimm_reaction;
@@ -70,10 +71,13 @@
 	factory(const boost::filesystem::path& wd, const boost::filesystem::path& mass_response, const boost::filesystem::path& ccd_response );
 	~factory();
 
-	const weight_type& weight( const std::string& spectype, const vector_type& inner_diam, const vector_type& outer_diam, double dimmbase, double dimmsize );
+	const weight_type& nullexposure_weight( const std::string& spectype, const boost::numeric::ublas::vector<double>& inner_diam, const boost::numeric::ublas::vector<double>& outer_diam, double dimmbase, double dimmsize );
+	const weight_type& shortexposure_weight( const std::string& spectype, const boost::numeric::ublas::vector<double>& inner_diam, const boost::numeric::ublas::vector<double>& outer_diam );
 
-	const const_iterator begin() const { return m_storage.begin(); }
-	const const_iterator end()   const { return m_storage.end(); }
+	const const_iterator nullexposure_begin() const { return m_nullexposure.begin(); }
+	const const_iterator nullexposure_end()   const { return m_nullexposure.end(); }
+	const const_iterator shortexposure_begin() const { return m_shortexposure.begin(); }
+	const const_iterator shortexposure_end()   const { return m_shortexposure.end(); }
 };
 
 }}

atmos-2.97.1.tar.gz/src/spline.h -> atmos-2.97.2.tar.gz/src/spline.h Changed

atmos-2.97.1.tar.gz/src/table_storage.h -> atmos-2.97.2.tar.gz/src/table_storage.h Changed

atmos-2.97.1.tar.gz/src/weif.cpp -> atmos-2.97.2.tar.gz/src/weif.cpp Changed

atmos-2.97.1.tar.gz/src/weif.h -> atmos-2.97.2.tar.gz/src/weif.h Changed

atmos-2.97.1.tar.gz/src/worksheet.cpp -> atmos-2.97.2.tar.gz/src/worksheet.cpp Changed

@@ -1,5 +1,5 @@
 /*
-    $Id: worksheet.cpp 180 2011-01-30 15:56:33Z matwey $
+    $Id: worksheet.cpp 203 2011-02-03 17:40:12Z matwey $
     Copyright (C) 2009  Sternberg Astronomical Institute, MSU
 
     This program is free software: you can redistribute it and/or modify
@@ -57,6 +57,7 @@
 
 	double airmass( const time_type& time ) const;
 	const weight_type& weight();
+	const weight_type& short_weight();
 	const std::string& sclass() const { return m_spectral; }
 };
 
@@ -120,6 +121,7 @@
 	const atmos::parameters::vector_type m_deadtime;
 	const atmos::parameters::vector_type m_scatter;
 	const atmos::parameters::vector_type m_npoisson;
+	const double m_exposition;	
 
 	atmos::indices::filter::flux_threshold m_flux_filter;
 public:
@@ -130,11 +132,14 @@
 		m_deadtime( ctx.parameters().nonlinearity() ),
 		m_scatter( ctx.parameters().scatter() ),
 		m_npoisson( ctx.parameters().nonpoisson() ),
+		m_exposition( ctx.parameters().exposition() ),
 		m_flux_filter( ctx.m_va"flux-threshold".as<double>() ) {
 
 	}
 	virtual ~normal_accumulator() {};
 	virtual void accumulate( const time_type& time, const flux_type& mean, const vars_type& vars, const rho1_type& rho1, const rho2_type& rho2 ) {
+		/* mean is in photons
+		   flux is in photons/msec */
 		using namespace atmos::indices;
 		
 		if( m_accumulation_begin == boost::posix_time::not_a_date_time )
@@ -144,20 +149,25 @@
 		m_num++;
 		m_mean_mz = m_mean_mz * ( double(m_num-1) / double(m_num) ) + context().curobject().airmass( time ) * ( 1.0/m_num );
 
-		flux_type flux_corrected = flux_nonlinearity_correction( mean, m_deadtime );
+		flux_type flux = mean / m_exposition;
+		flux_type flux_corrected = flux_nonlinearity_correction( flux, m_deadtime );
 
-		const boost::tuple< vars_type, rho1_type, rho2_type >& vars_corrected = vars_correction(
+		const boost::tuple< vars_type, rho1_type, rho2_type >& vars_corrected = vars_nonlinearity_correction(
 			flux_corrected,
-			element_prod( vars, outer_prod( mean, mean ) ),
-			element_prod( rho1, outer_prod( mean, mean ) ),
-			element_prod( rho2, element_prod( mean, mean ) ), m_deadtime, m_npoisson );
+			element_prod( vars, outer_prod( flux, flux ) ),
+			element_prod( rho1, outer_prod( flux, flux ) ),
+			element_prod( rho2, element_prod( flux, flux ) ), m_deadtime, m_npoisson, m_exposition );
 
 		flux_corrected = flux_background_correction( flux_corrected, m_scatter, context().background_factory().background( time ) );
 
 		if( ! m_flux_filter( flux_corrected ) ) {
-			const boost::tuple< indices_type, indesis_type>& inds = calculate( flux_corrected, vars_corrected.get<0>(), vars_corrected.get<1>(), vars_corrected.get<2>() );
+			const vars_type& vars_nullexp = vars_exposure_correction( flux_corrected, vars_corrected.get<0>(), vars_corrected.get<1>() );
+			
+			const indesis_type& indesis         = calculate_indesis( flux_corrected, vars_corrected.get<0>(), vars_corrected.get<1>(), vars_corrected.get<2>() );
+			const indices_type& indecis_nullexp = calculate_indices( flux_corrected, vars_nullexp );
+			const indices_type& indecis         = calculate_indices( flux_corrected, vars_corrected.get<0>() );
 
-			m_indices.push_back( time, flux_corrected, inds.get<0>(), inds.get<1>() );
+			m_indices.push_back( time, flux_corrected, indecis_nullexp, indesis, indecis - calculate_indices(flux_corrected, vars_corrected.get<1>()) );
 		}
 	}
 	virtual void dimm_accumulate( const time_type& time, unsigned int framenum, const dimm_flux_type& left, const dimm_flux_type& right, const dimm_vars_type& xvars, const dimm_flux_type& yvars ) {
@@ -224,11 +234,11 @@
 			context().m_file_output.indices( meantime, indices_avg, indices_err, indesis_avg, indesis_err );
 			context().m_file_output.fluxes( meantime, flux_avg, flux_err );
 
-			integrals_functor_type ifun = context().m_integrals_factory( context().curobject().sclass(), context().curobject().weight() );
+			integrals_functor_type ifun = context().m_integrals_factory( context().curobject().sclass(), context().curobject().weight(), context().curobject().short_weight() );
 
 			integrals::storage integrals;
 			for( indices::storage::const_iterator it = m_indices.begin(); it != m_indices.end(); ++it ) {
-				integrals.push_back( ifun( it->indices, it->indesis, m_mean_mz ) );
+				integrals.push_back( ifun( it->indices, it->indesis, it->wind_indices, m_mean_mz ) );
 			}
 
 			integrals.unstable_remove_if( integrals::filter::regular() );
@@ -236,12 +246,14 @@
 			double free_seeing_avg       = average( integrals.free_seeing().begin(),       integrals.free_seeing().end(),       mean< double >() );
 			double isoplanatic_angle_avg = average( integrals.isoplanatic_angle().begin(), integrals.isoplanatic_angle().end(), mean< double >() );
 			double altitude_eff_avg      = average( integrals.altitude_eff().begin(),      integrals.altitude_eff().end(),      mean< double >() );
+			double desi_time_constant_avg= average( integrals.desi_time_constant().begin(),integrals.desi_time_constant().end(),mean< double >() );
 			double time_constant_avg     = average( integrals.time_constant().begin(),     integrals.time_constant().end(),     mean< double >() );
 			double moment_0_avg          = average( integrals.moment_0().begin(),          integrals.moment_0().end(),          mean< double >() );
 			double moment_2_avg          = average( integrals.moment_2().begin(),          integrals.moment_2().end(),          mean< double >() );
 			double free_seeing_err       = average( integrals.free_seeing().begin(),       integrals.free_seeing().end(),       error_of_mean< double >(free_seeing_avg) );
 			double isoplanatic_angle_err = average( integrals.isoplanatic_angle().begin(), integrals.isoplanatic_angle().end(), error_of_mean< double >(isoplanatic_angle_avg) );
 			double altitude_eff_err      = average( integrals.altitude_eff().begin(),      integrals.altitude_eff().end(),      error_of_mean< double >(altitude_eff_avg) );
+			double desi_time_constant_err= average( integrals.desi_time_constant().begin(),integrals.desi_time_constant().end(),error_of_mean< double >(desi_time_constant_avg) );
 			double time_constant_err     = average( integrals.time_constant().begin(),     integrals.time_constant().end(),     error_of_mean< double >(time_constant_avg) );
 			double moment_0_err          = average( integrals.moment_0().begin(),          integrals.moment_0().end(),          error_of_mean< double >(moment_0_avg) );
 			double moment_2_err          = average( integrals.moment_2().begin(),          integrals.moment_2().end(),          error_of_mean< double >(moment_2_avg) );
@@ -252,6 +264,7 @@
 				dimm_turbulence, dimm_turbulence_err,
 				isoplanatic_angle_avg, isoplanatic_angle_err,
 				altitude_eff_avg, altitude_eff_err,
+				desi_time_constant_avg, desi_time_constant_err,
 				time_constant_avg, time_constant_err,
 				moment_0_avg, moment_0_err,
 				moment_2_avg, moment_2_err );
@@ -310,13 +323,19 @@
 	return secz1 - secz2*( 0.0018167 + secz2*( 0.002875 + 0.0008083*secz2 ) );
 }
 const worksheet::object::weight_type& worksheet::object::weight() {
-	return context().m_weight_factory.weight(
+	return context().m_weight_factory.nullexposure_weight(
 		m_spectral,
 		context().m_pproxy.inner_diam(),
 		context().m_pproxy.outer_diam(),
 		context().m_pproxy.dimm_aperture_base(),
 		context().m_pproxy.dimm_aperture_size() );
 }
+const worksheet::object::weight_type& worksheet::object::short_weight() {
+	return context().m_weight_factory.shortexposure_weight(
+		m_spectral,
+		context().m_pproxy.inner_diam(),
+		context().m_pproxy.outer_diam() );
+}
 
 void validate(boost::any& v, const std::vector<std::string>& values, worksheet::profile_grid* target_type, int) {
 	namespace po = boost::program_options;
@@ -370,6 +389,7 @@
 		( "accum-time,a",    po::value<double>()->default_value(60.0), "force accumulation time" )
 		( "background-threshold,b", po::value<double>()->default_value(4.0), "set background threshold" )
 		( "chi2-threshold,x", po::value<double>()->default_value(25.0), "set chi2 threshold" )
+		( "enable-simple-background", "do not have regard to the sun's altitude")
 		( "flux-threshold",  po::value<double>()->default_value(2.0), "set flux threshold" )
 		( "profile-grid,Z",  po::value<profile_grid>()->default_value( profile_grid( utils::ublas::vector_generator< utils::algo::lg_generator >( utils::algo::lg_generator( 0.5, 1.0, 2.0 ),6 ) ), "0.5,1,2,4,8,16" ), "set for profile grid generator" )
 		( "mass-response,r", po::value<std::string>()->default_value("mass_maid.crv"), "file with mass detector response" )
@@ -398,7 +418,7 @@
 	}
 	return 0;
 }
-atmos::background::model::abstract* worksheet::select_background_model( const boost::posix_time::ptime& time ) const {
+atmos::background::model::abstract* worksheet::select_background_model_smart( const boost::posix_time::ptime& time ) const {
 	using namespace atmos::background::model;
 	
 	double alt = sun_altitude( time ) * 180 / M_PI; // degrees
@@ -410,6 +430,10 @@
 	// else
 	return new night( time );
 }
+atmos::background::model::abstract* worksheet::select_background_model_simple( const boost::posix_time::ptime& time ) const {
+	using namespace atmos::background::model;
+	return new night( time );
+}
 double worksheet::sun_altitude( const boost::posix_time::ptime& time ) const {
 	yalpa::ecliptic<double> sun = yalpa::sun_position( time );
 	yalpa::horizon<double> sun_hz( make_horizon_helper( sun, m_pproxy.site(), yalpa::sidereal_time(time) ) );
@@ -420,7 +444,7 @@
 	m_pcontainer(),
 	m_pproxy( m_pcontainer ),
 	m_weight_factory( boost::filesystem::path( va"spectra-dir".as<std::string>() ), boost::filesystem::path( va"mass-response".as<std::string>() ), boost::filesystem::path( va"ccd-response".as<std::string>() ) ),
-	m_background_factory( boost::bind( std::mem_fun( &worksheet::select_background_model ), this, _1) ),
+	m_background_factory( boost::bind( std::mem_fun( va.count("enable-simple-background") ? &worksheet::select_background_model_simple : &worksheet::select_background_model_smart ), this, _1) ),
 	m_file_input( std::cin, *this ),
 	m_file_output( std::cout )
 {
@@ -442,9 +466,12 @@
 }
 void worksheet::post() {
 	if( m_va.count("dump-weight") ) {
-		for( weight_factory_type::const_iterator it = weight_factory().begin(); it != weight_factory().end(); ++it ) {
+		for( weight_factory_type::const_iterator it = weight_factory().nullexposure_begin(); it != weight_factory().nullexposure_end(); ++it ) {
 			atmos::weight::dump_weight( ( boost::format("weight-%1%.wf") % it->first ).str(), it->second );
 		}
+		for( weight_factory_type::const_iterator it = weight_factory().shortexposure_begin(); it != weight_factory().shortexposure_end(); ++it ) {
+			atmos::weight::dump_weight( ( boost::format("weight-short-%1%.wf") % it->first ).str(), it->second );
+		}
 	}
 }
 worksheet::~worksheet() {

atmos-2.97.1.tar.gz/src/worksheet.h -> atmos-2.97.2.tar.gz/src/worksheet.h Changed

Changes of Revision 11