##// END OF EJS Templates
Sync
paul -
r4:294fc10efc0e default
parent child
Show More
@@ -43,7 +43,8
43
43
44 float k_f0[NB_BINS_COMPRESSED_MATRIX_f0][32];
44 float k_f0[NB_BINS_COMPRESSED_MATRIX_f0][32];
45
45
46 void init_k_f0(){
46 void init_k_f0( void )
47 {
47 unsigned char i;
48 unsigned char i;
48
49
49 for(i=0; i<NB_BINS_COMPRESSED_MATRIX_f0; i++){
50 for(i=0; i<NB_BINS_COMPRESSED_MATRIX_f0; i++){
@@ -87,16 +88,32 void init_k_f0(){
87
88
88 float alpha_M = 45 * (3.1415927/180);
89 float alpha_M = 45 * (3.1415927/180);
89
90
90 void BP1_set(){
91 void BP1_set( float * compressed_spec_mat, unsigned char nb_bins_compressed_spec_mat, unsigned char * LFR_BP1 ){
91 int i, exponent;
92 int i, exponent;
92 float PSDB, PSDE, tmp, NVEC_V0, NVEC_V1, NVEC_V2, aux, tr_SB_SB,
93 float PSDB;
93 e_cross_b_re, e_cross_b_im,
94 float PSDE;
94 n_cross_e_scal_b_re, n_cross_e_scal_b_im,
95 float tmp;
95 ny, nz, bx_bx_star, vphi,
96 float NVEC_V0;
96 significand;
97 float NVEC_V1;
97 signed char nbitexp, nbitsig, expmin, expmax; // 8 bits
98 float NVEC_V2;
99 float aux;
100 float tr_SB_SB;
101 float e_cross_b_re;
102 float e_cross_b_im;
103 float n_cross_e_scal_b_re;
104 float n_cross_e_scal_b_im;
105 float ny;
106 float nz;
107 float bx_bx_star;
108 float vphi;
109 float significand;
110 signed char nbitexp;
111 signed char nbitsig;
112 signed char expmin;
113 signed char expmax; // 8 bits
98 short int rangesig; // 16 bits
114 short int rangesig; // 16 bits
99 unsigned short int psd, tmp_u_short_int; // 16 bits
115 unsigned short int psd;
116 unsigned short int tmp_u_short_int; // 16 bits
100 unsigned short int *pt_u_short_int; // pointer on unsigned 16-bit words
117 unsigned short int *pt_u_short_int; // pointer on unsigned 16-bit words
101 unsigned char tmp_u_char; // 8 bits
118 unsigned char tmp_u_char; // 8 bits
102 unsigned char *pt_u_char; // pointer on unsigned 8-bit bytes
119 unsigned char *pt_u_char; // pointer on unsigned 8-bit bytes
@@ -104,7 +121,7 void BP1_set(){
104 init_k_f0();
121 init_k_f0();
105
122
106 #ifdef DEBUG_TCH
123 #ifdef DEBUG_TCH
107 printf("Number of bins: %d\n", NB_BINS_COMPRESSED_MATRIX_f0);
124 printf("Number of bins: %d\n", nb_bins_compressed_spec_mat);
108 printf("BP1 : \n");
125 printf("BP1 : \n");
109 #endif
126 #endif
110
127
@@ -121,12 +138,12 void BP1_set(){
121 printf("nbitsig : %d, rangesig : %d\n", nbitsig, rangesig);
138 printf("nbitsig : %d, rangesig : %d\n", nbitsig, rangesig);
122 #endif
139 #endif
123
140
124 for(i=0; i<NB_BINS_COMPRESSED_MATRIX_f0; i++){
141 for(i=0; i<nb_bins_compressed_spec_mat; i++){
125 //==============================================
142 //==============================================
126 // BP1 PSDB == PA_LFR_SC_BP1_PB_F0 == 12 bits = 5 bits (exponent) + 7 bits (significand)
143 // BP1 PSDB == PA_LFR_SC_BP1_PB_F0 == 12 bits = 5 bits (exponent) + 7 bits (significand)
127 PSDB = compressed_spectral_matrix_f0[i*30] // S11
144 PSDB = compressed_spec_mat[i*30] // S11
128 + compressed_spectral_matrix_f0[i*30+9] // S22
145 + compressed_spec_mat[i*30+9] // S22
129 + compressed_spectral_matrix_f0[i*30+16]; // S33
146 + compressed_spec_mat[i*30+16]; // S33
130
147
131 significand = frexpf(PSDB/3, &exponent); // 0.5 <= significand < 1
148 significand = frexpf(PSDB/3, &exponent); // 0.5 <= significand < 1
132 // PSDB/3 = significand * 2^exponent
149 // PSDB/3 = significand * 2^exponent
@@ -149,7 +166,7 void BP1_set(){
149 // where just the first nbitsig bits are used (0, ..., 2^nbitsig-1)
166 // where just the first nbitsig bits are used (0, ..., 2^nbitsig-1)
150 tmp_u_short_int = (unsigned short int) (exponent-expmin); // Shift and cast into a 16-bit unsigned int
167 tmp_u_short_int = (unsigned short int) (exponent-expmin); // Shift and cast into a 16-bit unsigned int
151 // where just the first nbitexp bits are used (0, ..., 2^nbitexp-1)
168 // where just the first nbitexp bits are used (0, ..., 2^nbitexp-1)
152 pt_u_short_int = (unsigned short int*) &LFR_BP1_F0[i*9+2]; // Affect an unsigned short int pointer with the
169 pt_u_short_int = (unsigned short int*) &LFR_BP1[i*9+2]; // Affect an unsigned short int pointer with the
153 // adress where the 16-bit word result will be stored
170 // adress where the 16-bit word result will be stored
154 *pt_u_short_int = psd | (tmp_u_short_int << nbitsig); // Put the exponent bits (nbitexp) next to the
171 *pt_u_short_int = psd | (tmp_u_short_int << nbitsig); // Put the exponent bits (nbitexp) next to the
155 // left place of the significand bits (nbitsig), making
172 // left place of the significand bits (nbitsig), making
@@ -161,15 +178,15 void BP1_set(){
161 printf("psd for PSDB significand : %d\n",psd);
178 printf("psd for PSDB significand : %d\n",psd);
162 printf("tmp_u_short_int for PSDB exponent : %d\n",tmp_u_short_int);
179 printf("tmp_u_short_int for PSDB exponent : %d\n",tmp_u_short_int);
163 printf("*pt_u_short_int for PSDB exponent + significand: %.3d or %.4x\n",*pt_u_short_int, *pt_u_short_int);
180 printf("*pt_u_short_int for PSDB exponent + significand: %.3d or %.4x\n",*pt_u_short_int, *pt_u_short_int);
164 printf("LFR_BP1_F0[i*9+3] : %.3d or %.2x\n",LFR_BP1_F0[i*9+3], LFR_BP1_F0[i*9+3]);
181 printf("LFR_BP1[i*9+3] : %.3d or %.2x\n",LFR_BP1[i*9+3], LFR_BP1[i*9+3]);
165 printf("LFR_BP1_F0[i*9+2] : %.3d or %.2x\n",LFR_BP1_F0[i*9+2], LFR_BP1_F0[i*9+2]);
182 printf("LFR_BP1[i*9+2] : %.3d or %.2x\n",LFR_BP1[i*9+2], LFR_BP1[i*9+2]);
166 #endif
183 #endif
167 //==============================================
184 //==============================================
168 // BP1 PSDE == PA_LFR_SC_BP1_PE_F0 == 12 bits = 5 bits (exponent) + 7 bits (significand)
185 // BP1 PSDE == PA_LFR_SC_BP1_PE_F0 == 12 bits = 5 bits (exponent) + 7 bits (significand)
169 PSDE = compressed_spectral_matrix_f0[i*30+21] * k_f0[i][K44_PE] // S44
186 PSDE = compressed_spec_mat[i*30+21] * k_f0[i][K44_PE] // S44
170 + compressed_spectral_matrix_f0[i*30+24] * k_f0[i][K55_PE] // S55
187 + compressed_spec_mat[i*30+24] * k_f0[i][K55_PE] // S55
171 + compressed_spectral_matrix_f0[i*30+22] * k_f0[i][K45_PE_RE] // S45 Re
188 + compressed_spec_mat[i*30+22] * k_f0[i][K45_PE_RE] // S45 Re
172 - compressed_spectral_matrix_f0[i*30+23] * k_f0[i][K45_PE_IM]; // S45 Im
189 - compressed_spec_mat[i*30+23] * k_f0[i][K45_PE_IM]; // S45 Im
173
190
174 significand = frexpf(PSDE/2, &exponent); // 0.5 <= significand < 1
191 significand = frexpf(PSDE/2, &exponent); // 0.5 <= significand < 1
175 // PSDE/2 = significand * 2^exponent
192 // PSDE/2 = significand * 2^exponent
@@ -193,7 +210,7 void BP1_set(){
193 // where just the first nbitsig bits are used (0, ..., 2^nbitsig-1)
210 // where just the first nbitsig bits are used (0, ..., 2^nbitsig-1)
194 tmp_u_short_int = (unsigned short int) (exponent-expmin); // Shift and cast into a 16-bit unsigned int
211 tmp_u_short_int = (unsigned short int) (exponent-expmin); // Shift and cast into a 16-bit unsigned int
195 // where just the first nbitexp bits are used (0, ..., 2^nbitexp-1)
212 // where just the first nbitexp bits are used (0, ..., 2^nbitexp-1)
196 pt_u_short_int = (unsigned short int*) &LFR_BP1_F0[i*9+0]; // Affect an unsigned short int pointer with the
213 pt_u_short_int = (unsigned short int*) &LFR_BP1[i*9+0]; // Affect an unsigned short int pointer with the
197 // adress where the 16-bit word result will be stored
214 // adress where the 16-bit word result will be stored
198 *pt_u_short_int = psd | (tmp_u_short_int << nbitsig); // Put the exponent bits (nbitexp) next to the
215 *pt_u_short_int = psd | (tmp_u_short_int << nbitsig); // Put the exponent bits (nbitexp) next to the
199 // left place of the significand bits (nbitsig), making
216 // left place of the significand bits (nbitsig), making
@@ -205,39 +222,39 void BP1_set(){
205 printf("psd for PSDE significand : %d\n",psd);
222 printf("psd for PSDE significand : %d\n",psd);
206 printf("tmp_u_short_int for PSDE exponent : %d\n",tmp_u_short_int);
223 printf("tmp_u_short_int for PSDE exponent : %d\n",tmp_u_short_int);
207 printf("*pt_u_short_int for PSDE exponent + significand: %.3d or %.4x\n",*pt_u_short_int, *pt_u_short_int);
224 printf("*pt_u_short_int for PSDE exponent + significand: %.3d or %.4x\n",*pt_u_short_int, *pt_u_short_int);
208 printf("LFR_BP1_F0[i*9+1] : %.3d or %.2x\n",LFR_BP1_F0[i*9+1], LFR_BP1_F0[i*9+1]);
225 printf("LFR_BP1[i*9+1] : %.3d or %.2x\n",LFR_BP1[i*9+1], LFR_BP1[i*9+1]);
209 printf("LFR_BP1_F0[i*9+0] : %.3d or %.2x\n",LFR_BP1_F0[i*9+0], LFR_BP1_F0[i*9+0]);
226 printf("LFR_BP1[i*9+0] : %.3d or %.2x\n",LFR_BP1[i*9+0], LFR_BP1[i*9+0]);
210 #endif
227 #endif
211 //==============================================================================
228 //==============================================================================
212 // BP1 normal wave vector == PA_LFR_SC_BP1_NVEC_V0_F0 == 8 bits
229 // BP1 normal wave vector == PA_LFR_SC_BP1_NVEC_V0_F0 == 8 bits
213 // == PA_LFR_SC_BP1_NVEC_V1_F0 == 8 bits
230 // == PA_LFR_SC_BP1_NVEC_V1_F0 == 8 bits
214 // == PA_LFR_SC_BP1_NVEC_V2_F0 == 1 sign bit
231 // == PA_LFR_SC_BP1_NVEC_V2_F0 == 1 sign bit
215 tmp = sqrt( compressed_spectral_matrix_f0[i*30+2] *compressed_spectral_matrix_f0[i*30+2] //Im S12
232 tmp = sqrt( compressed_spec_mat[i*30+2] *compressed_spec_mat[i*30+2] //Im S12
216 +compressed_spectral_matrix_f0[i*30+4] *compressed_spectral_matrix_f0[i*30+4] //Im S13
233 +compressed_spec_mat[i*30+4] *compressed_spec_mat[i*30+4] //Im S13
217 +compressed_spectral_matrix_f0[i*30+11]*compressed_spectral_matrix_f0[i*30+11] //Im S23
234 +compressed_spec_mat[i*30+11]*compressed_spec_mat[i*30+11] //Im S23
218 );
235 );
219 NVEC_V0 = compressed_spectral_matrix_f0[i*30+11]/ tmp; // S23 Im => n1
236 NVEC_V0 = compressed_spec_mat[i*30+11]/ tmp; // S23 Im => n1
220 NVEC_V1 = -compressed_spectral_matrix_f0[i*30+4] / tmp; // S13 Im => n2
237 NVEC_V1 = -compressed_spec_mat[i*30+4] / tmp; // S13 Im => n2
221 NVEC_V2 = compressed_spectral_matrix_f0[i*30+2] / tmp; // S12 Im => n3
238 NVEC_V2 = compressed_spec_mat[i*30+2] / tmp; // S12 Im => n3
222
239
223 LFR_BP1_F0[i*9+4] = (unsigned char) (NVEC_V0*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
240 LFR_BP1[i*9+4] = (unsigned char) (NVEC_V0*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
224 LFR_BP1_F0[i*9+5] = (unsigned char) (NVEC_V1*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
241 LFR_BP1[i*9+5] = (unsigned char) (NVEC_V1*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
225 pt_u_char = (unsigned char*) &NVEC_V2; // affect an unsigned char pointer with the adress of NVEC_V2
242 pt_u_char = (unsigned char*) &NVEC_V2; // affect an unsigned char pointer with the adress of NVEC_V2
226 #ifdef LSB_FIRST_TCH
243 #ifdef LSB_FIRST_TCH
227 LFR_BP1_F0[i*9+6] = pt_u_char[3] & 0x80; // extract the sign bit of NVEC_V2 (32-bit float, sign bit in the 4th octet:PC convention)
244 LFR_BP1[i*9+6] = pt_u_char[3] & 0x80; // extract the sign bit of NVEC_V2 (32-bit float, sign bit in the 4th octet:PC convention)
228 // record it at the 8th bit position (from the right to the left) of LFR_BP1_F0[i*9+6]
245 // record it at the 8th bit position (from the right to the left) of LFR_BP1[i*9+6]
229 #endif
246 #endif
230 #ifdef MSB_FIRST_TCH
247 #ifdef MSB_FIRST_TCH
231 LFR_BP1_F0[i*9+6] = pt_u_char[0] & 0x80; // extract the sign bit of NVEC_V2 (32-bit float, sign bit in the 0th octet:SPARC convention)
248 LFR_BP1[i*9+6] = pt_u_char[0] & 0x80; // extract the sign bit of NVEC_V2 (32-bit float, sign bit in the 0th octet:SPARC convention)
232 // record it at the 8th bit position (from the right to the left) of LFR_BP1_F0[i*9+6]
249 // record it at the 8th bit position (from the right to the left) of LFR_BP1[i*9+6]
233 #endif
250 #endif
234 #ifdef DEBUG_TCH
251 #ifdef DEBUG_TCH
235 printf("NVEC_V0 : %16.8e\n",NVEC_V0);
252 printf("NVEC_V0 : %16.8e\n",NVEC_V0);
236 printf("NVEC_V1 : %16.8e\n",NVEC_V1);
253 printf("NVEC_V1 : %16.8e\n",NVEC_V1);
237 printf("NVEC_V2 : %16.8e\n",NVEC_V2);
254 printf("NVEC_V2 : %16.8e\n",NVEC_V2);
238 printf("LFR_BP1_F0[i*9+4] for NVEC_V0 : %u\n",LFR_BP1_F0[i*9+4]);
255 printf("LFR_BP1[i*9+4] for NVEC_V0 : %u\n",LFR_BP1[i*9+4]);
239 printf("LFR_BP1_F0[i*9+5] for NVEC_V1 : %u\n",LFR_BP1_F0[i*9+5]);
256 printf("LFR_BP1[i*9+5] for NVEC_V1 : %u\n",LFR_BP1[i*9+5]);
240 printf("LFR_BP1_F0[i*9+6] for NVEC_V2 : %u\n",LFR_BP1_F0[i*9+6]);
257 printf("LFR_BP1[i*9+6] for NVEC_V2 : %u\n",LFR_BP1[i*9+6]);
241 #endif
258 #endif
242 //=======================================================
259 //=======================================================
243 // BP1 ellipticity == PA_LFR_SC_BP1_ELLIP_F0 == 4 bits
260 // BP1 ellipticity == PA_LFR_SC_BP1_ELLIP_F0 == 4 bits
@@ -245,81 +262,81 void BP1_set(){
245
262
246 tmp_u_char = (unsigned char) (aux*15 + 0.5); // shift and cast into a 8-bit unsigned char with rounding
263 tmp_u_char = (unsigned char) (aux*15 + 0.5); // shift and cast into a 8-bit unsigned char with rounding
247 // where just the first 4 bits are used (0, ..., 15)
264 // where just the first 4 bits are used (0, ..., 15)
248 LFR_BP1_F0[i*9+6] = LFR_BP1_F0[i*9+6] | (tmp_u_char << 3); // put these 4 bits next to the right place
265 LFR_BP1[i*9+6] = LFR_BP1[i*9+6] | (tmp_u_char << 3); // put these 4 bits next to the right place
249 // of the sign bit of NVEC_V2 (recorded
266 // of the sign bit of NVEC_V2 (recorded
250 // previously in LFR_BP1_F0[i*9+6])
267 // previously in LFR_BP1[i*9+6])
251 #ifdef DEBUG_TCH
268 #ifdef DEBUG_TCH
252 printf("ellipticity : %16.8e\n",aux);
269 printf("ellipticity : %16.8e\n",aux);
253 printf("tmp_u_char for ellipticity : %u\n",tmp_u_char);
270 printf("tmp_u_char for ellipticity : %u\n",tmp_u_char);
254 printf("LFR_BP1_F0[i*9+6] for NVEC_V2 + ellipticity : %u\n",LFR_BP1_F0[i*9+6]);
271 printf("LFR_BP1[i*9+6] for NVEC_V2 + ellipticity : %u\n",LFR_BP1[i*9+6]);
255 #endif
272 #endif
256 //==============================================================
273 //==============================================================
257 // BP1 degree of polarization == PA_LFR_SC_BP1_DOP_F0 == 3 bits
274 // BP1 degree of polarization == PA_LFR_SC_BP1_DOP_F0 == 3 bits
258 tr_SB_SB = compressed_spectral_matrix_f0[i*30] *compressed_spectral_matrix_f0[i*30]
275 tr_SB_SB = compressed_spec_mat[i*30] *compressed_spec_mat[i*30]
259 + compressed_spectral_matrix_f0[i*30+9] *compressed_spectral_matrix_f0[i*30+9]
276 + compressed_spec_mat[i*30+9] *compressed_spec_mat[i*30+9]
260 + compressed_spectral_matrix_f0[i*30+16] *compressed_spectral_matrix_f0[i*30+16]
277 + compressed_spec_mat[i*30+16] *compressed_spec_mat[i*30+16]
261 + 2 * compressed_spectral_matrix_f0[i*30+1] *compressed_spectral_matrix_f0[i*30+1]
278 + 2 * compressed_spec_mat[i*30+1] *compressed_spec_mat[i*30+1]
262 + 2 * compressed_spectral_matrix_f0[i*30+2] *compressed_spectral_matrix_f0[i*30+2]
279 + 2 * compressed_spec_mat[i*30+2] *compressed_spec_mat[i*30+2]
263 + 2 * compressed_spectral_matrix_f0[i*30+3] *compressed_spectral_matrix_f0[i*30+3]
280 + 2 * compressed_spec_mat[i*30+3] *compressed_spec_mat[i*30+3]
264 + 2 * compressed_spectral_matrix_f0[i*30+4] *compressed_spectral_matrix_f0[i*30+4]
281 + 2 * compressed_spec_mat[i*30+4] *compressed_spec_mat[i*30+4]
265 + 2 * compressed_spectral_matrix_f0[i*30+10]*compressed_spectral_matrix_f0[i*30+10]
282 + 2 * compressed_spec_mat[i*30+10]*compressed_spec_mat[i*30+10]
266 + 2 * compressed_spectral_matrix_f0[i*30+11]*compressed_spectral_matrix_f0[i*30+11];
283 + 2 * compressed_spec_mat[i*30+11]*compressed_spec_mat[i*30+11];
267 aux = PSDB*PSDB;
284 aux = PSDB*PSDB;
268 tmp = ( 3*tr_SB_SB - aux ) / ( 2 * aux ); // compute the degree of polarisation
285 tmp = ( 3*tr_SB_SB - aux ) / ( 2 * aux ); // compute the degree of polarisation
269
286
270 tmp_u_char = (unsigned char) (tmp*7 + 0.5);// shift and cast into a 8-bit unsigned char with rounding
287 tmp_u_char = (unsigned char) (tmp*7 + 0.5);// shift and cast into a 8-bit unsigned char with rounding
271 // where just the first 3 bits are used (0, ..., 7)
288 // where just the first 3 bits are used (0, ..., 7)
272 LFR_BP1_F0[i*9+6] = LFR_BP1_F0[i*9+6] | tmp_u_char; // record these 3 bits at the 3 first bit positions
289 LFR_BP1[i*9+6] = LFR_BP1[i*9+6] | tmp_u_char; // record these 3 bits at the 3 first bit positions
273 // (from the right to the left) of LFR_BP1_F0[i*9+6]
290 // (from the right to the left) of LFR_BP1[i*9+6]
274 #ifdef DEBUG_TCH
291 #ifdef DEBUG_TCH
275 printf("DOP : %16.8e\n",tmp);
292 printf("DOP : %16.8e\n",tmp);
276 printf("tmp_u_char for DOP : %u\n",tmp_u_char);
293 printf("tmp_u_char for DOP : %u\n",tmp_u_char);
277 printf("LFR_BP1_F0[i*9+6] for NVEC_V2 + ellipticity + DOP : %u\n",LFR_BP1_F0[i*9+6]);
294 printf("LFR_BP1[i*9+6] for NVEC_V2 + ellipticity + DOP : %u\n",LFR_BP1[i*9+6]);
278 #endif
295 #endif
279 //=======================================================================================
296 //=======================================================================================
280 // BP1 X_SO-component of the Poynting flux == PA_LFR_SC_BP1_SX_F0 == 8 (+ 2) bits
297 // BP1 X_SO-component of the Poynting flux == PA_LFR_SC_BP1_SX_F0 == 8 (+ 2) bits
281 // = 5 bits (exponent) + 3 bits (significand)
298 // = 5 bits (exponent) + 3 bits (significand)
282 // + 1 sign bit + 1 argument bit (two sectors)
299 // + 1 sign bit + 1 argument bit (two sectors)
283 e_cross_b_re = compressed_spectral_matrix_f0[i*30+17]*k_f0[i][K34_SX_RE] //S34 Re
300 e_cross_b_re = compressed_spec_mat[i*30+17]*k_f0[i][K34_SX_RE] //S34 Re
284 + compressed_spectral_matrix_f0[i*30+19]*k_f0[i][K35_SX_RE] //S35 Re
301 + compressed_spec_mat[i*30+19]*k_f0[i][K35_SX_RE] //S35 Re
285 + compressed_spectral_matrix_f0[i*30+5] *k_f0[i][K14_SX_RE] //S14 Re
302 + compressed_spec_mat[i*30+5] *k_f0[i][K14_SX_RE] //S14 Re
286 + compressed_spectral_matrix_f0[i*30+7] *k_f0[i][K15_SX_RE] //S15 Re
303 + compressed_spec_mat[i*30+7] *k_f0[i][K15_SX_RE] //S15 Re
287 + compressed_spectral_matrix_f0[i*30+12]*k_f0[i][K24_SX_RE] //S24 Re
304 + compressed_spec_mat[i*30+12]*k_f0[i][K24_SX_RE] //S24 Re
288 + compressed_spectral_matrix_f0[i*30+14]*k_f0[i][K25_SX_RE] //S25 Re
305 + compressed_spec_mat[i*30+14]*k_f0[i][K25_SX_RE] //S25 Re
289 + compressed_spectral_matrix_f0[i*30+18]*k_f0[i][K34_SX_IM] //S34 Im
306 + compressed_spec_mat[i*30+18]*k_f0[i][K34_SX_IM] //S34 Im
290 + compressed_spectral_matrix_f0[i*30+20]*k_f0[i][K35_SX_IM] //S35 Im
307 + compressed_spec_mat[i*30+20]*k_f0[i][K35_SX_IM] //S35 Im
291 + compressed_spectral_matrix_f0[i*30+6] *k_f0[i][K14_SX_IM] //S14 Im
308 + compressed_spec_mat[i*30+6] *k_f0[i][K14_SX_IM] //S14 Im
292 + compressed_spectral_matrix_f0[i*30+8] *k_f0[i][K15_SX_IM] //S15 Im
309 + compressed_spec_mat[i*30+8] *k_f0[i][K15_SX_IM] //S15 Im
293 + compressed_spectral_matrix_f0[i*30+13]*k_f0[i][K24_SX_IM] //S24 Im
310 + compressed_spec_mat[i*30+13]*k_f0[i][K24_SX_IM] //S24 Im
294 + compressed_spectral_matrix_f0[i*30+15]*k_f0[i][K25_SX_IM]; //S25 Im
311 + compressed_spec_mat[i*30+15]*k_f0[i][K25_SX_IM]; //S25 Im
295 // Im(S_ji) = -Im(S_ij)
312 // Im(S_ji) = -Im(S_ij)
296 // k_ji = k_ij
313 // k_ji = k_ij
297 e_cross_b_im = compressed_spectral_matrix_f0[i*30+17]*k_f0[i][K34_SX_IM] //S34 Re
314 e_cross_b_im = compressed_spec_mat[i*30+17]*k_f0[i][K34_SX_IM] //S34 Re
298 + compressed_spectral_matrix_f0[i*30+19]*k_f0[i][K35_SX_IM] //S35 Re
315 + compressed_spec_mat[i*30+19]*k_f0[i][K35_SX_IM] //S35 Re
299 + compressed_spectral_matrix_f0[i*30+5] *k_f0[i][K14_SX_IM] //S14 Re
316 + compressed_spec_mat[i*30+5] *k_f0[i][K14_SX_IM] //S14 Re
300 + compressed_spectral_matrix_f0[i*30+7] *k_f0[i][K15_SX_IM] //S15 Re
317 + compressed_spec_mat[i*30+7] *k_f0[i][K15_SX_IM] //S15 Re
301 + compressed_spectral_matrix_f0[i*30+12]*k_f0[i][K24_SX_IM] //S24 Re
318 + compressed_spec_mat[i*30+12]*k_f0[i][K24_SX_IM] //S24 Re
302 + compressed_spectral_matrix_f0[i*30+14]*k_f0[i][K25_SX_IM] //S25 Re
319 + compressed_spec_mat[i*30+14]*k_f0[i][K25_SX_IM] //S25 Re
303 - compressed_spectral_matrix_f0[i*30+18]*k_f0[i][K34_SX_RE] //S34 Im
320 - compressed_spec_mat[i*30+18]*k_f0[i][K34_SX_RE] //S34 Im
304 - compressed_spectral_matrix_f0[i*30+20]*k_f0[i][K35_SX_RE] //S35 Im
321 - compressed_spec_mat[i*30+20]*k_f0[i][K35_SX_RE] //S35 Im
305 - compressed_spectral_matrix_f0[i*30+6] *k_f0[i][K14_SX_RE] //S14 Im
322 - compressed_spec_mat[i*30+6] *k_f0[i][K14_SX_RE] //S14 Im
306 - compressed_spectral_matrix_f0[i*30+8] *k_f0[i][K15_SX_RE] //S15 Im
323 - compressed_spec_mat[i*30+8] *k_f0[i][K15_SX_RE] //S15 Im
307 - compressed_spectral_matrix_f0[i*30+13]*k_f0[i][K24_SX_RE] //S24 Im
324 - compressed_spec_mat[i*30+13]*k_f0[i][K24_SX_RE] //S24 Im
308 - compressed_spectral_matrix_f0[i*30+15]*k_f0[i][K25_SX_RE]; //S25 Im
325 - compressed_spec_mat[i*30+15]*k_f0[i][K25_SX_RE]; //S25 Im
309 #ifdef DEBUG_TCH
326 #ifdef DEBUG_TCH
310 printf("ReaSX / 2 : %16.8e\n",e_cross_b_re/2);
327 printf("ReaSX / 2 : %16.8e\n",e_cross_b_re/2);
311 #endif
328 #endif
312 pt_u_char = (unsigned char*) &e_cross_b_re; // Affect an unsigned char pointer with the adress of e_cross_b_re
329 pt_u_char = (unsigned char*) &e_cross_b_re; // Affect an unsigned char pointer with the adress of e_cross_b_re
313 #ifdef LSB_FIRST_TCH
330 #ifdef LSB_FIRST_TCH
314 LFR_BP1_F0[i*9+1] = LFR_BP1_F0[i*9+1] | (pt_u_char[3] & 0x80); // Extract its sign bit (32-bit float, sign bit in the 4th octet:PC convention)
331 LFR_BP1[i*9+1] = LFR_BP1[i*9+1] | (pt_u_char[3] & 0x80); // Extract its sign bit (32-bit float, sign bit in the 4th octet:PC convention)
315 // Record it at the 8th bit position (from the right to the left)
332 // Record it at the 8th bit position (from the right to the left)
316 // of LFR_BP1_F0[i*9+1]
333 // of LFR_BP1[i*9+1]
317 pt_u_char[3] = (pt_u_char[3] & 0x7f); // Make e_cross_b_re be positive in any case: |ReaSX|
334 pt_u_char[3] = (pt_u_char[3] & 0x7f); // Make e_cross_b_re be positive in any case: |ReaSX|
318 #endif
335 #endif
319 #ifdef MSB_FIRST_TCH
336 #ifdef MSB_FIRST_TCH
320 LFR_BP1_F0[i*9+1] = LFR_BP1_F0[i*9+1] | (pt_u_char[0] & 0x80); // Extract its sign bit (32-bit float, sign bit in the 0th octet:SPARC convention)
337 LFR_BP1[i*9+1] = LFR_BP1[i*9+1] | (pt_u_char[0] & 0x80); // Extract its sign bit (32-bit float, sign bit in the 0th octet:SPARC convention)
321 // Record it at the 8th bit position (from the right to the left)
338 // Record it at the 8th bit position (from the right to the left)
322 // of LFR_BP1_F0[i*9+1]
339 // of LFR_BP1[i*9+1]
323 pt_u_char[0] = (pt_u_char[0] & 0x7f); // Make e_cross_b_re be positive in any case: |ReaSX|
340 pt_u_char[0] = (pt_u_char[0] & 0x7f); // Make e_cross_b_re be positive in any case: |ReaSX|
324 #endif
341 #endif
325 significand = frexpf(e_cross_b_re/2, &exponent);// 0.5 <= significand < 1
342 significand = frexpf(e_cross_b_re/2, &exponent);// 0.5 <= significand < 1
@@ -340,7 +357,7 void BP1_set(){
340 significand = 0.5; // min value that can be recorded
357 significand = 0.5; // min value that can be recorded
341 }
358 }
342
359
343 LFR_BP1_F0[i*9+7] = (unsigned char) ((significand*2-1)*7 + 0.5); // Shift and cast into a 8-bit unsigned char with rounding
360 LFR_BP1[i*9+7] = (unsigned char) ((significand*2-1)*7 + 0.5); // Shift and cast into a 8-bit unsigned char with rounding
344 // where just the first 3 bits are used (0, ..., 7)
361 // where just the first 3 bits are used (0, ..., 7)
345 tmp_u_char = (unsigned char) (exponent-expmin); // Shift and cast into a 8-bit unsigned char where
362 tmp_u_char = (unsigned char) (exponent-expmin); // Shift and cast into a 8-bit unsigned char where
346 // just the first 5 bits are used (0, ..., 2^5-1)
363 // just the first 5 bits are used (0, ..., 2^5-1)
@@ -348,14 +365,14 void BP1_set(){
348 printf("|ReaSX| / 2 : %16.8e\n",e_cross_b_re/2);
365 printf("|ReaSX| / 2 : %16.8e\n",e_cross_b_re/2);
349 printf("significand : %16.8e\n",significand);
366 printf("significand : %16.8e\n",significand);
350 printf("exponent : %d\n" ,exponent);
367 printf("exponent : %d\n" ,exponent);
351 printf("LFR_BP1_F0[i*9+7] for ReaSX significand : %u\n",LFR_BP1_F0[i*9+7]);
368 printf("LFR_BP1[i*9+7] for ReaSX significand : %u\n",LFR_BP1[i*9+7]);
352 printf("tmp_u_char for ReaSX exponent : %d\n",tmp_u_char);
369 printf("tmp_u_char for ReaSX exponent : %d\n",tmp_u_char);
353 #endif
370 #endif
354 LFR_BP1_F0[i*9+7] = LFR_BP1_F0[i*9+7] | (tmp_u_char << 3); // shift these 5 bits to the left before logical addition
371 LFR_BP1[i*9+7] = LFR_BP1[i*9+7] | (tmp_u_char << 3); // shift these 5 bits to the left before logical addition
355 // with LFR_BP1_F0[i*9+7]
372 // with LFR_BP1[i*9+7]
356 #ifdef DEBUG_TCH
373 #ifdef DEBUG_TCH
357 printf("LFR_BP1_F0[i*9+7] for ReaSX exponent + significand : %u\n",LFR_BP1_F0[i*9+7]);
374 printf("LFR_BP1[i*9+7] for ReaSX exponent + significand : %u\n",LFR_BP1[i*9+7]);
358 printf("LFR_BP1_F0[i*9+1] for ReaSX sign + PSDE 'exponent' : %u\n",LFR_BP1_F0[i*9+1]);
375 printf("LFR_BP1[i*9+1] for ReaSX sign + PSDE 'exponent' : %u\n",LFR_BP1[i*9+1]);
359 printf("ImaSX / 2 : %16.8e\n",e_cross_b_im/2);
376 printf("ImaSX / 2 : %16.8e\n",e_cross_b_im/2);
360 #endif
377 #endif
361 pt_u_char = (unsigned char*) &e_cross_b_im; // Affect an unsigned char pointer with the adress of e_cross_b_im
378 pt_u_char = (unsigned char*) &e_cross_b_im; // Affect an unsigned char pointer with the adress of e_cross_b_im
@@ -367,12 +384,12 void BP1_set(){
367 #endif
384 #endif
368 tmp_u_char = (e_cross_b_im > e_cross_b_re) ? 0x40 : 0x00; // Determine the sector argument of SX. If |Im| > |Re| affect
385 tmp_u_char = (e_cross_b_im > e_cross_b_re) ? 0x40 : 0x00; // Determine the sector argument of SX. If |Im| > |Re| affect
369 // an unsigned 8-bit char with 01000000; otherwise with null.
386 // an unsigned 8-bit char with 01000000; otherwise with null.
370 LFR_BP1_F0[i*9+1] = LFR_BP1_F0[i*9+1] | tmp_u_char; // Record it as a sign bit at the 7th bit position (from the right
387 LFR_BP1[i*9+1] = LFR_BP1[i*9+1] | tmp_u_char; // Record it as a sign bit at the 7th bit position (from the right
371 // to the left) of LFR_BP1_F0[i*9+1], by simple logical addition.
388 // to the left) of LFR_BP1[i*9+1], by simple logical addition.
372 #ifdef DEBUG_TCH
389 #ifdef DEBUG_TCH
373 printf("|ImaSX| / 2 : %16.8e\n",e_cross_b_im/2);
390 printf("|ImaSX| / 2 : %16.8e\n",e_cross_b_im/2);
374 printf("ArgSX sign : %u\n",tmp_u_char);
391 printf("ArgSX sign : %u\n",tmp_u_char);
375 printf("LFR_BP1_F0[i*9+1] for ReaSX & ArgSX signs + PSDE 'exponent' : %u\n",LFR_BP1_F0[i*9+1]);
392 printf("LFR_BP1[i*9+1] for ReaSX & ArgSX signs + PSDE 'exponent' : %u\n",LFR_BP1[i*9+1]);
376 #endif
393 #endif
377 //======================================================================
394 //======================================================================
378 // BP1 phase velocity estimator == PA_LFR_SC_BP1_VPHI_F0 == 8 (+ 2) bits
395 // BP1 phase velocity estimator == PA_LFR_SC_BP1_VPHI_F0 == 8 (+ 2) bits
@@ -380,44 +397,44 void BP1_set(){
380 // + 1 sign bit + 1 argument bit (two sectors)
397 // + 1 sign bit + 1 argument bit (two sectors)
381 ny = sin(alpha_M)*NVEC_V1 + cos(alpha_M)*NVEC_V2;
398 ny = sin(alpha_M)*NVEC_V1 + cos(alpha_M)*NVEC_V2;
382 nz = NVEC_V0;
399 nz = NVEC_V0;
383 bx_bx_star = cos(alpha_M)*cos(alpha_M)*compressed_spectral_matrix_f0[i*30+9] // S22 Re
400 bx_bx_star = cos(alpha_M)*cos(alpha_M)*compressed_spec_mat[i*30+9] // S22 Re
384 + sin(alpha_M)*sin(alpha_M)*compressed_spectral_matrix_f0[i*30+16] // S33 Re
401 + sin(alpha_M)*sin(alpha_M)*compressed_spec_mat[i*30+16] // S33 Re
385 - 2*sin(alpha_M)*cos(alpha_M)*compressed_spectral_matrix_f0[i*30+10]; // S23 Re
402 - 2*sin(alpha_M)*cos(alpha_M)*compressed_spec_mat[i*30+10]; // S23 Re
386
403
387 n_cross_e_scal_b_re = ny * (compressed_spectral_matrix_f0[i*30+12]*k_f0[i][K24_NY_RE] //S24 Re
404 n_cross_e_scal_b_re = ny * (compressed_spec_mat[i*30+12]*k_f0[i][K24_NY_RE] //S24 Re
388 +compressed_spectral_matrix_f0[i*30+14]*k_f0[i][K25_NY_RE] //S25 Re
405 +compressed_spec_mat[i*30+14]*k_f0[i][K25_NY_RE] //S25 Re
389 +compressed_spectral_matrix_f0[i*30+17]*k_f0[i][K34_NY_RE] //S34 Re
406 +compressed_spec_mat[i*30+17]*k_f0[i][K34_NY_RE] //S34 Re
390 +compressed_spectral_matrix_f0[i*30+19]*k_f0[i][K35_NY_RE] //S35 Re
407 +compressed_spec_mat[i*30+19]*k_f0[i][K35_NY_RE] //S35 Re
391 +compressed_spectral_matrix_f0[i*30+13]*k_f0[i][K24_NY_IM] //S24 Im
408 +compressed_spec_mat[i*30+13]*k_f0[i][K24_NY_IM] //S24 Im
392 +compressed_spectral_matrix_f0[i*30+15]*k_f0[i][K25_NY_IM] //S25 Im
409 +compressed_spec_mat[i*30+15]*k_f0[i][K25_NY_IM] //S25 Im
393 +compressed_spectral_matrix_f0[i*30+18]*k_f0[i][K34_NY_IM] //S34 Im
410 +compressed_spec_mat[i*30+18]*k_f0[i][K34_NY_IM] //S34 Im
394 +compressed_spectral_matrix_f0[i*30+20]*k_f0[i][K35_NY_IM]) //S35 Im
411 +compressed_spec_mat[i*30+20]*k_f0[i][K35_NY_IM]) //S35 Im
395 + nz * (compressed_spectral_matrix_f0[i*30+12]*k_f0[i][K24_NZ_RE] //S24 Re
412 + nz * (compressed_spec_mat[i*30+12]*k_f0[i][K24_NZ_RE] //S24 Re
396 +compressed_spectral_matrix_f0[i*30+14]*k_f0[i][K25_NZ_RE] //S25 Re
413 +compressed_spec_mat[i*30+14]*k_f0[i][K25_NZ_RE] //S25 Re
397 +compressed_spectral_matrix_f0[i*30+17]*k_f0[i][K34_NZ_RE] //S34 Re
414 +compressed_spec_mat[i*30+17]*k_f0[i][K34_NZ_RE] //S34 Re
398 +compressed_spectral_matrix_f0[i*30+19]*k_f0[i][K35_NZ_RE] //S35 Re
415 +compressed_spec_mat[i*30+19]*k_f0[i][K35_NZ_RE] //S35 Re
399 +compressed_spectral_matrix_f0[i*30+13]*k_f0[i][K24_NZ_IM] //S24 Im
416 +compressed_spec_mat[i*30+13]*k_f0[i][K24_NZ_IM] //S24 Im
400 +compressed_spectral_matrix_f0[i*30+15]*k_f0[i][K25_NZ_IM] //S25 Im
417 +compressed_spec_mat[i*30+15]*k_f0[i][K25_NZ_IM] //S25 Im
401 +compressed_spectral_matrix_f0[i*30+18]*k_f0[i][K34_NZ_IM] //S34 Im
418 +compressed_spec_mat[i*30+18]*k_f0[i][K34_NZ_IM] //S34 Im
402 +compressed_spectral_matrix_f0[i*30+20]*k_f0[i][K35_NZ_IM]);//S35 Im
419 +compressed_spec_mat[i*30+20]*k_f0[i][K35_NZ_IM]);//S35 Im
403 // Im(S_ji) = -Im(S_ij)
420 // Im(S_ji) = -Im(S_ij)
404 // k_ji = k_ij
421 // k_ji = k_ij
405 n_cross_e_scal_b_im = ny * (compressed_spectral_matrix_f0[i*30+12]*k_f0[i][K24_NY_IM] //S24 Re
422 n_cross_e_scal_b_im = ny * (compressed_spec_mat[i*30+12]*k_f0[i][K24_NY_IM] //S24 Re
406 +compressed_spectral_matrix_f0[i*30+14]*k_f0[i][K25_NY_IM] //S25 Re
423 +compressed_spec_mat[i*30+14]*k_f0[i][K25_NY_IM] //S25 Re
407 +compressed_spectral_matrix_f0[i*30+17]*k_f0[i][K34_NY_IM] //S34 Re
424 +compressed_spec_mat[i*30+17]*k_f0[i][K34_NY_IM] //S34 Re
408 +compressed_spectral_matrix_f0[i*30+19]*k_f0[i][K35_NY_IM] //S35 Re
425 +compressed_spec_mat[i*30+19]*k_f0[i][K35_NY_IM] //S35 Re
409 -compressed_spectral_matrix_f0[i*30+13]*k_f0[i][K24_NY_RE] //S24 Im
426 -compressed_spec_mat[i*30+13]*k_f0[i][K24_NY_RE] //S24 Im
410 -compressed_spectral_matrix_f0[i*30+15]*k_f0[i][K25_NY_RE] //S25 Im
427 -compressed_spec_mat[i*30+15]*k_f0[i][K25_NY_RE] //S25 Im
411 -compressed_spectral_matrix_f0[i*30+18]*k_f0[i][K34_NY_RE] //S34 Im
428 -compressed_spec_mat[i*30+18]*k_f0[i][K34_NY_RE] //S34 Im
412 -compressed_spectral_matrix_f0[i*30+20]*k_f0[i][K35_NY_RE]) //S35 Im
429 -compressed_spec_mat[i*30+20]*k_f0[i][K35_NY_RE]) //S35 Im
413 + nz * (compressed_spectral_matrix_f0[i*30+12]*k_f0[i][K24_NZ_IM] //S24 Re
430 + nz * (compressed_spec_mat[i*30+12]*k_f0[i][K24_NZ_IM] //S24 Re
414 +compressed_spectral_matrix_f0[i*30+14]*k_f0[i][K25_NZ_IM] //S25 Re
431 +compressed_spec_mat[i*30+14]*k_f0[i][K25_NZ_IM] //S25 Re
415 +compressed_spectral_matrix_f0[i*30+17]*k_f0[i][K34_NZ_IM] //S34 Re
432 +compressed_spec_mat[i*30+17]*k_f0[i][K34_NZ_IM] //S34 Re
416 +compressed_spectral_matrix_f0[i*30+19]*k_f0[i][K35_NZ_IM] //S35 Re
433 +compressed_spec_mat[i*30+19]*k_f0[i][K35_NZ_IM] //S35 Re
417 -compressed_spectral_matrix_f0[i*30+13]*k_f0[i][K24_NZ_RE] //S24 Im
434 -compressed_spec_mat[i*30+13]*k_f0[i][K24_NZ_RE] //S24 Im
418 -compressed_spectral_matrix_f0[i*30+15]*k_f0[i][K25_NZ_RE] //S25 Im
435 -compressed_spec_mat[i*30+15]*k_f0[i][K25_NZ_RE] //S25 Im
419 -compressed_spectral_matrix_f0[i*30+18]*k_f0[i][K34_NZ_RE] //S34 Im
436 -compressed_spec_mat[i*30+18]*k_f0[i][K34_NZ_RE] //S34 Im
420 -compressed_spectral_matrix_f0[i*30+20]*k_f0[i][K35_NZ_RE]);//S35 Im
437 -compressed_spec_mat[i*30+20]*k_f0[i][K35_NZ_RE]);//S35 Im
421 #ifdef DEBUG_TCH
438 #ifdef DEBUG_TCH
422 printf("n_cross_e_scal_b_re : %16.8e\n",n_cross_e_scal_b_re);
439 printf("n_cross_e_scal_b_re : %16.8e\n",n_cross_e_scal_b_re);
423 printf("n_cross_e_scal_b_im : %16.8e\n",n_cross_e_scal_b_im);
440 printf("n_cross_e_scal_b_im : %16.8e\n",n_cross_e_scal_b_im);
@@ -425,15 +442,15 void BP1_set(){
425 // vphi = n_cross_e_scal_b_re / bx_bx_star => sign(VPHI) = sign(n_cross_e_scal_b_re)
442 // vphi = n_cross_e_scal_b_re / bx_bx_star => sign(VPHI) = sign(n_cross_e_scal_b_re)
426 pt_u_char = (unsigned char*) &n_cross_e_scal_b_re; // Affect an unsigned char pointer with the adress of n_cross_e_scal_b_re
443 pt_u_char = (unsigned char*) &n_cross_e_scal_b_re; // Affect an unsigned char pointer with the adress of n_cross_e_scal_b_re
427 #ifdef LSB_FIRST_TCH
444 #ifdef LSB_FIRST_TCH
428 LFR_BP1_F0[i*9+3] = LFR_BP1_F0[i*9+3] | (pt_u_char[3] & 0x80); // Extract its sign bit (32-bit float, sign bit in the 4th octet:PC convention)
445 LFR_BP1[i*9+3] = LFR_BP1[i*9+3] | (pt_u_char[3] & 0x80); // Extract its sign bit (32-bit float, sign bit in the 4th octet:PC convention)
429 // Record it at the 8th bit position (from the right to the left)
446 // Record it at the 8th bit position (from the right to the left)
430 // of LFR_BP1_F0[i*9+3]
447 // of LFR_BP1[i*9+3]
431 pt_u_char[3] = (pt_u_char[3] & 0x7f); // Make n_cross_e_scal_b_re be positive in any case: |n_cross_e_scal_b_re|
448 pt_u_char[3] = (pt_u_char[3] & 0x7f); // Make n_cross_e_scal_b_re be positive in any case: |n_cross_e_scal_b_re|
432 #endif
449 #endif
433 #ifdef MSB_FIRST_TCH
450 #ifdef MSB_FIRST_TCH
434 LFR_BP1_F0[i*9+3] = LFR_BP1_F0[i*9+3] | (pt_u_char[0] & 0x80); // Extract its sign bit (32-bit float, sign bit in the 0th octet:SPARC convention)
451 LFR_BP1[i*9+3] = LFR_BP1[i*9+3] | (pt_u_char[0] & 0x80); // Extract its sign bit (32-bit float, sign bit in the 0th octet:SPARC convention)
435 // Record it at the 8th bit position (from the right to the left)
452 // Record it at the 8th bit position (from the right to the left)
436 // of LFR_BP1_F0[i*9+3]
453 // of LFR_BP1[i*9+3]
437 pt_u_char[0] = (pt_u_char[0] & 0x7f); // Make n_cross_e_scal_b_re be positive in any case: |n_cross_e_scal_b_re|
454 pt_u_char[0] = (pt_u_char[0] & 0x7f); // Make n_cross_e_scal_b_re be positive in any case: |n_cross_e_scal_b_re|
438 #endif
455 #endif
439 vphi = n_cross_e_scal_b_re / bx_bx_star; // Compute |VPHI|
456 vphi = n_cross_e_scal_b_re / bx_bx_star; // Compute |VPHI|
@@ -460,19 +477,19 void BP1_set(){
460 printf("significand : %16.8e\n",significand);
477 printf("significand : %16.8e\n",significand);
461 printf("exponent : %d\n" ,exponent);
478 printf("exponent : %d\n" ,exponent);
462 #endif
479 #endif
463 LFR_BP1_F0[i*9+8] = (unsigned char) ((significand*2-1)*7 + 0.5); // Shift and cast into a 8-bit unsigned char with rounding
480 LFR_BP1[i*9+8] = (unsigned char) ((significand*2-1)*7 + 0.5); // Shift and cast into a 8-bit unsigned char with rounding
464 // where just the first 3 bits are used (0, ..., 7)
481 // where just the first 3 bits are used (0, ..., 7)
465 tmp_u_char = (unsigned char) (exponent-expmin); // Shift and cast into a 8-bit unsigned char where
482 tmp_u_char = (unsigned char) (exponent-expmin); // Shift and cast into a 8-bit unsigned char where
466 // just the first 5 bits are used (0, ..., 2^5-1)
483 // just the first 5 bits are used (0, ..., 2^5-1)
467 #ifdef DEBUG_TCH
484 #ifdef DEBUG_TCH
468 printf("LFR_BP1_F0[i*9+8] for VPHI significand : %u\n",LFR_BP1_F0[i*9+8]);
485 printf("LFR_BP1[i*9+8] for VPHI significand : %u\n",LFR_BP1[i*9+8]);
469 printf("tmp_u_char for VPHI exponent : %d\n",tmp_u_char);
486 printf("tmp_u_char for VPHI exponent : %d\n",tmp_u_char);
470 #endif
487 #endif
471 LFR_BP1_F0[i*9+8] = LFR_BP1_F0[i*9+8] | (tmp_u_char << 3); // shift these 5 bits to the left before logical addition
488 LFR_BP1[i*9+8] = LFR_BP1[i*9+8] | (tmp_u_char << 3); // shift these 5 bits to the left before logical addition
472 // with LFR_BP1_F0[i*9+8]
489 // with LFR_BP1[i*9+8]
473 #ifdef DEBUG_TCH
490 #ifdef DEBUG_TCH
474 printf("LFR_BP1_F0[i*9+8] for VPHI exponent + significand : %u\n",LFR_BP1_F0[i*9+8]);
491 printf("LFR_BP1[i*9+8] for VPHI exponent + significand : %u\n",LFR_BP1[i*9+8]);
475 printf("LFR_BP1_F0[i*9+3] for VPHI sign + PSDB 'exponent' : %u\n",LFR_BP1_F0[i*9+3]);
492 printf("LFR_BP1[i*9+3] for VPHI sign + PSDB 'exponent' : %u\n",LFR_BP1[i*9+3]);
476 #endif
493 #endif
477 pt_u_char = (unsigned char*) &n_cross_e_scal_b_im; // Affect an unsigned char pointer with the adress of n_cross_e_scal_b_im
494 pt_u_char = (unsigned char*) &n_cross_e_scal_b_im; // Affect an unsigned char pointer with the adress of n_cross_e_scal_b_im
478 #ifdef LSB_FIRST_TCH
495 #ifdef LSB_FIRST_TCH
@@ -483,18 +500,19 void BP1_set(){
483 #endif
500 #endif
484 tmp_u_char = (n_cross_e_scal_b_im > n_cross_e_scal_b_re) ? 0x40 : 0x00; // Determine the sector argument of SX. If |Im| > |Re| affect
501 tmp_u_char = (n_cross_e_scal_b_im > n_cross_e_scal_b_re) ? 0x40 : 0x00; // Determine the sector argument of SX. If |Im| > |Re| affect
485 // an unsigned 8-bit char with 01000000; otherwise with null.
502 // an unsigned 8-bit char with 01000000; otherwise with null.
486 LFR_BP1_F0[i*9+3] = LFR_BP1_F0[i*9+3] | tmp_u_char; // Record it as a sign bit at the 7th bit position (from the right
503 LFR_BP1[i*9+3] = LFR_BP1[i*9+3] | tmp_u_char; // Record it as a sign bit at the 7th bit position (from the right
487 // to the left) of LFR_BP1_F0[i*9+3], by simple logical addition.
504 // to the left) of LFR_BP1[i*9+3], by simple logical addition.
488 #ifdef DEBUG_TCH
505 #ifdef DEBUG_TCH
489 printf("|n_cross_e_scal_b_im| : %16.8e\n",n_cross_e_scal_b_im);
506 printf("|n_cross_e_scal_b_im| : %16.8e\n",n_cross_e_scal_b_im);
490 printf("|n_cross_e_scal_b_im|/bx_bx_star/2: %16.8e\n",n_cross_e_scal_b_im/bx_bx_star/2);
507 printf("|n_cross_e_scal_b_im|/bx_bx_star/2: %16.8e\n",n_cross_e_scal_b_im/bx_bx_star/2);
491 printf("ArgNEBX sign : %u\n",tmp_u_char);
508 printf("ArgNEBX sign : %u\n",tmp_u_char);
492 printf("LFR_BP1_F0[i*9+3] for VPHI & ArgNEBX signs + PSDB 'exponent' : %u\n",LFR_BP1_F0[i*9+3]);
509 printf("LFR_BP1[i*9+3] for VPHI & ArgNEBX signs + PSDB 'exponent' : %u\n",LFR_BP1[i*9+3]);
493 #endif
510 #endif
494 }
511 }
495 }
512 }
496
513
497 void BP2_set(){
514 void BP2_set( float * compressed_spec_mat, unsigned char nb_bins_compressed_spec_mat, unsigned char * LFR_BP2 )
515 {
498 int i, exponent;
516 int i, exponent;
499 float aux, significand, cross_re, cross_im;
517 float aux, significand, cross_re, cross_im;
500 signed char nbitexp, nbitsig, expmin, expmax; // 8 bits
518 signed char nbitexp, nbitsig, expmin, expmax; // 8 bits
@@ -503,7 +521,7 void BP2_set(){
503 unsigned short int *pt_u_short_int; // pointer on unsigned 16-bit words
521 unsigned short int *pt_u_short_int; // pointer on unsigned 16-bit words
504
522
505 #ifdef DEBUG_TCH
523 #ifdef DEBUG_TCH
506 printf("Number of bins: %d\n", NB_BINS_COMPRESSED_MATRIX_f0);
524 printf("Number of bins: %d\n", nb_bins_compressed_spec_mat);
507 printf("BP2 : \n");
525 printf("BP2 : \n");
508 #endif
526 #endif
509
527
@@ -519,7 +537,7 void BP2_set(){
519 printf("expmin : %d, expmax : %d\n", expmin, expmax);
537 printf("expmin : %d, expmax : %d\n", expmin, expmax);
520 #endif
538 #endif
521
539
522 for(i = 0; i<NB_BINS_COMPRESSED_MATRIX_f0; i++){
540 for(i = 0; i<nb_bins_compressed_spec_mat; i++){
523 //==============================================
541 //==============================================
524 // BP2 normalized cross correlations == PA_LFR_SC_BP2_CROSS_F0 == 10 * (8+8) bits
542 // BP2 normalized cross correlations == PA_LFR_SC_BP2_CROSS_F0 == 10 * (8+8) bits
525 // == PA_LFR_SC_BP2_CROSS_RE_0_F0 == 8 bits
543 // == PA_LFR_SC_BP2_CROSS_RE_0_F0 == 8 bits
@@ -543,104 +561,104 void BP2_set(){
543 // == PA_LFR_SC_BP2_CROSS_RE_9_F0 == 8 bits
561 // == PA_LFR_SC_BP2_CROSS_RE_9_F0 == 8 bits
544 // == PA_LFR_SC_BP2_CROSS_IM_9_F0 == 8 bits
562 // == PA_LFR_SC_BP2_CROSS_IM_9_F0 == 8 bits
545 // S12
563 // S12
546 aux = sqrt(compressed_spectral_matrix_f0[i*30]*compressed_spectral_matrix_f0[i*30+9]);
564 aux = sqrt(compressed_spec_mat[i*30]*compressed_spec_mat[i*30+9]);
547 cross_re = compressed_spectral_matrix_f0[i*30+1] / aux;
565 cross_re = compressed_spec_mat[i*30+1] / aux;
548 cross_im = compressed_spectral_matrix_f0[i*30+2] / aux;
566 cross_im = compressed_spec_mat[i*30+2] / aux;
549 LFR_BP2_F0[i*30+10] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
567 LFR_BP2[i*30+10] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
550 LFR_BP2_F0[i*30+20] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
568 LFR_BP2[i*30+20] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
551 #ifdef DEBUG_TCH
569 #ifdef DEBUG_TCH
552 printf("LFR_BP2_F0[i*30+10] for cross12_re (%16.8e) : %.3u\n",cross_re, LFR_BP2_F0[i*30+10]);
570 printf("LFR_BP2[i*30+10] for cross12_re (%16.8e) : %.3u\n",cross_re, LFR_BP2[i*30+10]);
553 printf("LFR_BP2_F0[i*30+20] for cross12_im (%16.8e) : %.3u\n",cross_im, LFR_BP2_F0[i*30+20]);
571 printf("LFR_BP2[i*30+20] for cross12_im (%16.8e) : %.3u\n",cross_im, LFR_BP2[i*30+20]);
554 #endif
572 #endif
555 // S13
573 // S13
556 aux = sqrt(compressed_spectral_matrix_f0[i*30]*compressed_spectral_matrix_f0[i*30+16]);
574 aux = sqrt(compressed_spec_mat[i*30]*compressed_spec_mat[i*30+16]);
557 cross_re = compressed_spectral_matrix_f0[i*30+3] / aux;
575 cross_re = compressed_spec_mat[i*30+3] / aux;
558 cross_im = compressed_spectral_matrix_f0[i*30+4] / aux;
576 cross_im = compressed_spec_mat[i*30+4] / aux;
559 LFR_BP2_F0[i*30+11] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
577 LFR_BP2[i*30+11] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
560 LFR_BP2_F0[i*30+21] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
578 LFR_BP2[i*30+21] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
561 #ifdef DEBUG_TCH
579 #ifdef DEBUG_TCH
562 printf("LFR_BP2_F0[i*30+11] for cross13_re (%16.8e) : %.3u\n",cross_re, LFR_BP2_F0[i*30+11]);
580 printf("LFR_BP2[i*30+11] for cross13_re (%16.8e) : %.3u\n",cross_re, LFR_BP2[i*30+11]);
563 printf("LFR_BP2_F0[i*30+21] for cross13_im (%16.8e) : %.3u\n",cross_im, LFR_BP2_F0[i*30+21]);
581 printf("LFR_BP2[i*30+21] for cross13_im (%16.8e) : %.3u\n",cross_im, LFR_BP2[i*30+21]);
564 #endif
582 #endif
565 // S14
583 // S14
566 aux = sqrt(compressed_spectral_matrix_f0[i*30]*compressed_spectral_matrix_f0[i*30+21]);
584 aux = sqrt(compressed_spec_mat[i*30]*compressed_spec_mat[i*30+21]);
567 cross_re = compressed_spectral_matrix_f0[i*30+5] / aux;
585 cross_re = compressed_spec_mat[i*30+5] / aux;
568 cross_im = compressed_spectral_matrix_f0[i*30+6] / aux;
586 cross_im = compressed_spec_mat[i*30+6] / aux;
569 LFR_BP2_F0[i*30+12] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
587 LFR_BP2[i*30+12] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
570 LFR_BP2_F0[i*30+22] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
588 LFR_BP2[i*30+22] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
571 #ifdef DEBUG_TCH
589 #ifdef DEBUG_TCH
572 printf("LFR_BP2_F0[i*30+12] for cross14_re (%16.8e) : %.3u\n",cross_re, LFR_BP2_F0[i*30+12]);
590 printf("LFR_BP2[i*30+12] for cross14_re (%16.8e) : %.3u\n",cross_re, LFR_BP2[i*30+12]);
573 printf("LFR_BP2_F0[i*30+22] for cross14_im (%16.8e) : %.3u\n",cross_im, LFR_BP2_F0[i*30+22]);
591 printf("LFR_BP2[i*30+22] for cross14_im (%16.8e) : %.3u\n",cross_im, LFR_BP2[i*30+22]);
574 #endif
592 #endif
575 // S15
593 // S15
576 aux = sqrt(compressed_spectral_matrix_f0[i*30]*compressed_spectral_matrix_f0[i*30+24]);
594 aux = sqrt(compressed_spec_mat[i*30]*compressed_spec_mat[i*30+24]);
577 cross_re = compressed_spectral_matrix_f0[i*30+7] / aux;
595 cross_re = compressed_spec_mat[i*30+7] / aux;
578 cross_im = compressed_spectral_matrix_f0[i*30+8] / aux;
596 cross_im = compressed_spec_mat[i*30+8] / aux;
579 LFR_BP2_F0[i*30+13] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
597 LFR_BP2[i*30+13] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
580 LFR_BP2_F0[i*30+23] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
598 LFR_BP2[i*30+23] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
581 #ifdef DEBUG_TCH
599 #ifdef DEBUG_TCH
582 printf("LFR_BP2_F0[i*30+13] for cross15_re (%16.8e) : %.3u\n",cross_re, LFR_BP2_F0[i*30+13]);
600 printf("LFR_BP2[i*30+13] for cross15_re (%16.8e) : %.3u\n",cross_re, LFR_BP2[i*30+13]);
583 printf("LFR_BP2_F0[i*30+23] for cross15_im (%16.8e) : %.3u\n",cross_im, LFR_BP2_F0[i*30+23]);
601 printf("LFR_BP2[i*30+23] for cross15_im (%16.8e) : %.3u\n",cross_im, LFR_BP2[i*30+23]);
584 #endif
602 #endif
585 // S23
603 // S23
586 aux = sqrt(compressed_spectral_matrix_f0[i*30+9]*compressed_spectral_matrix_f0[i*30+16]);
604 aux = sqrt(compressed_spec_mat[i*30+9]*compressed_spec_mat[i*30+16]);
587 cross_re = compressed_spectral_matrix_f0[i*30+10] / aux;
605 cross_re = compressed_spec_mat[i*30+10] / aux;
588 cross_im = compressed_spectral_matrix_f0[i*30+11] / aux;
606 cross_im = compressed_spec_mat[i*30+11] / aux;
589 LFR_BP2_F0[i*30+14] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
607 LFR_BP2[i*30+14] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
590 LFR_BP2_F0[i*30+24] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
608 LFR_BP2[i*30+24] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
591 #ifdef DEBUG_TCH
609 #ifdef DEBUG_TCH
592 printf("LFR_BP2_F0[i*30+14] for cross23_re (%16.8e) : %.3u\n",cross_re, LFR_BP2_F0[i*30+14]);
610 printf("LFR_BP2[i*30+14] for cross23_re (%16.8e) : %.3u\n",cross_re, LFR_BP2[i*30+14]);
593 printf("LFR_BP2_F0[i*30+24] for cross23_im (%16.8e) : %.3u\n",cross_im, LFR_BP2_F0[i*30+24]);
611 printf("LFR_BP2[i*30+24] for cross23_im (%16.8e) : %.3u\n",cross_im, LFR_BP2[i*30+24]);
594 #endif
612 #endif
595 // S24
613 // S24
596 aux = sqrt(compressed_spectral_matrix_f0[i*30+9]*compressed_spectral_matrix_f0[i*30+21]);
614 aux = sqrt(compressed_spec_mat[i*30+9]*compressed_spec_mat[i*30+21]);
597 cross_re = compressed_spectral_matrix_f0[i*30+12] / aux;
615 cross_re = compressed_spec_mat[i*30+12] / aux;
598 cross_im = compressed_spectral_matrix_f0[i*30+13] / aux;
616 cross_im = compressed_spec_mat[i*30+13] / aux;
599 LFR_BP2_F0[i*30+15] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
617 LFR_BP2[i*30+15] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
600 LFR_BP2_F0[i*30+25] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
618 LFR_BP2[i*30+25] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
601 #ifdef DEBUG_TCH
619 #ifdef DEBUG_TCH
602 printf("LFR_BP2_F0[i*30+15] for cross24_re (%16.8e) : %.3u\n",cross_re, LFR_BP2_F0[i*30+15]);
620 printf("LFR_BP2[i*30+15] for cross24_re (%16.8e) : %.3u\n",cross_re, LFR_BP2[i*30+15]);
603 printf("LFR_BP2_F0[i*30+25] for cross24_im (%16.8e) : %.3u\n",cross_im, LFR_BP2_F0[i*30+25]);
621 printf("LFR_BP2[i*30+25] for cross24_im (%16.8e) : %.3u\n",cross_im, LFR_BP2[i*30+25]);
604 #endif
622 #endif
605 // S25
623 // S25
606 aux = sqrt(compressed_spectral_matrix_f0[i*30+9]*compressed_spectral_matrix_f0[i*30+24]);
624 aux = sqrt(compressed_spec_mat[i*30+9]*compressed_spec_mat[i*30+24]);
607 cross_re = compressed_spectral_matrix_f0[i*30+14] / aux;
625 cross_re = compressed_spec_mat[i*30+14] / aux;
608 cross_im = compressed_spectral_matrix_f0[i*30+15] / aux;
626 cross_im = compressed_spec_mat[i*30+15] / aux;
609 LFR_BP2_F0[i*30+16] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
627 LFR_BP2[i*30+16] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
610 LFR_BP2_F0[i*30+26] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
628 LFR_BP2[i*30+26] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
611 #ifdef DEBUG_TCH
629 #ifdef DEBUG_TCH
612 printf("LFR_BP2_F0[i*30+16] for cross25_re (%16.8e) : %.3u\n",cross_re, LFR_BP2_F0[i*30+16]);
630 printf("LFR_BP2[i*30+16] for cross25_re (%16.8e) : %.3u\n",cross_re, LFR_BP2[i*30+16]);
613 printf("LFR_BP2_F0[i*30+26] for cross25_im (%16.8e) : %.3u\n",cross_im, LFR_BP2_F0[i*30+26]);
631 printf("LFR_BP2[i*30+26] for cross25_im (%16.8e) : %.3u\n",cross_im, LFR_BP2[i*30+26]);
614 #endif
632 #endif
615 // S34
633 // S34
616 aux = sqrt(compressed_spectral_matrix_f0[i*30+16]*compressed_spectral_matrix_f0[i*30+21]);
634 aux = sqrt(compressed_spec_mat[i*30+16]*compressed_spec_mat[i*30+21]);
617 cross_re = compressed_spectral_matrix_f0[i*30+17] / aux;
635 cross_re = compressed_spec_mat[i*30+17] / aux;
618 cross_im = compressed_spectral_matrix_f0[i*30+18] / aux;
636 cross_im = compressed_spec_mat[i*30+18] / aux;
619 LFR_BP2_F0[i*30+17] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
637 LFR_BP2[i*30+17] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
620 LFR_BP2_F0[i*30+27] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
638 LFR_BP2[i*30+27] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
621 #ifdef DEBUG_TCH
639 #ifdef DEBUG_TCH
622 printf("LFR_BP2_F0[i*30+17] for cross34_re (%16.8e) : %.3u\n",cross_re, LFR_BP2_F0[i*30+17]);
640 printf("LFR_BP2[i*30+17] for cross34_re (%16.8e) : %.3u\n",cross_re, LFR_BP2[i*30+17]);
623 printf("LFR_BP2_F0[i*30+27] for cross34_im (%16.8e) : %.3u\n",cross_im, LFR_BP2_F0[i*30+27]);
641 printf("LFR_BP2[i*30+27] for cross34_im (%16.8e) : %.3u\n",cross_im, LFR_BP2[i*30+27]);
624 #endif
642 #endif
625 // S35
643 // S35
626 aux = sqrt(compressed_spectral_matrix_f0[i*30+16]*compressed_spectral_matrix_f0[i*30+24]);
644 aux = sqrt(compressed_spec_mat[i*30+16]*compressed_spec_mat[i*30+24]);
627 cross_re = compressed_spectral_matrix_f0[i*30+19] / aux;
645 cross_re = compressed_spec_mat[i*30+19] / aux;
628 cross_im = compressed_spectral_matrix_f0[i*30+20] / aux;
646 cross_im = compressed_spec_mat[i*30+20] / aux;
629 LFR_BP2_F0[i*30+18] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
647 LFR_BP2[i*30+18] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
630 LFR_BP2_F0[i*30+28] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
648 LFR_BP2[i*30+28] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
631 #ifdef DEBUG_TCH
649 #ifdef DEBUG_TCH
632 printf("LFR_BP2_F0[i*30+18] for cross35_re (%16.8e) : %.3u\n",cross_re, LFR_BP2_F0[i*30+18]);
650 printf("LFR_BP2[i*30+18] for cross35_re (%16.8e) : %.3u\n",cross_re, LFR_BP2[i*30+18]);
633 printf("LFR_BP2_F0[i*30+28] for cross35_im (%16.8e) : %.3u\n",cross_im, LFR_BP2_F0[i*30+28]);
651 printf("LFR_BP2[i*30+28] for cross35_im (%16.8e) : %.3u\n",cross_im, LFR_BP2[i*30+28]);
634 #endif
652 #endif
635 // S45
653 // S45
636 aux = sqrt(compressed_spectral_matrix_f0[i*30+21]*compressed_spectral_matrix_f0[i*30+24]);
654 aux = sqrt(compressed_spec_mat[i*30+21]*compressed_spec_mat[i*30+24]);
637 cross_re = compressed_spectral_matrix_f0[i*30+22] / aux;
655 cross_re = compressed_spec_mat[i*30+22] / aux;
638 cross_im = compressed_spectral_matrix_f0[i*30+23] / aux;
656 cross_im = compressed_spec_mat[i*30+23] / aux;
639 LFR_BP2_F0[i*30+19] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
657 LFR_BP2[i*30+19] = (unsigned char) (cross_re*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
640 LFR_BP2_F0[i*30+29] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
658 LFR_BP2[i*30+29] = (unsigned char) (cross_im*127.5 + 128); // shift and cast into a 8-bit unsigned char (0, ..., 255) with rounding
641 #ifdef DEBUG_TCH
659 #ifdef DEBUG_TCH
642 printf("LFR_BP2_F0[i*30+19] for cross45_re (%16.8e) : %.3u\n",cross_re, LFR_BP2_F0[i*30+19]);
660 printf("LFR_BP2[i*30+19] for cross45_re (%16.8e) : %.3u\n",cross_re, LFR_BP2[i*30+19]);
643 printf("LFR_BP2_F0[i*30+29] for cross45_im (%16.8e) : %.3u\n",cross_im, LFR_BP2_F0[i*30+29]);
661 printf("LFR_BP2[i*30+29] for cross45_im (%16.8e) : %.3u\n",cross_im, LFR_BP2[i*30+29]);
644 #endif
662 #endif
645 //==============================================
663 //==============================================
646 // BP2 auto correlations == PA_LFR_SC_BP2_AUTO_F0 == 5*16 bits = 5*[6 bits (exponent) + 10 bits (significand)]
664 // BP2 auto correlations == PA_LFR_SC_BP2_AUTO_F0 == 5*16 bits = 5*[6 bits (exponent) + 10 bits (significand)]
@@ -650,10 +668,10 void BP2_set(){
650 // == PA_LFR_SC_BP2_AUTO_A3_F0 == 16 bits
668 // == PA_LFR_SC_BP2_AUTO_A3_F0 == 16 bits
651 // == PA_LFR_SC_BP2_AUTO_A4_F0 == 16 bits
669 // == PA_LFR_SC_BP2_AUTO_A4_F0 == 16 bits
652 // S11
670 // S11
653 significand = frexpf(compressed_spectral_matrix_f0[i*30], &exponent); // 0.5 <= significand < 1
671 significand = frexpf(compressed_spec_mat[i*30], &exponent); // 0.5 <= significand < 1
654 // S11 = significand * 2^exponent
672 // S11 = significand * 2^exponent
655 #ifdef DEBUG_TCH
673 #ifdef DEBUG_TCH
656 printf("S11 : %16.8e\n",compressed_spectral_matrix_f0[i*30]);
674 printf("S11 : %16.8e\n",compressed_spec_mat[i*30]);
657 printf("significand : %16.8e\n",significand);
675 printf("significand : %16.8e\n",significand);
658 printf("exponent : %d\n" ,exponent);
676 printf("exponent : %d\n" ,exponent);
659 #endif
677 #endif
@@ -674,7 +692,7 void BP2_set(){
674 // where just the first nbitsig bits are used (0, ..., 2^nbitsig-1)
692 // where just the first nbitsig bits are used (0, ..., 2^nbitsig-1)
675 tmp_u_short_int = (unsigned short int) (exponent-expmin); // Shift and cast into a 16-bit unsigned int
693 tmp_u_short_int = (unsigned short int) (exponent-expmin); // Shift and cast into a 16-bit unsigned int
676 // where just the first nbitexp bits are used (0, ..., 2^nbitexp-1)
694 // where just the first nbitexp bits are used (0, ..., 2^nbitexp-1)
677 pt_u_short_int = (unsigned short int*) &LFR_BP2_F0[i*30+0]; // Affect an unsigned short int pointer with the
695 pt_u_short_int = (unsigned short int*) &LFR_BP2[i*30+0]; // Affect an unsigned short int pointer with the
678 // adress where the 16-bit word result will be stored
696 // adress where the 16-bit word result will be stored
679 *pt_u_short_int = autocor | (tmp_u_short_int << nbitsig); // Put the exponent bits (nbitexp) next to the
697 *pt_u_short_int = autocor | (tmp_u_short_int << nbitsig); // Put the exponent bits (nbitexp) next to the
680 // left place of the significand bits (nbitsig), making
698 // left place of the significand bits (nbitsig), making
@@ -683,14 +701,14 void BP2_set(){
683 printf("autocor for S11 significand : %u\n",autocor );
701 printf("autocor for S11 significand : %u\n",autocor );
684 printf("tmp_u_char for S11 exponent : %u\n",tmp_u_short_int );
702 printf("tmp_u_char for S11 exponent : %u\n",tmp_u_short_int );
685 printf("*pt_u_short_int for S11 exponent + significand : %.3d or %x\n",*pt_u_short_int, *pt_u_short_int);
703 printf("*pt_u_short_int for S11 exponent + significand : %.3d or %x\n",*pt_u_short_int, *pt_u_short_int);
686 printf("LFR_BP2_F0[i*30+1] : %u or %x\n",LFR_BP2_F0[i*30+1], LFR_BP2_F0[i*30+1]);
704 printf("LFR_BP2[i*30+1] : %u or %x\n",LFR_BP2[i*30+1], LFR_BP2[i*30+1]);
687 printf("LFR_BP2_F0[i*30+0] : %u or %x\n",LFR_BP2_F0[i*30+0], LFR_BP2_F0[i*30+0]);
705 printf("LFR_BP2[i*30+0] : %u or %x\n",LFR_BP2[i*30+0], LFR_BP2[i*30+0]);
688 #endif
706 #endif
689 // S22
707 // S22
690 significand = frexpf(compressed_spectral_matrix_f0[i*30+9], &exponent); // 0.5 <= significand < 1
708 significand = frexpf(compressed_spec_mat[i*30+9], &exponent); // 0.5 <= significand < 1
691 // S22 = significand * 2^exponent
709 // S22 = significand * 2^exponent
692 #ifdef DEBUG_TCH
710 #ifdef DEBUG_TCH
693 printf("S22 : %16.8e\n",compressed_spectral_matrix_f0[i*30+9]);
711 printf("S22 : %16.8e\n",