00001
00002
00003
00004
00005
00076
00077
00078
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128 #include <stdio.h>
00129 #include <stdlib.h>
00130 #include <string.h>
00131 #include <math.h>
00132 #include <complex.h>
00133 #include <omp.h>
00134 #include "interpol_code.h"
00135 #include "polcal.h"
00136 #include "HMIparam.h"
00137 #include "fstats.h"
00138 #include "drms_defs.h"
00139
00140 #undef I //I is the complex number (0,1) in complex.h. We un-define it to avoid confusion with the loop iterative variable i
00141
00142 char *module_name= "HMI_observables2";
00143
00144 #define kRecSetIn "begin" //beginning time for which an output is wanted. MANDATORY PARAMETER.
00145 #define kRecSetIn2 "end" //end time for which an output is wanted. MANDATORY PARAMETER.
00146
00147 #define kTypeSetIn "levin" //data level of the input filtergrams (lev1,lev1d,lev1p) LEV1 BY DEFAULT
00148 #define kTypeSetOut "levout" //data level of the output series (lev1d,lev1p, or lev1.5) LEV1.5 BY DEFAULT
00149 #define WaveLengthIn "wavelength" //filtergram name Ii starting the framelist (i ranges from 0 to 5). MANDATORY PARAMETER.
00150 #define QuickLookIn "quicklook" //quicklook data (yes=1,no=0)? 0 BY DEFAULT
00151 #define CamIDIn "camid" //front camera (camid=1) or side camera (camid=0)?
00152
00153
00154
00155 #define DataCadenceIn "cadence" //cadence
00156 #define SeriesIn "lev1" //series name for the lev1 filtergrams
00157 #define SmoothTables "smooth" //Use smooth look-up tables for the MDI-like algorithm?
00158 #define RotationalFlat "rotational" //force the use of rotational flat fields?
00159
00160 #define minval(x,y) (((x) < (y)) ? (x) : (y))
00161 #define maxval(x,y) (((x) < (y)) ? (y) : (x))
00162
00163
00164 #define LIGHT_SIDE 2 //SIDE CAMERA
00165 #define LIGHT_FRONT 3 //FRONT CAMERA
00166 #define DARK_SIDE 0 //SIDE CAMERA
00167 #define DARK_FRONT 1 //FRONT CAMERA
00168
00169 #define Q_ACS_ECLP 0x2000 //eclipse keyword for the lev1 data
00170
00171
00172
00173
00174 #define QUAL_NODATA (0x80000000) //no l.o.s. observables (Dopplergram, magnetogram, etc...) image was produced (record created, but NO DATA SEGMENT. Most keywords have default value)
00175 #define QUAL_TARGETFILTERGRAMMISSING (0x40000000) //no target filtergram was found near target time (the target filtergram is the filtergram used to identify the framelist): could be due to missing filtergrams or because no observable sequence was run at that time
00176 #define QUAL_NOINTERPOLATEDKEYWORDS (0x20000000) //could not interpolate some keywords at target time (CROTA2, DSUN_OBS, and CRLT_OBS are required by do_interpolate()), because some lev 1 records are missing/corrupted
00177 #define QUAL_NOFRAMELISTINFO (0x10000000) //could not figure out which observables framelist was used, or the framelist run for the required time range is not an observables framelist
00178 #define QUAL_WRONGCADENCE (0x8000000) //the cadence corresponding to the framelist run at the required time does not match the expected value provided by user (could be an error from user, or an issue with the framelist)
00179 #define QUAL_WRONGTARGET (0x4000000) //the target filtergram does not belong to the current framelist (there is something wrong either with the framelist or the target filtergram)
00180 #define QUAL_MISSINGLEV1D (0x2000000) //not enough lev1d filtergrams to produce an observable (too many lev 1 records were missing or corrupted)
00181 #define QUAL_MISSINGKEYWORDLEV1D (0x1000000) //could not read some needed keywords (like FID) in the lev1d data (too many lev 1 records were missing or corrupted and the corresponding lev 1d record is unusable)
00182 #define QUAL_WRONGWAVELENGTHNUM (0x800000) //the number of wavelengths in the lev1d records is not correct (issue with the framelist, or too many lev 1 records were missing or corrupted)
00183 #define QUAL_MISSINGKEYWORDLEV1P (0x400000) //could not read some needed keywords in the lev1p data (too many lev 1 records were missing or corrupted and the corresponding lev 1p record is unusable)
00184 #define QUAL_NOLOOKUPRECORD (0x200000) //could not find a record for the look-up tables for the MDI-like algorithm (the MDI-like algorithm cannot be used)
00185 #define QUAL_NOLOOKUPKEYWORD (0x100000) //could not read the keywords of the look-up tables for the MDI-like algorithm (the MDI-like algorithm cannot be used)
00186 #define QUAL_NOTENOUGHINTERPOLANTS (0x80000) //not enough interpolation points for the temporal interpolation at a given wavelength and polarization (too many lev 1 records were missing or corrupted)
00187 #define QUAL_INTERPOLATIONFAILED (0x40000) //the temporal interpolation routine failed (no lev1d record was produced)
00188 #define QUAL_MISSINGLEV1P (0x20000) //not enough lev1p records to produce an observable (too many lev 1 records were missing or corrupted)
00189 #define QUAL_NOCOEFFKEYWORD (0x200) //could not read the keywords of the polynomial coefficient series
00190 #define QUAL_NOCOEFFPRECORD (0x80) //could not find a record for the polynomial coefficient for the correction of the MDI-like algorithm, or could not access the keywords of a specific record
00191
00192
00193 #define QUAL_LOWINTERPNUM (0x10000) //the number of interpolation points (for temporal interpolation) was not always as high as expected, AND/OR 2 interpolation points were separated by more than the cadence
00194 #define QUAL_LOWKEYWORDNUM (0x8000) //some keywords (especially CROTA2, DSUN_OBS, and CRLT_OBS) could not be interpolated properly at target time, but a closest-neighbour approximation or an extrapolation was used
00195 #define QUAL_ISSTARGET (0x4000) //the ISS loop was OPEN for one or several filtergrams used to produce the observable
00196 #define QUAL_NOTEMP (0x2000) //cannot read the temperatures needed for polarization calibration (default temperature used)
00197 #define QUAL_NOGAPFILL (0x1000) //the code could not properly gap-fill all the lev 1 filtergrams
00198 #define QUAL_LIMBFITISSUE (0x800) //some lev1 records were discarded because R_SUN, and/or CRPIX1/CRPIX2 were missing or too different from the median value of the other lev 1 records (too much jitter for instance)
00199 #define QUAL_NOCOSMICRAY (0x400) //some cosmic-ray hit lists could not be read for the level 1 filtergrams
00200 #define QUAL_ECLIPSE (0x100) //at least one lev1 record was taken during an eclipse
00201
00202
00203
00204
00205
00206
00207
00208
00209
00210
00211
00212
00213 ModuleArgs_t module_args[] =
00214 {
00215 {ARG_STRING, kRecSetIn, "" , "beginning time for which an output is wanted"},
00216 {ARG_STRING, kRecSetIn2, "" , "end time for which an output is wanted"},
00217 {ARG_STRING, kTypeSetIn, "lev1.0", "Level of input series: 1.0,1d,1p"},
00218 {ARG_STRING, kTypeSetOut,"lev1.5", "Level of output series, combination of: 1d,1p,1.5. For example: 1p,1.5"},
00219 {ARG_INT , WaveLengthIn,"3" , "Index (0 to 5) of the wavelength starting the framelist"},
00220 {ARG_INT , QuickLookIn, "0" , "Quicklook data? No=0; Yes=1"},
00221 {ARG_INT , CamIDIn , "1" , "Front (1) or side (0) camera?"},
00222 {ARG_DOUBLE, DataCadenceIn,"45.0" ,"Cadence (in seconds)"},
00223 {ARG_STRING, SeriesIn, "hmi.lev1", "Name of the lev1 series"},
00224 {ARG_INT , SmoothTables, "0", "Use smooth look-up tables? yes=1, no=0 (default)"},
00225 {ARG_INT , RotationalFlat, "0", "Use rotational flat fields? yes=1, no=0 (default)"},
00226 {ARG_END}
00227 };
00228
00229
00230
00231
00232
00233
00234
00235
00236
00237
00238 int needtochangeFID(int HFLID)
00239 {
00240 int need;
00241
00242 switch(HFLID)
00243 {
00244 case 58312: need=1;
00245 break;
00246 default: need=0;
00247 }
00248
00249 return need;
00250 }
00251
00252
00253
00254
00255
00256
00257
00258
00259
00260 int signj(int v)
00261 {
00262 return v > 0 ? 1 : (v < 0 ? -1 : 0);
00263 }
00264
00265
00266
00267
00268
00269
00270
00271
00272
00273
00274
00275
00276
00277 int WhichWavelength(int FID)
00278 {
00279 int result;
00280 int temp;
00281 if(FID < 100000) temp=(FID/10) % 20;
00282 else
00283 {
00284 temp = ((FID-100000)/10) % 20;
00285 }
00286 switch (temp)
00287 {
00288 case 19: result=9;
00289 break;
00290 case 17: result=7;
00291 break;
00292 case 15: result=0;
00293 break;
00294 case 14: result=0;
00295 break;
00296 case 13: result=1;
00297 break;
00298 case 12: result=1;
00299 break;
00300 case 11: result=2;
00301 break;
00302 case 10: result=2;
00303 break;
00304 case 9: result=3;
00305 break;
00306 case 8: result=3;
00307 break;
00308 case 7: result=4;
00309 break;
00310 case 6: result=4;
00311 break;
00312 case 5: result=5;
00313 break;
00314 case 3: result=6;
00315 break;
00316 case 1: result=8;
00317 break;
00318 default: result=-101;
00319 }
00320
00321 return result;
00322
00323 }
00324
00325
00326
00327
00328
00329
00330
00331
00332
00333
00334
00335
00336
00337
00338
00339
00340
00341
00342
00343
00344
00345
00346
00347
00348
00349
00350
00351
00352
00353 int framelistInfo(int HFLID,int HPLTID,int HWLTID,int WavelengthID,int *PHWPLPOS,int *WavelengthIndex,int *WavelengthLocation, int *PPolarizationType,int CamIdIn,int *Pcombine,int *Pnpol,int MaxNumFiltergrams,TIME *PDataCadence,int *CameraValues,int *FID)
00354 {
00355 int framelistSize=0,HFLIDread,FIDread,i,j,compteur;
00356 int PLINDEX,WLINDEX;
00357 FILE *sequencefile;
00358
00359
00360
00361
00362
00363
00364
00365 char *filename=NULL;
00366 char *filename2=NULL;
00367 char *filename3=NULL;
00368
00369 filename=strdup(DEFS_MKPATH("/../Sequences3.txt"));
00370 filename2=strdup(DEFS_MKPATH("/../std.w"));
00371 filename3=strdup(DEFS_MKPATH("/../std.p"));
00372
00373 char line[256];
00374 int PL_Index[MaxNumFiltergrams],WL_Index[MaxNumFiltergrams];
00375 int WT1P,WT2P,WT3P,WT4P,PS1P,PS2P,PS3P;
00376 int found=0;
00377
00378 int combinef,combines,npolf,npols,framelistSizef,framelistSizes,PolarizationTypef,PolarizationTypes,CAMERA;
00379 float DataCadencef,DataCadences;
00380
00381
00382
00383
00384
00385 printf("HFLID OF FRAMELIST = %d\n",HFLID);
00386 sequencefile = fopen(filename,"r");
00387 if(sequencefile == NULL)
00388 {
00389 printf("The file %s does not exist or cannot be read\n",filename);
00390 free(filename);
00391 free(filename2);
00392 free(filename3);
00393 filename=NULL;
00394 filename2=NULL;
00395 filename3=NULL;
00396 return 1;
00397
00398 }
00399
00400 i=0;
00401 compteur=0;
00402 while (fgets(line,256,sequencefile) != NULL)
00403 {
00404 sscanf(line,"%d %d %d %f %f %d %d %d %d %d %d %d %d",&HFLIDread,&PolarizationTypef,&PolarizationTypes,&DataCadencef,&DataCadences,&npolf,&npols,&combinef,&combines,&framelistSizef,&framelistSizes,&j,&CAMERA);
00405 if(HFLIDread == HFLID)
00406 {
00407
00408 if(CamIdIn == LIGHT_FRONT)
00409 {
00410 *Pcombine=combinef;
00411 *Pnpol=npolf;
00412 *PDataCadence=DataCadencef;
00413 *PPolarizationType=PolarizationTypef;
00414 framelistSize=framelistSizef;
00415 if(combinef == 0)
00416 {
00417 if(CAMERA == 3)
00418 {
00419 FID[i]=j;
00420 WavelengthLocation[i]=compteur;
00421 CameraValues[i]=LIGHT_FRONT;
00422 i+=1;
00423 }
00424 }
00425 else
00426 {
00427 FID[i]=j;
00428 WavelengthLocation[i]=compteur;
00429 if(CAMERA == 3) CameraValues[i]=LIGHT_FRONT;
00430 else CameraValues[i]=LIGHT_SIDE;
00431 i+=1;
00432 }
00433 }
00434 if(CamIdIn == LIGHT_SIDE)
00435 {
00436 *Pcombine=combines;
00437 *Pnpol=npols;
00438 *PDataCadence=DataCadences;
00439 *PPolarizationType=PolarizationTypes;
00440 framelistSize=framelistSizes;
00441 if(combines == 0)
00442 {
00443 if(CAMERA == 2)
00444 {
00445 FID[i]=j;
00446 WavelengthLocation[i]=compteur;
00447 CameraValues[i]=LIGHT_SIDE;
00448 i+=1;
00449 }
00450 }
00451 else
00452 {
00453 FID[i]=j;
00454 WavelengthLocation[i]=compteur;
00455 if(CAMERA == 3) CameraValues[i]=LIGHT_FRONT;
00456 else CameraValues[i]=LIGHT_SIDE;
00457 i+=1;
00458 }
00459
00460 }
00461 found=1;
00462 compteur+=1;
00463 }
00464 }
00465 fclose(sequencefile);
00466
00467
00468
00469
00470 int baseindexW,baseindexP;
00471
00472 if((framelistSize == 12) || (framelistSize == 24) || (framelistSize == 36) || (framelistSize == 48) || (framelistSize == 16) || (framelistSize == 32) || (framelistSize == 20) || (framelistSize == 40))
00473 {
00474 switch(WavelengthID)
00475 {
00476 case 9: baseindexW=HWLTID-19;
00477 break;
00478 case 7: baseindexW=HWLTID-17;
00479 break;
00480 case 0: baseindexW=HWLTID-15;
00481 break;
00482 case 1: baseindexW=HWLTID-13;
00483 break;
00484 case 2: baseindexW=HWLTID-11;
00485 break;
00486 case 3: baseindexW=HWLTID-9;
00487 break;
00488 case 4: baseindexW=HWLTID-7;
00489 break;
00490 case 5: baseindexW=HWLTID-5;
00491 break;
00492 case 6: baseindexW=HWLTID-3;
00493 break;
00494 case 8: baseindexW=HWLTID-1;
00495 break;
00496 }
00497 }
00498 if( (framelistSize == 10) )
00499 {
00500 switch(WavelengthID)
00501 {
00502 case 0: baseindexW=HWLTID-14;
00503 break;
00504 case 1: baseindexW=HWLTID-12;
00505 break;
00506 case 2: baseindexW=HWLTID-10;
00507 break;
00508 case 3: baseindexW=HWLTID-8;
00509 break;
00510 case 4: baseindexW=HWLTID-6;
00511 break;
00512 }
00513 }
00514
00515 baseindexP=(HPLTID/10)*10;
00516
00517
00518 for(i=0;i<framelistSize;++i)
00519 {
00520 WLINDEX=(FID[i]/10) % 20;
00521 WL_Index[i]=WLINDEX+baseindexW;
00522 PLINDEX=FID[i]%10;
00523 PL_Index[i]=PLINDEX+baseindexP;
00524 switch(WLINDEX)
00525 {
00526 case 19: WavelengthIndex[i]=9;
00527 break;
00528 case 17: WavelengthIndex[i]=7;
00529 break;
00530 case 15: WavelengthIndex[i]=0;
00531 break;
00532 case 14: WavelengthIndex[i]=0;
00533 break;
00534 case 13: WavelengthIndex[i]=1;
00535 break;
00536 case 12: WavelengthIndex[i]=1;
00537 break;
00538 case 11: WavelengthIndex[i]=2;
00539 break;
00540 case 10: WavelengthIndex[i]=2;
00541 break;
00542 case 9: WavelengthIndex[i]=3;
00543 break;
00544 case 8: WavelengthIndex[i]=3;
00545 break;
00546 case 7: WavelengthIndex[i]=4;
00547 break;
00548 case 6: WavelengthIndex[i]=4;
00549 break;
00550 case 5: WavelengthIndex[i]=5;
00551 break;
00552 case 3: WavelengthIndex[i]=6;
00553 break;
00554 case 1: WavelengthIndex[i]=8;
00555 break;
00556 default: WavelengthIndex[i]=-101;
00557 }
00558 if(WavelengthIndex[i] == -101)
00559 {
00560 printf("Error: WavelengthIndex[i]=-1 \n");
00561 free(filename);
00562 free(filename2);
00563 free(filename3);
00564 filename=NULL;
00565 filename2=NULL;
00566 filename3=NULL;
00567 return 1;
00568
00569 }
00570 }
00571
00572
00573
00574
00575 sequencefile = fopen(filename2,"r");
00576 if(sequencefile == NULL)
00577 {
00578 printf("The file %s does not exist or cannot be read\n",filename2);
00579 free(filename);
00580 free(filename2);
00581 free(filename3);
00582 filename=NULL;
00583 filename2=NULL;
00584 filename3=NULL;
00585 return 1;
00586
00587 }
00588
00589 for(j=0;j<6;++j) fgets(line,256,sequencefile);
00590 while (fgets(line,256,sequencefile) != NULL)
00591 {
00592 sscanf(line,"%d %d %d %d %d",&j,&WT1P,&WT2P,&WT3P,&WT4P);
00593 for(i=0;i<framelistSize;++i)
00594 {
00595 if(j == WL_Index[i])
00596 {
00597 PHWPLPOS[i*7 ]=WT1P;
00598 PHWPLPOS[i*7+1]=WT2P;
00599 PHWPLPOS[i*7+2]=WT3P;
00600 PHWPLPOS[i*7+3]=WT4P;
00601 }
00602 }
00603 }
00604 fclose(sequencefile);
00605
00606 sequencefile = fopen(filename3,"r");
00607 if(sequencefile == NULL)
00608 {
00609 printf("The file %s does not exist or cannot be read\n",filename3);
00610 free(filename);
00611 free(filename2);
00612 free(filename3);
00613 filename=NULL;
00614 filename2=NULL;
00615 filename3=NULL;
00616 return 1;
00617
00618 }
00619
00620 for(j=0;j<6;++j) fgets(line,256,sequencefile);
00621 while (fgets(line,256,sequencefile) != NULL)
00622 {
00623 sscanf(line,"%d %d %d %d",&j,&PS1P,&PS2P,&PS3P);
00624 for(i=0;i<framelistSize;++i)
00625 {
00626 if(j == PL_Index[i])
00627 {
00628 PHWPLPOS[i*7+4]=PS1P;
00629 PHWPLPOS[i*7+5]=PS2P;
00630 PHWPLPOS[i*7+6]=PS3P;
00631 }
00632 }
00633 }
00634 fclose(sequencefile);
00635
00636
00637 if(found == 0)
00638 {
00639 framelistSize=0;
00640 }
00641
00642 free(filename);
00643 free(filename2);
00644 free(filename3);
00645 filename=NULL;
00646 filename2=NULL;
00647 filename3=NULL;
00648 return framelistSize;
00649 }
00650
00651
00652
00653
00654
00655
00656
00657
00658
00659
00660 int StringLocate(const char *s1,const char *s2)
00661 {
00662 int i=0,result=0;
00663 char *ptr;
00664 char temp[2]="a\0";
00665
00666 while(s2[i] != '\0')
00667 {
00668 temp[0]=s2[i];
00669 ptr = strpbrk(s1,temp);
00670 if(ptr == NULL) result=1;
00671 ++i;
00672 }
00673
00674 return result;
00675 }
00676
00677
00678
00679
00680
00681
00682
00683
00684
00685
00686
00687
00688
00689 int MaskCreation(unsigned char *Mask, int nx, int ny, DRMS_Array_t *BadPixels, int HIMGCFID, float *image, DRMS_Array_t *CosmicRays, int nbadperm)
00690 {
00691
00692 int status=2;
00693
00694 if(nx != 4096 || ny != 4096) return status;
00695
00696 char *filename="/home/production/img_cnfg_ids";
00697 const int minid=80;
00698 const int maxid=124;
00699 const int n_tables=1;
00700 const int maxtab=256;
00701
00702 int *badpixellist=BadPixels->data;
00703 int nBadPixels=BadPixels->axis[0];
00704 int *cosmicraylist=NULL;
00705 int ncosmic=0;
00706
00707 if(CosmicRays != NULL)
00708 {
00709 cosmicraylist=CosmicRays->data;
00710 ncosmic=CosmicRays->axis[0];
00711 }
00712 else ncosmic = -1;
00713
00714 int x_orig,y_orig,x_dir,y_dir;
00715 int skip, take, skip_x, skip_y;
00716 int datum, nsx, nss, nlin,nq,nqq,nsy;
00717 int idn=-1;
00718
00719 int *id=NULL, *tab=NULL, *nrows=NULL, *ncols=NULL, *rowstr=NULL, *colstr=NULL, *config=NULL;
00720 id =(int *)(malloc(maxtab*sizeof(int)));
00721 if(id == NULL)
00722 {
00723 printf("Error: memory could not be allocated to id\n");
00724 return 1;
00725
00726 }
00727 tab =(int *)(malloc(maxtab*sizeof(int)));
00728 if(tab == NULL)
00729 {
00730 printf("Error: memory could not be allocated to tab\n");
00731 return 1;
00732
00733 }
00734 nrows =(int *)(malloc(maxtab*sizeof(int)));
00735 if(nrows == NULL)
00736 {
00737 printf("Error: memory could not be allocated to nrows\n");
00738 return 1;
00739
00740 }
00741 ncols =(int *)(malloc(maxtab*sizeof(int)));
00742 if(ncols == NULL)
00743 {
00744 printf("Error: memory could not be allocated to ncols\n");
00745 return 1;
00746
00747 }
00748 rowstr=(int *)(malloc(maxtab*sizeof(int)));
00749 if(rowstr == NULL)
00750 {
00751 printf("Error: memory could not be allocated to rowstr\n");
00752 return 1;
00753
00754 }
00755 colstr=(int *)(malloc(maxtab*sizeof(int)));
00756 if(colstr == NULL)
00757 {
00758 printf("Error: memory could not be allocated to colstr\n");
00759 return 1;
00760
00761 }
00762 config=(int *)(malloc(maxtab*sizeof(int)));
00763 if(config == NULL)
00764 {
00765 printf("Error: memory could not be allocated to config\n");
00766 return 1;
00767
00768 }
00769
00770
00771 int skipt[4*2048];
00772 int taket[4*2048];
00773
00774 int **kx=NULL;
00775 int *kkx=NULL;
00776
00777 kx=(int **)(malloc(9*sizeof(int *)));
00778 if(kx == NULL)
00779 {
00780 printf("Error: memory could not be allocated to kx\n");
00781 return 1;
00782
00783 }
00784
00785
00786 int kx0[11]={4,0,0,1,0,1,1,0,1,1,1};
00787 kx[0]=kx0;
00788 int kx1[7]={2,1,0,1,1,1,2};
00789 kx[2]=kx1;
00790 int kx2[7]={2,0,1,0,0,1,2} ;
00791 kx[4]=kx2;
00792 int kx3[7]={2,0,0,1,0,2,1};
00793 kx[1]=kx3;
00794 int kx4[7]={2,1,1,0,1,2,1};
00795 kx[3]=kx4;
00796 int kx5[7]={1,0,0,2,2};
00797 kx[5]=kx5;
00798 int kx6[5]={1,1,0,2,2};
00799 kx[6]=kx6;
00800 int kx7[5]={1,1,1,2,2};
00801 kx[7]=kx7;
00802 int kx8[5]={1,0,1,2,2};
00803 kx[8]=kx8;
00804
00805 char **filename_table=NULL;
00806 filename_table=(char **)(malloc(n_tables*sizeof(char *)));
00807 if(filename_table == NULL)
00808 {
00809 printf("Error: memory could not be allocated to filename_table\n");
00810 return 1;
00811
00812 }
00813
00814 char *filename_croptable="/home/cvsuser/cvsroot/EGSE/tables/crop/crop6";
00815 filename_table[0]=filename_croptable;
00816
00817 int i,j,k;
00818 int ix, jx;
00819 int n_config;
00820 char string[256];
00821 char strng[5];
00822
00823
00824 FILE *config_file=NULL;
00825 FILE *crop_table=NULL;
00826
00827
00828 config_file=fopen(filename, "r");
00829
00830 if (config_file != NULL){
00831 fgets(string, 256, config_file);
00832 fgets(string, 256, config_file);
00833 fgets(string, 256, config_file);
00834
00835
00836 n_config=0;
00837 fscanf(config_file, "%d", &datum);
00838 id[n_config]=datum;
00839
00840
00841 do {
00842
00843 fscanf(config_file, "%s", string);
00844
00845 config[n_config]=0;
00846 if (string[0] == '4')config[n_config]=0;
00847 if (string[0] == '2') if (string[7] == 'E')config[n_config]=1;
00848 if (string[0] == '2') if (string[7] == 'F')config[n_config]=2;
00849 if (string[0] == '2') if (string[7] == 'G')config[n_config]=3;
00850 if (string[0] == '2') if (string[7] == 'H')config[n_config]=4;
00851 if (string[0] == '1') if (string[7] == 'E')config[n_config]=5;
00852 if (string[0] == '1') if (string[7] == 'F')config[n_config]=6;
00853 if (string[0] == '1') if (string[7] == 'G')config[n_config]=7;
00854 if (string[0] == '1') if (string[7] == 'H')config[n_config]=8;
00855
00856
00857 fscanf(config_file, "%s", string);
00858
00859 fscanf(config_file, "%s", string);
00860
00861 fscanf(config_file, "%s", string);
00862
00863 if (string[0] == 'N') tab[n_config]=-1;
00864 if (string[0] == 'c') tab[n_config]=0;
00865
00866 fscanf(config_file, "%s", string);
00867
00868 fscanf(config_file, "%d", &datum);
00869
00870 fscanf(config_file, "%d", &datum);
00871
00872 fscanf(config_file, "%s", string);
00873
00874 fscanf(config_file, "%d", &datum);
00875 nrows[n_config]=datum;
00876
00877 fscanf(config_file, "%d", &datum);
00878 ncols[n_config]=datum;
00879
00880 fscanf(config_file, "%d", &datum);
00881 rowstr[n_config]=datum;
00882
00883 fscanf(config_file, "%d", &datum);
00884 colstr[n_config]=datum;
00885
00886 fscanf(config_file, "%d", &datum);
00887
00888 fscanf(config_file, "%d", &datum);
00889
00890
00891
00892 ++n_config;
00893 fscanf(config_file, "%d", &datum);
00894 id[n_config]=datum;
00895
00896
00897
00898
00899 }
00900 while (!feof(config_file) && n_config < maxtab);
00901
00902 fclose(config_file);
00903 }
00904
00905
00906
00907
00908 for (i=0; i<n_config; ++i) if (id[i] == HIMGCFID) idn=i;
00909
00910 skip_x=ncols[idn]/2;
00911 skip_y=nrows[idn]/2;
00912
00913
00914
00915 kkx=kx[config[idn]];
00916 nq=kkx[0];
00917 nlin=kkx[2*nq+3-2]*ny/2;
00918 nss=kkx[2*nq+3-1]*nx/2;
00919
00920
00921
00922 if (idn == -1)
00923 {printf("Error: invalid HIMGCFID\n"); status=2;}
00924 else
00925 {
00926 if (tab[idn] != -1)
00927 {
00928 filename_croptable=filename_table[tab[idn]];
00929 crop_table=fopen(filename_croptable, "r");
00930
00931
00932 fscanf(crop_table, "%d", &datum);
00933
00934 fscanf(crop_table, "%d", &datum);
00935 fscanf(crop_table, "%d", &nqq);
00936
00937 for (j=0; j<nqq; ++j)
00938 {
00939 fscanf(crop_table, "%d", &skip);
00940 fscanf(crop_table, "%d", &take);
00941 skipt[j]=skip;
00942 taket[j]=take;
00943 }
00944 fclose(crop_table);
00945 }
00946 else
00947 {
00948 printf("no crop table\n");
00949
00950 for (k=0; k<nq; ++k)
00951 for (j=0; j<nlin; ++j)
00952 {
00953 skipt[k*nss+j]=0;
00954 taket[k*nss+j]=nss;
00955 }
00956 }
00957
00958 nsx=nss-skip_x;
00959 nsy=nlin-skip_y;
00960
00961
00962
00963 for (k=0; k<nq; ++k)
00964 {
00965 x_orig=kkx[k*2+1]*nx - kkx[k*2+1];
00966 y_orig=kkx[k*2+2]*ny - kkx[k*2+2];
00967 x_dir=-(kkx[k*2+1])*2+1;
00968 y_dir=-(kkx[k*2+2])*2+1;
00969
00970
00971 for (j=0; j<skip_y; ++j) for (i=0; i<nss; ++i) Mask[(y_orig+y_dir*j)*nx+x_orig+x_dir*i]=2;
00972 for (j=0; j<nlin; ++j) for (i=0; i<skip_x; ++i) Mask[(y_orig+y_dir*j)*nx+x_orig+x_dir*i]=2;
00973
00974
00975 for (j=0; j<nsy; ++j)
00976 {
00977 jx=j+skip_y;
00978
00979 for (i=0; i<minval(skipt[k*nss+j],nss); ++i){ix=i+skip_x; Mask[(y_orig+y_dir*jx)*nx+x_orig+x_dir*ix]=2;}
00980 for (i=skipt[k*nss+j]; i<minval(skipt[k*nss+j]+taket[k*nss+j],nss); ++i){ix=i+skip_x; Mask[(y_orig+y_dir*jx)*nx+x_orig+x_dir*ix]=0;}
00981 for (i=(skipt[k*nss+j]+taket[k*nss+j]); i<nss; ++i){ix=i+skip_x; Mask[(y_orig+y_dir*jx)*nx+x_orig+x_dir*ix]=2;}
00982 }
00983 }
00984
00985
00986 for(k=0;k<nx*ny;++k)
00987 {
00988 if(Mask[k] == 0)
00989 {
00990 if(isnan(image[k])) Mask[k] = 1;
00991 }
00992 }
00993
00994
00995 if(ncosmic != -1 && nbadperm != -1) nBadPixels = nbadperm;
00996 if(nBadPixels > 0)
00997 {
00998 for (k=0;k<nBadPixels;++k)
00999 {
01000 if(Mask[badpixellist[k]] == 0) Mask[badpixellist[k]] = 1;
01001 }
01002 }
01003
01004
01005
01006 if(ncosmic > 0)
01007 {
01008 for (k=0;k<ncosmic;++k)
01009 {
01010 if(Mask[cosmicraylist[k]] == 0) Mask[cosmicraylist[k]] = 1;
01011 }
01012 }
01013
01014 status=0;
01015
01016 }
01017
01018 free(id);
01019 free(nrows);
01020 free(ncols);
01021 free(rowstr);
01022 free(colstr);
01023 free(config);
01024 free(kx);
01025 free(filename_table);
01026 id=NULL;
01027 nrows=NULL;
01028 ncols=NULL;
01029 rowstr=NULL;
01030 colstr=NULL;
01031 config=NULL;
01032 kx=NULL;
01033 filename_table=NULL;
01034
01035 return status;
01036 }
01037
01038
01039
01040
01041 int heightformation(int FID, double OBSVR, float *CDELT1, float *RSUN, float *CRPIX1, float *CRPIX2, float CROTA2)
01042 {
01043 int wl=0;
01044 int status=0;
01045 float correction=0.0,correction2=0.0;
01046
01047 wl = (FID/10)%20;
01048
01049 if( (wl >= 0) && (wl < 20) )
01050 {
01051 correction = 0.445*exp(-(wl-10.-(float)OBSVR/(0.690/6173.*3.e8/20.)-0.25)*(wl-10.-(float)OBSVR/(0.690/6173.*3.e8/20.)-0.25)/7.1);
01052 correction2 = 0.39*(-2.0*(wl-10.- (float)OBSVR/(0.690/6173.*3.e8/20.)-0.35)/6.15)*exp(-(wl-10.-(float)OBSVR/(0.690/6173.*3.e8/20.)-0.35)*(wl-10.-(float)OBSVR/(0.690/6173.*3.e8/20.)-0.35)/6.15);
01053
01054 *CDELT1 = *CDELT1*(*RSUN)/((*RSUN)-correction);
01055 *RSUN = *RSUN-correction;
01056 *CRPIX1 = *CRPIX1-cos(M_PI-CROTA2*M_PI/180.)*correction2;
01057 *CRPIX2 = *CRPIX2-sin(M_PI-CROTA2*M_PI/180.)*correction2;
01058 }
01059 else status=1;
01060
01061 return status;
01062
01063 }
01064
01065
01066
01067
01068 char *observables_version()
01069 {
01070 return strdup("$Id: HMI_observables2.c,v 1.10 2016/10/03 17:43:52 arta Exp $");
01071 }
01072
01073
01074
01075
01076
01077
01078
01079
01080
01081
01082
01083
01084
01085
01086
01087
01088 int DoIt(void)
01089 {
01090 #define MaxNString 256 //maximum length of strings in character number
01091 double tstart=dsecnd();
01092
01093
01094 char *CODEVERSION =NULL;
01095 CODEVERSION=observables_version();
01096 char *CODEVERSION1=NULL;
01097 CODEVERSION1=interpol_version();
01098 char *CODEVERSION2=NULL;
01099 CODEVERSION2=interpol_version();
01100 char *CODEVERSION3=NULL;
01101 CODEVERSION3=polcal_version();
01102
01103
01104 char *DISTCOEFPATH=NULL;
01105 char *ROTCOEFPATH =NULL;
01106
01107 char HISTORY[MaxNString];
01108 char COMMENT[]="De-rotation: ON; Un-distortion: ON; Re-centering: ON; Re-sizing: OFF; RSUNerr=0.5; correction for cosmic-ray hits; support for rotational flat fields; support for smooth versions of the look-up tables";
01109 struct init_files initfiles;
01110
01111
01112
01113
01114
01115 char *DISTCOEFFILEF=NULL;
01116 char *DISTCOEFFILES=NULL;
01117 char *ROTCOEFFILE=NULL;
01118 DISTCOEFFILEF=strdup(DEFS_MKPATH("/../libs/lev15/distmodel_front_o6_100624.txt"));
01119 DISTCOEFFILES=strdup(DEFS_MKPATH("/../libs/lev15/distmodel_side_o6_100624.txt"));
01120 ROTCOEFFILE =strdup(DEFS_MKPATH("/../libs/lev15/rotcoef_file.txt"));
01121 DISTCOEFPATH =strdup(DEFS_MKPATH("/../libs/lev15/"));
01122 ROTCOEFPATH =strdup(DEFS_MKPATH("/../libs/lev15/"));
01123
01124
01125
01126
01127
01128 initfiles.dist_file_front=DISTCOEFFILEF;
01129 initfiles.dist_file_side =DISTCOEFFILES;
01130 initfiles.diffrot_coef =ROTCOEFFILE;
01131
01132
01133
01134
01135 char *inRecQuery = cmdparams_get_str(&cmdparams, kRecSetIn, NULL);
01136 char *inRecQuery2 = cmdparams_get_str(&cmdparams, kRecSetIn2, NULL);
01137 char *inLev = cmdparams_get_str(&cmdparams, kTypeSetIn, NULL);
01138 char *outLev = cmdparams_get_str(&cmdparams, kTypeSetOut, NULL);
01139 int WavelengthID = cmdparams_get_int(&cmdparams,WaveLengthIn , NULL);
01140 int QuickLook = cmdparams_get_int(&cmdparams,QuickLookIn, NULL);
01141 int CamId = cmdparams_get_int(&cmdparams,CamIDIn, NULL);
01142 TIME DataCadence = cmdparams_get_double(&cmdparams,DataCadenceIn,NULL);
01143 char *inLev1Series = cmdparams_get_str(&cmdparams,SeriesIn, NULL);
01144 int inSmoothTables = cmdparams_get_int(&cmdparams,SmoothTables, NULL);
01145 int inRotationalFlat = cmdparams_get_int(&cmdparams,RotationalFlat, NULL);
01146
01147 if(CamId == 0) CamId = LIGHT_SIDE;
01148 else CamId = LIGHT_FRONT;
01149
01150
01151 if(QuickLook != 0 && QuickLook != 1)
01152 {
01153 printf("The parameter quicklook must be either 0 or 1\n");
01154 return 1;
01155
01156 }
01157
01158 printf("COMMAND LINE PARAMETERS:\n inRecquery = %s \n inRecquery2 = %s \ninLev = %s \n outLev = %s \n WavelengthID = %d \n QuickLook = %d \n CamId = %d \n DataCadence = %f smooth= %d rotational = %d \n",inRecQuery,inRecQuery2,inLev,outLev,WavelengthID,QuickLook,CamId,DataCadence,inSmoothTables,inRotationalFlat);
01159
01160
01161
01162
01163 int NumWavelengths=10;
01164 int MaxNumFiltergrams=72;
01165 int TempIntNum;
01166 int nRecmax = 23040;
01167 char HMISeriesLev1[MaxNString];
01168 char HMISeriesLev10[MaxNString];
01169 char HMISeriesLev1d[MaxNString];
01170 char HMISeriesLev1pa[MaxNString];
01171 char HMISeriesLev1pb[MaxNString];
01172 char HMISeriesLev15a[MaxNString];
01173 char HMISeriesLev15b[MaxNString];
01174 char HMISeriesLev15c[MaxNString];
01175 char HMISeriesLev15d[MaxNString];
01176 char HMISeriesLev15e[MaxNString];
01177 char HMISeriesLookup[MaxNString];
01178 if(inSmoothTables == 1) strcpy(HMISeriesLookup,"su_couvidat.lookup"); else strcpy(HMISeriesLookup,"hmi.lookup");
01179 printf("Series used for the look-up tables: %s\n",HMISeriesLookup);
01180
01181 char CosmicRaySeries[MaxNString]= "hmi.cosmic_rays";
01182 char HMISeriesTemperature[MaxNString]= "hmi.temperature_summary_300s";
01183 char HMISeriesCoeffs[MaxNString]= "hmi.coefficients";
01184 char HMIRotationalFlats[MaxNString]= "su_richard.flatfield_update";
01185
01186 if(QuickLook == 0) TempIntNum=6; else TempIntNum=2;
01187
01188
01189
01190
01191 TIME CadenceRead;
01192 TIME TimeCaution = DataCadence;
01193 TIME MaxSearchDistanceL,MaxSearchDistanceR;
01194 TIME TREC_EPOCH = sscan_time("1993.01.01_00:00:00_TAI");
01195 TIME TREC_STEP = 0.;
01196 TIME TREC_EPOCH0= sscan_time("1993.01.01_00:00:00_TAI");
01197 TIME temptime=0.0, temptime2=0.0;
01198 TIME TimeBegin,TimeEnd,TimeBegin2,TimeEnd2,TargetTime,PreviousTargetTime;
01199 TIME *internTOBS=NULL;
01200 TIME trec;
01201 TIME tobs;
01202 TIME *timeL;
01203
01204 char HMISeries[MaxNString];
01205 char HMILookup[MaxNString];
01206 char HMICoeffs[MaxNString];
01207 char HMISeriesTemp[MaxNString];
01208 char DATEOBS[MaxNString];
01209 char AcceptedLev[16][4] = {"0","d","p","5","dp","pd","d5","5d","dp5","d5p","pd5","p5d","5dp","5pd","p5","5p"};
01210 char timeBegin[MaxNString] ="2000.12.25_00:00:00";
01211 char timeEnd[MaxNString] ="2000.12.25_00:00:00";
01212 char timeBegin2[MaxNString]="2000.12.25_00:00:00";
01213 char timeEnd2[MaxNString] ="2000.12.25_00:00:00";
01214 char **IMGTYPE=NULL;
01215 char CamIds[]="000";
01216 char FSNtemps[]="00000000000000";
01217 char jsocverss[MaxNString];
01218 char **HWLTNSET=NULL;
01219 char TargetISS[]="CLOSED";
01220 char source[64000];
01221 char recnums[MaxNString];
01222 char HMIFlatField0[MaxNString];
01223 char *HMIFlatField;
01224
01225 int *keyL=NULL;
01226 int TestLevIn[3], TestLevOut[15];
01227 int i,temp,temp2,TotalIn,TotalOut,observables,k,ii,iii,j;
01228 int status = DRMS_SUCCESS, status2 = DRMS_SUCCESS, statusA[51], TotalStatus, CreateEmptyRecord=0;
01229 int PolarizationType=-1;
01230 int nRecs1,nRecs1d,nRecs1p,nRecs15;
01231 int ActualTempIntNum;
01232 int framelistSize=0;
01233 int *HWL1POS=NULL;
01234 int *HWL2POS=NULL;
01235 int *HWL3POS=NULL;
01236 int *HWL4POS=NULL;
01237 int *HPL1POS=NULL;
01238 int *HPL2POS=NULL;
01239 int *HPL3POS=NULL;
01240 int *FID=NULL;
01241 int *HFLID=NULL;
01242 int *CFINDEX=NULL;
01243 int *FSN=NULL;
01244 int *HWLTID=NULL;
01245 int *HPLTID=NULL;
01246 int *HIMGCFID=NULL;
01247 int *KeywordMissing=NULL;
01248 int *SegmentRead=NULL;
01249
01250
01251
01252 int *HCAMID=NULL;
01253 int TargetWavelength=0;
01254 int *IndexFiltergram=NULL;
01255 int nIndexFiltergram;
01256 int TargetHFLID=0;
01257 int TargetHWLPOS[4];
01258 int TargetHPLPOS[3];
01259 int TargetHPLTID;
01260 int TargetHWLTID;
01261 int TargetCFINDEX;
01262 int *CARROT=NULL;
01263 int *NBADPERM=NULL;
01264 int axisin[2] ;
01265 int axisout[2];
01266 int Nelem=0;
01267 int PHWPLPOS[MaxNumFiltergrams*7];
01268 int WavelengthIndex[MaxNumFiltergrams], WavelengthLocation[MaxNumFiltergrams], *OrganizeFramelist=NULL, *OrganizeFramelist2=NULL, *FramelistArray=NULL, *SegmentStatus=NULL;
01269 int FIDValues[MaxNumFiltergrams];
01270 int CameraValues[MaxNumFiltergrams];
01271 int FiltergramLocation, Needed;
01272 int Lev1pWanted=0;
01273 int Lev1dWanted=0;
01274 int Lev15Wanted=0;
01275 int Segments1d=0;
01276 int Segments1p=0;
01277 int *ps1=NULL,*ps2=NULL,*ps3=NULL,*fid=NULL,*Wavelengths=NULL;
01278 int method=1;
01279 int npol;
01280 int npolout;
01281 int nSegs1p;
01282 int Lev1pOffset;
01283 int combine;
01284 int ThresholdPol;
01285 int nthreads;
01286 int CARROTint;
01287 int camera,fidfilt;
01288 int ngood;
01289 int MISSVALS[5];
01290 int SATVALS=0;
01291 int WavelengthID2;
01292 int QUALITY=0;
01293 int QUALITYLEV1=0;
01294 int wl=0;
01295 int row,column;
01296 int *FSNL=NULL;
01297 int *QUALITYin=NULL;
01298 int *QUALITYlev1=NULL;
01299 int COSMICCOUNT=0;
01300 int initialrun=0;
01301
01302 DRMS_RecordSet_t *recLev1 = NULL;
01303 DRMS_RecordSet_t *recLev1d = NULL;
01304 DRMS_RecordSet_t *recLev1p = NULL;
01305 DRMS_RecordSet_t *recLev15a= NULL;
01306 DRMS_RecordSet_t *recLev15b= NULL;
01307 DRMS_RecordSet_t *recLev15c= NULL;
01308 DRMS_RecordSet_t *recLev15d= NULL;
01309 DRMS_RecordSet_t *recLev15e= NULL;
01310 DRMS_RecordSet_t *lookup = NULL;
01311 DRMS_RecordSet_t *rectemp = NULL;
01312 DRMS_RecordSet_t *recpoly = NULL;
01313 DRMS_RecordSet_t *recpoly2 = NULL;
01314 DRMS_RecordSet_t *recflat = NULL;
01315 DRMS_RecordSet_t *recflatrot= NULL;
01316
01317 DRMS_Array_t *arrayL0=NULL;
01318 DRMS_Array_t *arrayL1=NULL;
01319 DRMS_Array_t *arrayL2=NULL;
01320 DRMS_Array_t *rotationalflats=NULL;
01321
01322 double diftime=0.0;
01323 double coeff[4],coeff2[4];
01324 double *count=NULL;
01325
01326 float *RSUN=NULL;
01327 float *CROTA2=NULL;
01328 float *CRLTOBS=NULL;
01329 float *X0=NULL;
01330 float *Y0=NULL;
01331 float *CDELT1=NULL;
01332 double *DSUNOBS=NULL;
01333 float *image =NULL;
01334 float *image0=NULL,*image1=NULL,*image2=NULL,*image3=NULL,*image4=NULL;
01335 float **images =NULL;
01336 float **imagesout=NULL;
01337 double *OBSVR=NULL;
01338 double *OBSVW=NULL;
01339 double *OBSVN=NULL;
01340 float *CRLNOBS=NULL;
01341
01342 float TSEL;
01343 float TFRONT;
01344 float X0AVG,Y0AVG,RSUNAVG;
01345 float X0RMS,Y0RMS,RSUNRMS;
01346 float X0AVGF,Y0AVGF,RSUNAVGF,X0AVGS,Y0AVGS,RSUNAVGS;
01347 float CRLNOBSint,CRLTOBSint,CROTA2int,RSUNint,ctime1,ctime2;
01348 double OBSVRint,OBSVWint,OBSVNint,DSUNOBSint;
01349 float cdelt1;
01350 float *X0ARR=NULL, *Y0ARR=NULL, *RSUNARR=NULL;
01351 float RSUNerr=0.5;
01352 float correction,correction2;
01353 float *temparr1=NULL,*temparr2=NULL,*temparr3=NULL,*temparr4=NULL,*LCP=NULL,*RCP=NULL;
01354 float distance;
01355 float obsvr;
01356 float X0LF=0.0,Y0LF=0.0;
01357 float *pztflat;
01358 float *rotflat;
01359
01360
01361 char *FSNS = "FSN";
01362 char *TOBSS = "T_OBS";
01363 char *HWL1POSS = "HWL1POS";
01364 char *HWL2POSS = "HWL2POS";
01365 char *HWL3POSS = "HWL3POS";
01366 char *HWL4POSS = "HWL4POS";
01367 char *HPL1POSS = "HPL1POS";
01368 char *HPL2POSS = "HPL2POS";
01369 char *HPL3POSS = "HPL3POS";
01370 char *FIDS = "FID";
01371 char *HFLIDS = "HFLID";
01372 char *HIMGCFIDS = "HIMGCFID";
01373 char *IMGTYPES = "IMG_TYPE";
01374 char *RSUNS = "R_SUN";
01375 char *CROTA2S = "CROTA2";
01376 char *CRLTOBSS = "CRLT_OBS";
01377 char *DSUNOBSS = "DSUN_OBS";
01378 char *CRPIX1S = "CRPIX1";
01379 char *CRPIX2S = "CRPIX2";
01380 char *X0LFS = "X0_LF";
01381 char *Y0LFS = "Y0_LF";
01382 char *HCAMIDS = "HCAMID";
01383 char *HCFTIDS = "HCFTID";
01384 char *CDELT1S = "CDELT1";
01385 char *CDELT2S = "CDELT2";
01386 char *CSYSER1S = "CSYSER1";
01387 char *CSYSER2S = "CSYSER2";
01388 char *WCSNAMES = "WCSNAME";
01389 char *OBSVRS = "OBS_VR";
01390 char *OBSVWS = "OBS_VW";
01391 char *OBSVNS = "OBS_VN";
01392 char *CRLNOBSS = "CRLN_OBS";
01393 char *CARROTS = "CAR_ROT";
01394
01395 char *TS08 = "T08_PSASM_MEAN";
01396 char *TS01 = "T01_FWMR1_MEAN";
01397 char *TS02 = "T02_FWMR2_MEAN";
01398 char *RSUNOBSS = "RSUN_OBS";
01399 char *HWLTIDS = "HWLTID";
01400 char *HPLTIDS = "HPLTID";
01401 char *WavelengthIDS = "WAVELNID";
01402 char *HWLTNSETS = "HWLTNSET";
01403 char *NBADPERMS = "NBADPERM";
01404 char *COEFF0S = "COEFF0";
01405 char *COEFF1S = "COEFF1";
01406 char *COEFF2S = "COEFF2";
01407 char *COEFF3S = "COEFF3";
01408 char *COUNTS = "COUNT";
01409 char *FLATREC = "FLAT_REC";
01410
01411
01412 char *TRECS = "T_REC";
01413 char *TRECEPOCHS = "T_REC_epoch";
01414 char *TRECSTEPS = "T_REC_step";
01415 char *CADENCES = "CADENCE";
01416 char *DATES = "DATE";
01417 char *DATEOBSS = "DATE__OBS";
01418 char *INSTRUMES = "INSTRUME";
01419 char *CAMERAS = "CAMERA";
01420 char *QUALITYS = "QUALITY";
01421 char *HISTORYS = "HISTORY";
01422 char *COMMENTS = "COMMENT";
01423 char *BLDVERSS = "BLD_VERS";
01424 char *TOTVALSS = "TOTVALS";
01425 char *DATAVALSS = "DATAVALS";
01426 char *MISSVALSS = "MISSVALS";
01427 char *DATAMINS = "DATAMIN2";
01428 char *DATAMAXS = "DATAMAX2";
01429 char *DATAMEDNS = "DATAMED2";
01430 char *DATAMEANS = "DATAMEA2";
01431 char *DATARMSS = "DATARMS2";
01432 char *DATASKEWS = "DATASKE2";
01433 char *DATAKURTS = "DATAKUR2";
01434 char *DATAMINS2 = "DATAMIN";
01435 char *DATAMAXS2 = "DATAMAX";
01436 char *DATAMEDNS2 = "DATAMEDN";
01437 char *DATAMEANS2 = "DATAMEAN";
01438 char *DATARMSS2 = "DATARMS";
01439 char *DATASKEWS2 = "DATASKEW";
01440 char *DATAKURTS2 = "DATAKURT";
01441 char *QLOOKS = "QLOOK";
01442 char *CALFSNS = "CAL_FSN";
01443 char *LUTQUERYS = "LUTQUERY";
01444 char *TSELS = "TSEL";
01445 char *TFRONTS = "TFRONT";
01446 char *TINTNUMS = "TINTNUM";
01447 char *SINTNUMS = "SINTNUM";
01448 char *DISTCOEFS = "DISTCOEF";
01449 char *ROTCOEFS = "ROTCOEF";
01450 char *ODICOEFFS = "ODICOEFF";
01451 char *OROCOEFFS = "OROCOEFF";
01452 char *POLCALMS = "POLCALM";
01453 char *CODEVER0S = "CODEVER0";
01454 char *CODEVER1S = "CODEVER1";
01455 char *CODEVER2S = "CODEVER2";
01456 char *CODEVER3S = "CODEVER3";
01457 char *SATVALSS = "SATVALS";
01458 char *SOURCES = "SOURCE";
01459 char *ierror = NULL;
01460 char **ierrors = NULL;
01461 char *RAWMEDNS = "RAWMEDN";
01462 char *QUALLEV1S = "QUALLEV1";
01463 char TOTVALSSS[80][13] ;
01464 char DATAVALSSS[80][13] ;
01465 char MISSVALSSS[80][13] ;
01466 char DATAMINSS[80][13] ;
01467 char DATAMAXSS[80][13] ;
01468 char DATAMEDNSS[80][13] ;
01469 char DATAMEANSS[80][13] ;
01470 char DATARMSSS[80][13] ;
01471 char DATASKEWSS[80][13] ;
01472 char DATAKURTSS[80][13] ;
01473 char DATAMINS2S[80][13] ;
01474 char DATAMAXS2S[80][13] ;
01475 char DATAMEDNS2S[80][13] ;
01476 char DATAMEANS2S[80][13] ;
01477 char DATARMSS2S[80][13] ;
01478 char DATASKEWS2S[80][13] ;
01479 char DATAKURTS2S[80][13] ;
01480 char *ROTFLAT = "ROT_FLAT";
01481 char query[MaxNString]="QUERY";
01482 char QueryFlatField[MaxNString];
01483 strcpy(QueryFlatField,"");
01484
01485 double minimum,maximum,median,mean,sigma,skewness,kurtosis;
01486
01487 struct initial const_param;
01488 struct keyword *KeyInterp=NULL;
01489 struct keyword KeyInterpOut;
01490 struct parameterDoppler DopplerParameters;
01491 struct polcal_struct pars;
01492
01493 unsigned char *Mask =NULL;
01494
01495 DRMS_Array_t *BadPixels= NULL;
01496 DRMS_Array_t *CosmicRays= NULL;
01497 DRMS_Array_t *arrin[TempIntNum];
01498 DRMS_Array_t *arrerrors[TempIntNum];
01499 DRMS_Array_t **Segments=NULL;
01500 DRMS_Array_t **Ierror=NULL;
01501 DRMS_Array_t **arrLev1d= NULL;
01502 DRMS_Array_t **arrLev1p= NULL;
01503 DRMS_Array_t **arrLev15= NULL;
01504 DRMS_Array_t *arrintable= NULL;
01505 DRMS_Array_t *flatfield=NULL;
01506 DRMS_Array_t *flatfieldrot=NULL;
01507
01508
01509 DRMS_Type_t type1d = DRMS_TYPE_FLOAT;
01510 DRMS_Type_t type1p = DRMS_TYPE_FLOAT;
01511 DRMS_Type_t type15 = DRMS_TYPE_FLOAT;
01512 DRMS_Type_t typeLO = DRMS_TYPE_INT;
01513 DRMS_Type_t typet = DRMS_TYPE_TIME;
01514 DRMS_Type_t typeEr = DRMS_TYPE_CHAR;
01515
01516 DRMS_Segment_t *segin = NULL;
01517 DRMS_Segment_t *segout = NULL;
01518
01519 DRMS_Record_t *rec = NULL;
01520
01521 double t0,t1;
01522 double *keyF=NULL;
01523 double TSTARTFLAT=0.0, TSTOPFLAT=0.0;
01524
01525
01526
01527 char Lev1pSegName[60][5]={"I0","Q0","U0","V0","I1","Q1","U1","V1","I2","Q2","U2","V2","I3","Q3","U3","V3","I4","Q4","U4","V4","I5","Q5","U5","V5","I6","Q6","U6","V6","I7","Q7","U7","V7","I8","Q8","U8","V8","I9","Q9","U9","V9","LCP0","RCP0","LCP1","RCP1","LCP2","RCP2","LCP3","RCP3","LCP4","RCP4","LCP5","RCP5","LCP6","RCP6","LCP7","RCP7","LCP8","RCP8","LCP9","RCP9"};
01528
01529
01530
01531
01532
01533
01534
01535 nthreads=omp_get_max_threads();
01536 printf("NUMBER OF THREADS USED BY OPEN MP= %d\n",nthreads);
01537
01538
01539
01540
01541 TotalIn=0;
01542 for(i=0;i<=2;++i)
01543 {
01544 temp = StringLocate(inLev,AcceptedLev[i]);
01545
01546 if(temp == 0) TestLevIn[i]=1; else TestLevIn[i]=0;
01547 TotalIn+=TestLevIn[i];
01548 }
01549
01550 TotalOut=0;
01551 for(i=0;i<=14;++i)
01552 {
01553 temp = StringLocate(outLev,AcceptedLev[i+1]);
01554 if(temp == 0) TestLevOut[i]=1; else TestLevOut[i]=0;
01555 TotalOut+=TestLevOut[i];
01556 }
01557
01558
01559
01560
01561
01562 if(TotalIn != 1)
01563 {
01564 if(StringLocate(inLev,"1") == 0)
01565 {
01566 TestLevIn[0]=1;
01567 TotalIn=1;
01568 }
01569 else
01570 {
01571 printf("The parameter levin must be one of the following strings (select only one): lev1, lev1d, or lev1p\n");
01572 return 1;
01573
01574 }
01575 }
01576
01577 if(TotalOut == 0)
01578 {
01579 printf("The parameter levout must be a combination of the following strings: lev1d, lev1p, or lev1.5\n");
01580 return 1;
01581 }
01582
01583
01584 if(TestLevIn[1]==1 && TestLevOut[13]==0 && TestLevOut[14]==0 && TestLevOut[1]==0 && TestLevOut[2]==0)
01585 {
01586 printf("The parameter levin is the string lev1d, therefore levout can only be a combination of the strings lev1p and lev1.5\n");
01587 return 1;
01588 }
01589 if(TestLevIn[2]==1 && TestLevOut[2]==0)
01590 {
01591 printf("The parameter levin is the string lev1p, therefore levout can only be the string lev1.5\n");
01592 return 1;
01593 }
01594
01595
01596
01597
01598
01599
01600
01601
01602
01603 if(WavelengthID > NumWavelengths-1 || WavelengthID < 0)
01604 {
01605 printf("The parameter WaveLengthIn is not in the range 0-9\n");
01606 return 1;
01607 }
01608
01609
01610
01611 if(DataCadence != 22.5 && DataCadence != 45.0 && DataCadence != 60.0 && DataCadence != 90.0 && DataCadence != 135.0 && DataCadence != 75.0 && DataCadence != 150.0 && DataCadence != 720.0)
01612 {
01613 printf("The command-line parameter cadence is not an accepted value\n");
01614 return 1;
01615 }
01616
01617
01618
01619
01620
01621
01622
01623 if(TestLevIn[0]==1) printf("Input data are level 1.0\n");
01624 if(TestLevIn[1]==1) printf("Input data are level 1d\n");
01625 if(TestLevIn[2]==1) printf("Input data are level 1p\n");
01626 if(TestLevOut[0]==1 || TestLevOut[3]==1 || TestLevOut[4]==1 || TestLevOut[5]==1 || TestLevOut[6]==1 || TestLevOut[7]==1 || TestLevOut[8]==1 || TestLevOut[9]==1 || TestLevOut[10]==1 || TestLevOut[11]==1 || TestLevOut[12]==1)
01627 {
01628 printf("Output data are level 1d\n");
01629 Lev1dWanted=1;
01630 }
01631 if(TestLevOut[1]==1 || TestLevOut[3]==1 || TestLevOut[4]==1 || TestLevOut[7]==1 || TestLevOut[8]==1 || TestLevOut[9]==1 || TestLevOut[10]==1 || TestLevOut[11]==1 || TestLevOut[12]==1 || TestLevOut[13]==1 || TestLevOut[14]==1)
01632 {
01633 printf("Output data are level 1p\n");
01634 Lev1pWanted=1;
01635 }
01636 if(TestLevOut[2]==1 || TestLevOut[5]==1 || TestLevOut[6]==1 || TestLevOut[7]==1 || TestLevOut[8]==1 || TestLevOut[9]==1 || TestLevOut[10]==1 || TestLevOut[11]==1 || TestLevOut[12]==1 || TestLevOut[13]==1 || TestLevOut[14]==1)
01637 {
01638 printf("Output data are level 1.5\n");
01639 Lev15Wanted=1;
01640 }
01641
01642
01643
01644
01645
01646
01647
01648
01649
01650
01651
01652
01653
01654
01655 printf("BEGINNING TIME: %s\n",inRecQuery);
01656 printf("ENDING TIME: %s\n",inRecQuery2);
01657
01658
01659
01660
01661 TimeBegin=sscan_time(inRecQuery);
01662 TimeEnd =sscan_time(inRecQuery2);
01663
01664
01665 if(TimeBegin > TimeEnd)
01666 {
01667 printf("Error: the ending time must be later than the beginning time!\n");
01668 return 1;
01669 }
01670
01671
01672
01673
01674
01675 if(TestLevIn[0]==1)
01676 {
01677 status = initialize_interpol(&const_param,&initfiles,4096,4096);
01678 if(status != 0)
01679 {
01680 printf("Error: could not initialize the gapfilling, derotation, and temporal interpolation routines\n");
01681 return 1;
01682 }
01683
01684
01685 }
01686 if(Lev1pWanted || (Lev15Wanted && TestLevIn[2]==0))
01687 {
01688
01689
01690 status = init_polcal(&pars, method);
01691 if(status != 0)
01692 {
01693 printf("Error: could not initialize the polarization calibration routine\n");
01694 return 1;
01695 }
01696 }
01697
01698
01699
01700 strcpy(HMISeriesLev1,inLev1Series);
01701
01702
01703 strcpy(HMISeriesLev10,HMISeriesLev1);
01704
01705
01706
01707
01708 if( QuickLook == 1)
01709 {
01710 if(DataCadence == 22.5) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d22Q");
01711 if(DataCadence == 45.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d45Q");
01712 if(DataCadence == 60.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d60Q");
01713 if(DataCadence == 75.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d75Q");
01714 if(DataCadence == 90.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d90Q");
01715 if(DataCadence == 120.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d120Q");
01716 if(DataCadence == 135.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d135Q");
01717 if(DataCadence == 150.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d150Q");
01718 }
01719 else
01720 {
01721 if(DataCadence == 22.5) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d22");
01722 if(DataCadence == 45.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d45");
01723 if(DataCadence == 60.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d60");
01724 if(DataCadence == 75.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d75");
01725 if(DataCadence == 90.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d90");
01726 if(DataCadence == 120.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d120");
01727 if(DataCadence == 135.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d135");
01728 if(DataCadence == 150.0) strcpy(HMISeriesLev1d,"hmi.HMISeriesLev1d150");
01729 }
01730
01731
01732
01733
01734 if(QuickLook == 1)
01735 {
01736 if(DataCadence == 45.0)
01737 {
01738 strcpy(HMISeriesLev15a,"hmi_test.V2_45s_nrt" );
01739 strcpy(HMISeriesLev15b,"hmi_test.M2_45s_nrt" );
01740 strcpy(HMISeriesLev15c,"hmi_test.Ld2_45s_nrt" );
01741 strcpy(HMISeriesLev15d,"hmi_test.Lw2_45s_nrt" );
01742 strcpy(HMISeriesLev15e,"hmi_test.Ic2_45s_nrt" );
01743 }
01744 if(DataCadence == 22.5)
01745 {
01746 strcpy(HMISeriesLev15a,"hmi_test.V_22s_nrt" );
01747 strcpy(HMISeriesLev15b,"hmi_test.M_22s_nrt" );
01748 strcpy(HMISeriesLev15c,"hmi_test.Ld_22s_nrt" );
01749 strcpy(HMISeriesLev15d,"hmi_test.Lw_22s_nrt" );
01750 strcpy(HMISeriesLev15e,"hmi_test.Ic_22s_nrt" );
01751 }
01752 if(DataCadence == 60.0)
01753 {
01754 strcpy(HMISeriesLev15a,"hmi.V_60s_nrt");
01755 strcpy(HMISeriesLev15b,"hmi.M_60s_nrt");
01756 strcpy(HMISeriesLev15c,"hmi.Ld_60s_nrt");
01757 strcpy(HMISeriesLev15d,"hmi.Lw_60s_nrt");
01758 strcpy(HMISeriesLev15e,"hmi.Ic_60s_nrt");
01759 }
01760 if(DataCadence == 75.0)
01761 {
01762 strcpy(HMISeriesLev15a,"hmi.V_75s_nrt" );
01763 strcpy(HMISeriesLev15b,"hmi.M_75s_nrt" );
01764 strcpy(HMISeriesLev15c,"hmi.Ld_75s_nrt" );
01765 strcpy(HMISeriesLev15d,"hmi.Lw_75s_nrt" );
01766 strcpy(HMISeriesLev15e,"hmi.Ic_75s_nrt" );
01767 }
01768 if(DataCadence == 720.0)
01769 {
01770 strcpy(HMISeriesLev15a,"hmi.V_720s_nrt");
01771 strcpy(HMISeriesLev15b,"hmi.M_720s_nrt");
01772 strcpy(HMISeriesLev15c,"hmi.Ld_720s_nrt");
01773 strcpy(HMISeriesLev15d,"hmi.Lw_720s_nrt");
01774 strcpy(HMISeriesLev15e,"hmi.Ic_720s_nrt");
01775 }
01776 }
01777 else
01778 {
01779 if(DataCadence == 45.0)
01780 {
01781 strcpy(HMISeriesLev15a,"hmi_test.V2_45s" );
01782 strcpy(HMISeriesLev15b,"hmi_test.M2_45s" );
01783 strcpy(HMISeriesLev15c,"hmi_test.Ld2_45s" );
01784 strcpy(HMISeriesLev15d,"hmi_test.Lw2_45s" );
01785 strcpy(HMISeriesLev15e,"hmi_test.Ic2_45s" );
01786 }
01787 if(DataCadence == 22.5)
01788 {
01789 strcpy(HMISeriesLev15a,"hmi_test.V_22s" );
01790 strcpy(HMISeriesLev15b,"hmi_test.M_22s" );
01791 strcpy(HMISeriesLev15c,"hmi_test.Ld_22s" );
01792 strcpy(HMISeriesLev15d,"hmi_test.Lw_22s" );
01793 strcpy(HMISeriesLev15e,"hmi_test.Ic_22s" );
01794 }
01795 if(DataCadence == 60.0)
01796 {
01797 strcpy(HMISeriesLev15a,"hmi.V_60s");
01798 strcpy(HMISeriesLev15b,"hmi.M_60s");
01799 strcpy(HMISeriesLev15c,"hmi.Ld_60s");
01800 strcpy(HMISeriesLev15d,"hmi.Lw_60s");
01801 strcpy(HMISeriesLev15e,"hmi.Ic_60s");
01802 }
01803 if(DataCadence == 75.0)
01804 {
01805 strcpy(HMISeriesLev15a,"hmi.V_75s" );
01806 strcpy(HMISeriesLev15b,"hmi.M_75s" );
01807 strcpy(HMISeriesLev15c,"hmi.Ld_75s" );
01808 strcpy(HMISeriesLev15d,"hmi.Lw_75s" );
01809 strcpy(HMISeriesLev15e,"hmi.Ic_75s" );
01810 }
01811 if(DataCadence == 720.0)
01812 {
01813 strcpy(HMISeriesLev15a,"hmi.V_720s" );
01814 strcpy(HMISeriesLev15b,"hmi.M_720s" );
01815 strcpy(HMISeriesLev15c,"hmi.Ld_720s" );
01816 strcpy(HMISeriesLev15d,"hmi.Lw_720s" );
01817 strcpy(HMISeriesLev15e,"hmi.Ic_720s" );
01818 }
01819 }
01820
01821
01822
01823
01824
01825 if(QuickLook == 1)
01826 {
01827 if(DataCadence == 22.5) strcpy(HMISeriesLev1pb,"hmi.HMISeriesLev1pb22Q");
01828 if(DataCadence == 45.0) strcpy(HMISeriesLev1pb,"hmi.HMISeriesLev1pb45Q");
01829 if(DataCadence == 60.0) strcpy(HMISeriesLev1pb,"hmi.HMISeriesLev1pb60Q");
01830 if(DataCadence == 75.0) strcpy(HMISeriesLev1pb,"hmi.HMISeriesLev1pb75Q");
01831 }
01832 else
01833 {
01834 if(DataCadence == 22.5) strcpy(HMISeriesLev1pb,"hmi.HMISeriesLev1pb22");
01835 if(DataCadence == 45.0) strcpy(HMISeriesLev1pb,"hmi.HMISeriesLev1pb45");
01836 if(DataCadence == 60.0) strcpy(HMISeriesLev1pb,"hmi.HMISeriesLev1pb60");
01837 if(DataCadence == 75.0) strcpy(HMISeriesLev1pb,"hmi.HMISeriesLev1pb75");
01838 }
01839
01840
01841
01842
01843 if(QuickLook == 1)
01844 {
01845 if( DataCadence == 45.0) strcpy(HMISeriesLev1pa,"hmi.HMISeriesLev1pa45Q");
01846 if( DataCadence == 90.0) strcpy(HMISeriesLev1pa,"hmi.HMISeriesLev1pa90Q");
01847 if( DataCadence == 120.0) strcpy(HMISeriesLev1pa,"hmi.HMISeriesLev1pa120Q");
01848 if( DataCadence == 135.0) strcpy(HMISeriesLev1pa,"hmi.HMISeriesLev1pa135Q");
01849 if( DataCadence == 150.0) strcpy(HMISeriesLev1pa,"hmi.HMISeriesLev1pa150Q");
01850 if( DataCadence == 720.0) strcpy(HMISeriesLev1pa,"hmi.S_720s_nrt");
01851 }
01852 else
01853 {
01854 if( DataCadence == 45.0) strcpy(HMISeriesLev1pa,"hmi.HMISeriesLev1pa45");
01855 if( DataCadence == 90.0) strcpy(HMISeriesLev1pa,"hmi.HMISeriesLev1pa90");
01856 if( DataCadence == 120.0) strcpy(HMISeriesLev1pa,"hmi.HMISeriesLev1pa120");
01857 if( DataCadence == 135.0) strcpy(HMISeriesLev1pa,"hmi.HMISeriesLev1pa135");
01858 if( DataCadence == 150.0) strcpy(HMISeriesLev1pa,"hmi.HMISeriesLev1pa150");
01859 if( DataCadence == 720.0) strcpy(HMISeriesLev1pa,"hmi.S_720s");
01860 }
01861
01862
01863
01864
01865
01866
01867
01868
01869
01870
01871
01872 if ( TestLevIn[0]==1 )
01873 {
01874
01875
01876
01877
01878
01879
01880
01881 TimeBegin2=TimeBegin-(TIME)TempIntNum*DataCadence/2.-TimeCaution;
01882 TimeEnd2 =TimeEnd +(TIME)TempIntNum*DataCadence/2.+TimeCaution;
01883 sprint_time(timeBegin2,TimeBegin2,"TAI",0);
01884 sprint_time(timeEnd2,TimeEnd2,"TAI",0);
01885 strcat(HMISeriesLev1,"[");
01886 strcat(HMISeriesLev1,timeBegin2);
01887 strcat(HMISeriesLev1,"-");
01888 strcat(HMISeriesLev1,timeEnd2);
01889 strcat(HMISeriesLev1,"]");
01890
01891
01892
01893
01894 printf("LEVEL 1 SERIES QUERY = %s\n",HMISeriesLev1);
01895 recLev1 = drms_open_records(drms_env,HMISeriesLev1,&status);
01896 if (status == DRMS_SUCCESS && recLev1 != NULL && recLev1->n > 0)
01897 {
01898 nRecs1 = recLev1->n;
01899
01900 if(nRecs1 >= nRecmax)
01901 {
01902 printf("Too many records requested\n");
01903 return 1;
01904 }
01905
01906 printf("Number of level 1 records opened: %d\n",nRecs1);
01907
01908
01909 FSN = (int *)malloc(nRecs1*sizeof(int));
01910 if(FSN == NULL)
01911 {
01912 printf("Error: memory could not be allocated to FSN\n");
01913 return 1;
01914 }
01915 internTOBS = (TIME *)malloc(nRecs1*sizeof(TIME));
01916 if(internTOBS == NULL)
01917 {
01918 printf("Error: memory could not be allocated to internTOBS\n");
01919 return 1;
01920 }
01921 HWL1POS = (int *)malloc(nRecs1*sizeof(int));
01922 if(HWL1POS == NULL)
01923 {
01924 printf("Error: memory could not be allocated to HWL1POS\n");
01925 return 1;
01926 }
01927 HWL2POS = (int *)malloc(nRecs1*sizeof(int));
01928 if(HWL2POS == NULL)
01929 {
01930 printf("Error: memory could not be allocated to HWL2POS\n");
01931 return 1;
01932 }
01933 HWL3POS = (int *)malloc(nRecs1*sizeof(int));
01934 if(HWL3POS == NULL)
01935 {
01936 printf("Error: memory could not be allocated to HWL3POS\n");
01937 return 1;
01938 }
01939 HWL4POS = (int *)malloc(nRecs1*sizeof(int));
01940 if(HWL4POS == NULL)
01941 {
01942 printf("Error: memory could not be allocated to HWL4POS\n");
01943 return 1;
01944 }
01945 HPL1POS = (int *)malloc(nRecs1*sizeof(int));
01946 if(HPL1POS == NULL)
01947 {
01948 printf("Error: memory could not be allocated to HPL1POS\n");
01949 return 1;
01950 }
01951 HPL2POS = (int *)malloc(nRecs1*sizeof(int));
01952 if(HPL2POS == NULL)
01953 {
01954 printf("Error: memory could not be allocated to HPL2POS\n");
01955 return 1;
01956 }
01957 HPL3POS = (int *)malloc(nRecs1*sizeof(int));
01958 if(HPL3POS == NULL)
01959 {
01960 printf("Error: memory could not be allocated to HPL3POS\n");
01961 return 1;
01962 }
01963 FID = (int *)malloc(nRecs1*sizeof(int));
01964 if(FID == NULL)
01965 {
01966 printf("Error: memory could not be allocated to FID\n");
01967 return 1;
01968 }
01969 HFLID = (int *)malloc(nRecs1*sizeof(int));
01970 if(HFLID == NULL)
01971 {
01972 printf("Error: memory could not be allocated to HFLID\n");
01973 return 1;
01974 }
01975 HCAMID = (int *)malloc(nRecs1*sizeof(int));
01976 if(HCAMID == NULL)
01977 {
01978 printf("Error: memory could not be allocated to HCAMID\n");
01979 return 1;
01980 }
01981 RSUN = (float *)malloc(nRecs1*sizeof(float));
01982 if(RSUN == NULL)
01983 {
01984 printf("Error: memory could not be allocated to RSUN\n");
01985 return 1;
01986 }
01987 CROTA2 = (float *)malloc(nRecs1*sizeof(float));
01988 if(CROTA2 == NULL)
01989 {
01990 printf("Error: memory could not be allocated to CROTA2\n");
01991 return 1;
01992 }
01993 CRLTOBS = (float *)malloc(nRecs1*sizeof(float));
01994 if(CRLTOBS == NULL)
01995 {
01996 printf("Error: memory could not be allocated to CRLTOBS\n");
01997 return 1;
01998 }
01999 DSUNOBS = (double *)malloc(nRecs1*sizeof(double));
02000 if(DSUNOBS == NULL)
02001 {
02002 printf("Error: memory could not be allocated to DSUNOBS\n");
02003 return 1;
02004 }
02005 X0 = (float *)malloc(nRecs1*sizeof(float));
02006 if(X0 == NULL)
02007 {
02008 printf("Error: memory could not be allocated to X0\n");
02009 return 1;
02010 }
02011 Y0 = (float *)malloc(nRecs1*sizeof(float));
02012 if(Y0 == NULL)
02013 {
02014 printf("Error: memory could not be allocated to Y0\n");
02015 return 1;
02016 }
02017 SegmentRead= (int *)malloc(nRecs1*sizeof(int));
02018 if(SegmentRead == NULL)
02019 {
02020 printf("Error: memory could not be allocated to SegmentRead\n");
02021 return 1;
02022 }
02023 KeywordMissing= (int *)malloc(nRecs1*sizeof(int));
02024 if(KeywordMissing == NULL)
02025 {
02026 printf("Error: memory could not be allocated to KeywordMissing\n");
02027 return 1;
02028 }
02029 Segments = (DRMS_Array_t **)malloc(nRecs1*sizeof(DRMS_Array_t *));
02030 if(Segments == NULL)
02031 {
02032 printf("Error: memory could not be allocated to Segments\n");
02033 return 1;
02034 }
02035 for(i=0;i<nRecs1;i++) Segments[i]=NULL;
02036 Ierror = (DRMS_Array_t **)malloc(nRecs1*sizeof(DRMS_Array_t *));
02037 if(Ierror == NULL)
02038 {
02039 printf("Error: memory could not be allocated to Ierror\n");
02040 return 1;
02041 }
02042 for(i=0;i<nRecs1;i++) Ierror[i]=NULL;
02043 IndexFiltergram = (int *)malloc(nRecs1*sizeof(int));
02044 if(IndexFiltergram == NULL)
02045 {
02046 printf("Error: memory could not be allocated to IndexFiltergram\n");
02047 return 1;
02048 }
02049 CFINDEX = (int *)malloc(nRecs1*sizeof(int));
02050 if(CFINDEX == NULL)
02051 {
02052 printf("Error: memory could not be allocated to CFINDEX\n");
02053 return 1;
02054 }
02055 HIMGCFID = (int *)malloc(nRecs1*sizeof(int));
02056 if(HIMGCFID == NULL)
02057 {
02058 printf("Error: memory could not be allocated to HIMGCFID\n");
02059 return 1;
02060 }
02061 IMGTYPE = (char **)malloc(nRecs1*sizeof(char *));
02062 if(IMGTYPE == NULL)
02063 {
02064 printf("Error: memory could not be allocated to IMGTYPE\n");
02065 return 1;
02066 }
02067 for(i=0;i<nRecs1;i++) IMGTYPE[i]=NULL;
02068 CDELT1 = (float *)malloc(nRecs1*sizeof(float));
02069 if(CDELT1 == NULL)
02070 {
02071 printf("Error: memory could not be allocated to CDELT1\n");
02072 return 1;
02073 }
02074 OBSVR = (double *)malloc(nRecs1*sizeof(double));
02075 if(OBSVR == NULL)
02076 {
02077 printf("Error: memory could not be allocated to OBSVR\n");
02078 return 1;
02079 }
02080 OBSVW = (double *)malloc(nRecs1*sizeof(double));
02081 if(OBSVW == NULL)
02082 {
02083 printf("Error: memory could not be allocated to OBSVW\n");
02084 return 1;
02085 }
02086 OBSVN = (double *)malloc(nRecs1*sizeof(double));
02087 if(OBSVN == NULL)
02088 {
02089 printf("Error: memory could not be allocated to OBSVN\n");
02090 return 1;
02091 }
02092 CRLNOBS = (float *)malloc(nRecs1*sizeof(float));
02093 if(CRLNOBS == NULL)
02094 {
02095 printf("Error: memory could not be allocated to CRLNOBS\n");
02096 return 1;
02097 }
02098
02099
02100
02101
02102
02103
02104 CARROT = (int *)malloc(nRecs1*sizeof(int));
02105 if(CARROT == NULL)
02106 {
02107 printf("Error: memory could not be allocated to CARROT\n");
02108 return 1;
02109 }
02110 HWLTID = (int *)malloc(nRecs1*sizeof(int));
02111 if(HWLTID == NULL)
02112 {
02113 printf("Error: memory could not be allocated to HWLTID\n");
02114 return 1;
02115 }
02116 HPLTID = (int *)malloc(nRecs1*sizeof(int));
02117 if(HPLTID == NULL)
02118 {
02119 printf("Error: memory could not be allocated to HPLTID\n");
02120 return 1;
02121 }
02122 HWLTNSET = (char **)malloc(nRecs1*sizeof(char *));
02123 if(HWLTNSET == NULL)
02124 {
02125 printf("Error: memory could not be allocated to HWLTNSET\n");
02126 return 1;
02127 }
02128 NBADPERM = (int *)malloc(nRecs1*sizeof(int));
02129 if(NBADPERM == NULL)
02130 {
02131 printf("Error: memory could not be allocated to NBADPERM\n");
02132 return 1;
02133 }
02134 QUALITYlev1 = (int *)malloc(nRecs1*sizeof(int));
02135 if(QUALITYlev1 == NULL)
02136 {
02137 printf("Error: memory could not be allocated to QUALITYlev1\n");
02138 return 1;
02139 }
02140 for(i=0;i<nRecs1;++i) QUALITYlev1[i]=0;
02141 QUALITYin = (int *)malloc(nRecs1*sizeof(int));
02142 if(QUALITYin == NULL)
02143 {
02144 printf("Error: memory could not be allocated to QUALITYin\n");
02145 return 1;
02146 }
02147
02148
02149
02150
02151
02152
02153
02154 t0=dsecnd();
02155
02156 k=0;
02157 for(i=0;i<nRecs1;++i)
02158 {
02159 FSN[i] = drms_getkey_int(recLev1->records[i] ,FSNS ,&statusA[0]);
02160
02161
02162 internTOBS[i] = drms_getkey_time(recLev1->records[i],TOBSS ,&statusA[1]);
02163 HWL1POS[i] = drms_getkey_int(recLev1->records[i] ,HWL1POSS ,&statusA[2]);
02164 HWL2POS[i] = drms_getkey_int(recLev1->records[i] ,HWL2POSS ,&statusA[3]);
02165 HWL3POS[i] = drms_getkey_int(recLev1->records[i] ,HWL3POSS ,&statusA[4]);
02166 HWL4POS[i] = drms_getkey_int(recLev1->records[i] ,HWL4POSS ,&statusA[5]);
02167 HPL1POS[i] = drms_getkey_int(recLev1->records[i] ,HPL1POSS ,&statusA[6]);
02168 HPL2POS[i] = drms_getkey_int(recLev1->records[i] ,HPL2POSS ,&statusA[7]);
02169 HPL3POS[i] = drms_getkey_int(recLev1->records[i] ,HPL3POSS ,&statusA[8]);
02170 FID[i] = drms_getkey_int(recLev1->records[i] ,FIDS ,&statusA[9]);
02171 HCAMID[i] = drms_getkey_int(recLev1->records[i] ,HCAMIDS ,&statusA[10]);
02172 if(HCAMID[i] != LIGHT_SIDE && HCAMID[i] != LIGHT_FRONT) statusA[11]=1;
02173 CFINDEX[i] = drms_getkey_int(recLev1->records[i] ,HCFTIDS ,&statusA[11]);
02174
02175
02176
02177 HFLID[i] = drms_getkey_int(recLev1->records[i] ,"HFTSACID" ,&statusA[12]);
02178
02179
02180
02181
02182 if(CamId == LIGHT_SIDE)
02183 {
02184 iii=needtochangeFID(HFLID[i]);
02185 if(iii == 1)
02186 {
02187 ii = drms_getkey_int(recLev1->records[i],"HFLPSITN",&status);
02188 if( (ii%72)/24 == 1)
02189 {
02190 FID[i]=FID[i]+100000;
02191 }
02192 }
02193
02194 }
02195
02196
02197 IMGTYPE[i] = (char *)malloc(6*sizeof(char *));
02198 if(IMGTYPE[i] == NULL)
02199 {
02200 printf("Error: memory could not be allocated to IMGTYPE[%d]\n",i);
02201 return 1;
02202 }
02203 IMGTYPE[i] = drms_getkey_string(recLev1->records[i],IMGTYPES ,&statusA[13]);
02204
02205 X0[i] = (float)drms_getkey_double(recLev1->records[i],CRPIX1S, &statusA[14]);
02206 if(statusA[14] == DRMS_SUCCESS && !isnan(X0[i])) X0[i]=X0[i]-1.0;
02207 else statusA[14] = 1;
02208
02209 Y0[i] = (float)drms_getkey_double(recLev1->records[i],CRPIX2S,&statusA[15]);
02210 if(statusA[15] == DRMS_SUCCESS && !isnan(Y0[i])) Y0[i]=Y0[i]-1.0;
02211 else statusA[15] = 1;
02212
02213 X0LF = (float)drms_getkey_double(recLev1->records[i],X0LFS, &status);
02214 Y0LF = (float)drms_getkey_double(recLev1->records[i],Y0LFS, &status2);
02215 if(status != DRMS_SUCCESS || status2 != DRMS_SUCCESS || isnan(X0LF) || isnan(Y0LF))
02216 {
02217 statusA[14]=1;
02218 statusA[15]=1;
02219 X0[i]=sqrt(-1);
02220 Y0[i]=sqrt(-1);
02221 }
02222
02223
02224 RSUN[i] = (float)drms_getkey_double(recLev1->records[i],RSUNS ,&statusA[16]);
02225
02226 CROTA2[i] = (float)drms_getkey_double(recLev1->records[i],CROTA2S ,&statusA[17]);
02227
02228
02229
02230 if(statusA[17] == DRMS_SUCCESS && !isnan(CROTA2[i]))
02231 {
02232 if(CROTA2[i] > 362. || CROTA2[i] < -362.) statusA[17] = 1;
02233 }
02234
02235
02236
02237 if(statusA[17] == DRMS_SUCCESS && !isnan(CROTA2[i])) CROTA2[i]=-CROTA2[i];
02238 else statusA[17] = 1;
02239
02240 CRLTOBS[i] = (float)drms_getkey_double(recLev1->records[i],CRLTOBSS,&statusA[18]);
02241 if(isnan(CRLTOBS[i])) statusA[18] = 1;
02242 DSUNOBS[i] = drms_getkey_double(recLev1->records[i],DSUNOBSS,&statusA[19]);
02243 if(isnan(DSUNOBS[i])) statusA[19] = 1;
02244 HIMGCFID[i] = drms_getkey_int(recLev1->records[i] ,HIMGCFIDS ,&statusA[20]);
02245 if(isnan(HIMGCFID[i])) statusA[20] = 1;
02246 CDELT1[i] = (float)drms_getkey_double(recLev1->records[i] ,CDELT1S,&statusA[21]);
02247 if(isnan(CDELT1[i])) statusA[21] = 1;
02248 OBSVR[i] = drms_getkey_double(recLev1->records[i] ,OBSVRS ,&statusA[22]);
02249 if(isnan(OBSVR[i])) statusA[22] = 1;
02250 OBSVW[i] = drms_getkey_double(recLev1->records[i] ,OBSVWS ,&statusA[23]);
02251 if(isnan(OBSVW[i])) statusA[23] = 1;
02252 OBSVN[i] = drms_getkey_double(recLev1->records[i] ,OBSVNS ,&statusA[24]);
02253 if(isnan(OBSVN[i])) statusA[24] = 1;
02254 CARROT[i] = drms_getkey_int(recLev1->records[i] ,CARROTS ,&statusA[25]);
02255 SegmentRead[i]= 0;
02256 KeywordMissing[i]=0;
02257 CRLNOBS[i] = (float)drms_getkey_double(recLev1->records[i] ,CRLNOBSS,&statusA[26]);
02258 if(isnan(CRLNOBS[i])) statusA[26] = 1;
02259
02260
02261
02262 statusA[27]=0;
02263
02264
02265 HWLTID[i] = drms_getkey_int(recLev1->records[i] ,HWLTIDS ,&statusA[28]);
02266 HPLTID[i] = drms_getkey_int(recLev1->records[i] ,HPLTIDS ,&statusA[29]);
02267 HWLTNSET[i] = (char *)malloc(7*sizeof(char *));
02268 if(HWLTNSET[i] == NULL)
02269 {
02270 printf("Error: memory could not be allocated to HWLTNSET[%d]\n",i);
02271 return 1;
02272 }
02273 HWLTNSET[i] = drms_getkey_string(recLev1->records[i] ,HWLTNSETS ,&statusA[30]);
02274 NBADPERM[i] = drms_getkey_int(recLev1->records[i] ,NBADPERMS ,&statusA[31]);
02275 if(statusA[31] != DRMS_SUCCESS) NBADPERM[i]=-1;
02276 QUALITYin[i] = drms_getkey_int(recLev1->records[i] ,QUALITYS ,&statusA[32]);
02277
02278
02279
02280
02281
02282
02283
02284
02285
02286
02287
02288
02289
02290
02291
02292
02293
02294
02295
02296
02297
02298
02299
02300 TotalStatus=0;
02301 for(ii=0;ii<=30;++ii) TotalStatus+=statusA[ii];
02302 if(TotalStatus != 0 || !strcmp(IMGTYPE[i],"DARK"))
02303 {
02304 printf("Error: the level 1 filtergram index %d is missing at least one keyword\n",i);
02305
02306 internTOBS[i] = MISSINGKEYWORD;
02307 HWL1POS[i] = MISSINGKEYWORDINT;
02308 HWL2POS[i] = MISSINGKEYWORDINT;
02309 HWL3POS[i] = MISSINGKEYWORDINT;
02310 HWL4POS[i] = MISSINGKEYWORDINT;
02311 HPL1POS[i] = MISSINGKEYWORDINT;
02312 HPL2POS[i] = MISSINGKEYWORDINT;
02313 HPL3POS[i] = MISSINGKEYWORDINT;
02314 FID[i] = MISSINGKEYWORDINT;
02315 HCAMID[i] = MISSINGKEYWORDINT;
02316 CFINDEX[i] = MISSINGKEYWORDINT;
02317
02318 FSN[i] = MISSINGKEYWORDINT;
02319 X0[i] = MISSINGKEYWORD;
02320 Y0[i] = MISSINGKEYWORD;
02321 RSUN[i] = MISSINGKEYWORD;
02322 CROTA2[i] = MISSINGKEYWORD;
02323 CRLTOBS[i] = MISSINGKEYWORD;
02324 DSUNOBS[i] = MISSINGKEYWORD;
02325 HIMGCFID[i] = MISSINGKEYWORD;
02326 CDELT1[i] = MISSINGKEYWORD;
02327 OBSVR[i] = MISSINGKEYWORD;
02328 OBSVW[i] = MISSINGKEYWORD;
02329 OBSVN[i] = MISSINGKEYWORD;
02330 CARROT[i] = MISSINGKEYWORDINT;
02331 CRLNOBS[i] = MISSINGKEYWORD;
02332
02333 HWLTID[i] = MISSINGKEYWORD;
02334 strcpy(HWLTNSET[i],"NONE");
02335
02336
02337 KeywordMissing[i]=1;
02338 }
02339 else
02340 {
02341 if(WhichWavelength(FID[i]) == WavelengthID && HCAMID[i] == CamId)
02342 {
02343 IndexFiltergram[k]=i;
02344 ++k;
02345 }
02346 }
02347
02348 }
02349 t1=dsecnd();
02350 printf("TIME ELAPSED TO READ THE KEYWORDS OF ALL LEVEL 1 RECORDS: %f\n",t1-t0);
02351
02352 nIndexFiltergram=k;
02353 if(nIndexFiltergram == 0)
02354 {
02355 printf("Error: no filtergram was found with the wavelength %d in the requested level 1 records %s\n",WavelengthID,HMISeriesLev1);
02356
02357 }
02358
02359 }
02360 else
02361 {
02362 printf("Error: unable to open the requested level 1 records %s\n",HMISeriesLev1);
02363 return 1;
02364 }
02365
02366 }
02367
02368
02369
02370
02371
02372
02373
02374
02375
02376
02377
02378
02379
02380
02381
02382 TREC_STEP = DataCadence;
02383
02384 TargetTime = (TIME)floor((TimeBegin-TREC_EPOCH0+TREC_STEP/2.0)/TREC_STEP)*TREC_STEP+TREC_EPOCH0;
02385 if(TargetTime < TimeBegin) TargetTime+=TREC_STEP;
02386 initialrun=1;
02387 PreviousTargetTime=TargetTime;
02388
02389
02390
02391 while(TargetTime <= TimeEnd)
02392 {
02393
02394
02395 if(inRotationalFlat == 1)
02396 {
02397
02398
02399 if(floor(TargetTime/86400.0) != floor(PreviousTargetTime/86400.0) )
02400 {
02401 printf("Error: the new target time is for a different day than the previous target time: you are not allowed to change day when applying a rotational flat field\n");
02402 return 1;
02403 }
02404
02405
02406 if(initialrun != 1)
02407 {
02408 if(TargetTime < TSTARTFLAT || TargetTime > TSTOPFLAT)
02409 {
02410 printf("Error: the target time is not within the time range for which the rotation flat field used is valid\n");
02411 return 1;
02412 }
02413 }
02414
02415 }
02416
02417
02418 sprint_time(timeBegin2,TargetTime,"TAI",0);
02419 printf("\n TARGET TIME = %s\n",timeBegin2);
02420 printf("-----------------------------------------------------------------------------------\n");
02421
02422 if(nIndexFiltergram == 0 && TestLevIn[0] == 1 )
02423 {
02424 QUALITY = QUALITY | QUAL_TARGETFILTERGRAMMISSING;
02425 CreateEmptyRecord = 1;
02426 goto NextTargetTime;
02427 }
02428
02429 QUALITY = 0;
02430 QUALITYLEV1 = 0;
02431 strcpy(HISTORY,"");
02432
02433
02434
02435
02436
02437
02438
02439
02440
02441 if (TestLevIn[1]==1)
02442 {
02443
02444 strcpy(source,"[RECORDS USED: NOT REPORTED");
02445
02446 strcpy(HMISeries,HMISeriesLev1d);
02447 strcat(HMISeries,"[");
02448 strcat(HMISeries,timeBegin2);
02449 strcat(HMISeries,"][][");
02450 sprintf(CamIds,"%d",CamId);
02451 strcat(HMISeries,CamIds);
02452 strcat(HMISeries,"]");
02453
02454 printf("LEVEL 1d QUERY= %s\n",HMISeries);
02455 recLev1d = drms_open_records(drms_env,HMISeries,&status);
02456
02457 if (status == DRMS_SUCCESS && recLev1d != NULL && recLev1d->n > 0)
02458 {
02459 nRecs1d = recLev1d->n;
02460
02461 if(nRecs1d >= MaxNumFiltergrams)
02462 {
02463 printf("Number of open record is larger than %d\n",MaxNumFiltergrams);
02464 return 1;
02465 }
02466 printf("Number of level 1d records opened= %d\n",nRecs1d);
02467
02468 trec = drms_getkey_time(recLev1d->records[0],TRECS,&status);
02469 if(status != DRMS_SUCCESS)
02470 {
02471 printf("Error: unable to read the %s keyword of level 1d data at target time %s\n",TRECS,timeBegin2);
02472 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
02473 }
02474 if(trec != TargetTime)
02475 {
02476 printf("Error: %s of a level 1d record is not equal to the target time %s\n",TRECS,timeBegin2);
02477 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
02478 }
02479 tobs = drms_getkey_time(recLev1d->records[0],TOBSS,&status);
02480 if(status != DRMS_SUCCESS)
02481 {
02482 printf("Error: unable to read the %s keyword of level 1d data at target time %s\n",TOBSS,timeBegin2);
02483 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
02484 }
02485
02486
02487
02488
02489
02490 TREC_EPOCH=drms_getkey_time(recLev1d->records[0],TRECEPOCHS,&status);
02491 if(TREC_EPOCH != TREC_EPOCH0)
02492 {
02493 printf("Error: TREC_EPOCH of level 1d data is not equal to the expected TREC_EPOCH: %f, at target time %s\n",TREC_EPOCH0,timeBegin2);
02494 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
02495 }
02496 TREC_STEP= drms_getkey_time(recLev1d->records[0],TRECSTEPS,&status);
02497 if(TREC_STEP != DataCadence)
02498 {
02499 printf("Error: the cadence is not equal to the T_REC_step keyword of the level 1d data, at target time %s\n",timeBegin2);
02500 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
02501 }
02502
02503 DRMS_SegmentDimInfo_t di;
02504 segin = drms_segment_lookupnum(recLev1d->records[0],0);
02505 status = drms_segment_getdims(segin,&di);
02506 if(status != DRMS_SUCCESS)
02507 {
02508 printf("Error: unable to read the dimensions of the data segment of level 1d data at target time %s\n",timeBegin2);
02509 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
02510 }
02511 axisin[0]= di.axis[0];
02512 axisin[1]= di.axis[1];
02513 axisout[0]=axisin[0];
02514 axisout[1]=axisin[1];
02515
02516 Segments1d=0;
02517 Segments1p=0;
02518 }
02519 else
02520 {
02521 printf("Unable to open the series %s for target time %s\n",HMISeries,timeBegin2);
02522 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
02523 }
02524 }
02525
02526
02527
02528
02529
02530
02531
02532
02533
02534
02535
02536 if ( TestLevIn[0]==1 )
02537 {
02538
02539 strcpy(source,HMISeriesLev10);
02540 strcat(source,"[:");
02541
02542
02543
02544
02545
02546 temptime = 100000000.0;
02547 i=TargetWavelength;
02548 while(fabs(internTOBS[IndexFiltergram[i]]-TargetTime) <= temptime)
02549 {
02550
02551 temptime=fabs(internTOBS[IndexFiltergram[i]]-TargetTime);
02552 if(i <= nIndexFiltergram-2) ++i;
02553 else break;
02554 }
02555 if(temptime > DataCadence/2.0 )
02556 {
02557 printf("Error: could not find a filtergram with the correct wavelength and close to the target time\n");
02558 QUALITY = QUALITY | QUAL_TARGETFILTERGRAMMISSING;
02559 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
02560 }
02561 TargetWavelength= i-1;
02562 temp = IndexFiltergram[TargetWavelength];
02563
02564 if(SegmentRead[temp] == 0)
02565 {
02566
02567
02568
02569
02570
02571 if(inRotationalFlat == 1)
02572 {
02573
02574 HMIFlatField = (char *)malloc(MaxNString*sizeof(char *));
02575 if(HMIFlatField == NULL)
02576 {
02577 printf("Error: memory could not be allocated to HMIFlatField\n");
02578 return 1;
02579 }
02580
02581 HMIFlatField = drms_getkey_string(recLev1->records[temp],FLATREC,&status);
02582 if (status != DRMS_SUCCESS)
02583 {
02584 printf("Error: could not read the FLAT_REC keyword for the target filtergram FSN= %d",FSN[temp]);
02585 return 1;
02586 }
02587 else
02588 {
02589 printf("PZT FLAT FIELD USED ON TARGET FILTERGRAM= %s\n",HMIFlatField);
02590 if(initialrun == 1)
02591 {
02592 strcpy(HMIFlatField0,HMIFlatField);
02593
02594
02595 recflat = drms_open_records(drms_env,HMIFlatField,&status);
02596 if (status != DRMS_SUCCESS || recflat == NULL || recflat->n == 0)
02597 {
02598 printf("Error: record missing or corrupt for the flat field query %s\n",HMIFlatField);
02599 return 1;
02600 }
02601 segin = drms_segment_lookup(recflat->records[0],"flatfield");
02602 flatfield = drms_segment_read(segin,type1d,&status);
02603 if (status != DRMS_SUCCESS || flatfield == NULL)
02604 {
02605 printf("Error: could not read the data segment for the flat field query %s\n",HMIFlatField);
02606 return 1;
02607 }
02608 pztflat = flatfield->data;
02609 status=drms_close_records(recflat,DRMS_FREE_RECORD);
02610
02611
02612 char keylist[]="T_START,T_STOP";
02613 int unique = 0;
02614 rotationalflats = drms_record_getvector(drms_env,HMIRotationalFlats,keylist,DRMS_TYPE_DOUBLE,unique,&status);
02615 if(status != DRMS_SUCCESS)
02616 {
02617 printf("Error: cannot read a list of keywords in the rotation flat-field series\n");
02618 return 1;
02619 }
02620 printf("DIMENSIONS OF ROTATIONAL FLAT-FIELD SERIES= %d %d\n",rotationalflats->axis[0],rotationalflats->axis[1]);
02621 int n0,n1;
02622 n1 =rotationalflats->axis[1];
02623 n0 =rotationalflats->axis[0];
02624 keyF=rotationalflats->data;
02625
02626 i=0;
02627 while(TargetTime > keyF[i] && TargetTime > keyF[n1+i])
02628 {
02629 i++;
02630 }
02631
02632 if(TargetTime >= keyF[i] && TargetTime <= keyF[n1+i])
02633 {
02634
02635 TSTARTFLAT=keyF[i];
02636 TSTOPFLAT =keyF[n1+i];
02637 sprint_time(query,keyF[i],"TAI",1);
02638 strcpy(QueryFlatField,HMIRotationalFlats);
02639 strcat(QueryFlatField,"[");
02640 if(CamId == LIGHT_SIDE) strcat(QueryFlatField,"1");
02641 if(CamId == LIGHT_FRONT) strcat(QueryFlatField,"2");
02642 strcat(QueryFlatField,"][");
02643 strcat(QueryFlatField,query);
02644 strcat(QueryFlatField,"]");
02645 printf("QUERY FOR ROTATIONAL FLAT FIELD= %s\n",QueryFlatField);
02646
02647 recflatrot = drms_open_records(drms_env,QueryFlatField,&status);
02648 if (status != DRMS_SUCCESS || recflatrot == NULL || recflatrot->n == 0)
02649 {
02650 printf("Error: record missing or corrupt for the rotational flat field query %s\n",QueryFlatField);
02651 return 1;
02652 }
02653 segin = drms_segment_lookup(recflatrot->records[0],"flatfield");
02654 flatfieldrot = drms_segment_read(segin,type1d,&status);
02655 if (status != DRMS_SUCCESS || flatfieldrot == NULL)
02656 {
02657 printf("Error: could not read the data segment for the rotational flat field query %s\n",QueryFlatField);
02658 return 1;
02659 }
02660 rotflat = flatfieldrot->data;
02661
02662 }
02663 else
02664 {
02665 printf("Error: no rotational flat field record exists for the target time %s\n",timeBegin2);
02666 return 1;
02667 }
02668
02669 drms_free_array(rotationalflats);
02670 rotationalflats=NULL;
02671 }
02672 else if(strcmp(HMIFlatField,HMIFlatField0) != 0)
02673 {
02674 printf("Warning: the hmi.flatfield record used to produce the level 1 records changed during the run of the observables code\n");
02675 printf("The new hmi,flatfield record used is: %s\n",HMIFlatField);
02676
02677 recflat = drms_open_records(drms_env,HMIFlatField,&status);
02678 if (status != DRMS_SUCCESS || recflat == NULL || recflat->n == 0)
02679 {
02680 printf("Error: record missing or corrupt for the flat field query %s\n",HMIFlatField);
02681 return 1;
02682 }
02683 drms_free_array(flatfield);
02684 strcpy(HMIFlatField0,HMIFlatField);
02685 segin = drms_segment_lookup(recflat->records[0],"flatfield");
02686 flatfield = drms_segment_read(segin,type1d,&status);
02687 if (status != DRMS_SUCCESS || flatfield == NULL)
02688 {
02689 printf("Error: could not read the data segment for the flat field query %s\n",HMIFlatField);
02690 return 1;
02691 }
02692 pztflat = flatfield->data;
02693 status=drms_close_records(recflat,DRMS_FREE_RECORD);
02694 }
02695
02696 }
02697 free(HMIFlatField);
02698 }
02699
02700
02701
02702
02703
02704 printf("READ SEGMENT OF TARGET FILTERGRAM\n");
02705 segin = drms_segment_lookupnum(recLev1->records[temp],0);
02706 Segments[temp] = drms_segment_read(segin,type1d, &status);
02707 if (status != DRMS_SUCCESS || Segments[temp] == NULL)
02708 {
02709 Segments[temp]=NULL;
02710 SegmentRead[temp]=-1;
02711 image = NULL;
02712
02713 }
02714 }
02715
02716 if(SegmentRead[temp] == -1)
02717 {
02718
02719 image = NULL;
02720 }
02721 else
02722 {
02723 image = Segments[temp]->data;
02724 }
02725
02726 printf("FSN OF TARGET FILTERGRAM = %d %d\n",FSN[temp],HCAMID[temp]);
02727 if(internTOBS[temp] > TargetTime) printf("Target filtergram is after target time\n");
02728 else printf("Target filtergram is before target time\n");
02729
02730 TargetHFLID = HFLID[temp];
02731 TargetHWLPOS[0] = HWL1POS[temp];
02732 TargetHWLPOS[1] = HWL2POS[temp];
02733 TargetHWLPOS[2] = HWL3POS[temp];
02734 TargetHWLPOS[3] = HWL4POS[temp];
02735 TargetHPLPOS[0] = HPL1POS[temp];
02736 TargetHPLPOS[1] = HPL2POS[temp];
02737 TargetHPLPOS[2] = HPL3POS[temp];
02738 TargetCFINDEX = CFINDEX[temp];
02739 TargetHWLTID = HWLTID[temp];
02740 TargetHPLTID = HPLTID[temp];
02741 strcpy(TargetISS,HWLTNSET[temp]);
02742 printf("ISS STATUS OF TARGET FILTERGRAM: %s\n",TargetISS);
02743 if(!strcmp(TargetISS,"OPEN")) QUALITY = QUALITY | QUAL_ISSTARGET;
02744 if( (QUALITYin[temp] & Q_ACS_ECLP) == Q_ACS_ECLP) QUALITY = QUALITY | QUAL_ECLIPSE;
02745
02746 axisin[0] = 4096;
02747 axisin[1] = 4096;
02748 axisout[0] = axisin[0];
02749 axisout[1] = axisin[1];
02750 Nelem = axisin[0]*axisin[1];
02751
02752
02753
02754
02755
02756
02757 i=IndexFiltergram[TargetWavelength];
02758 j=i;
02759
02760 if(internTOBS[IndexFiltergram[TargetWavelength]] < TargetTime && TargetWavelength < nIndexFiltergram-1)
02761 {
02762 j=TargetWavelength;
02763 while(internTOBS[IndexFiltergram[j]] < TargetTime && j < nIndexFiltergram-1 && fabs(internTOBS[IndexFiltergram[j]]-internTOBS[i]) <= 2.1*DataCadence) j++;
02764 if(fabs(internTOBS[IndexFiltergram[j]]-internTOBS[i]) > 2.1*DataCadence) j=i;
02765 else j=IndexFiltergram[j];
02766 }
02767 if(internTOBS[IndexFiltergram[TargetWavelength]] >= TargetTime && TargetWavelength > 0)
02768 {
02769 j=TargetWavelength;
02770 while(internTOBS[IndexFiltergram[j]] >= TargetTime && j > 0 && fabs(internTOBS[IndexFiltergram[j]]-internTOBS[i]) <= 2.1*DataCadence) --j;
02771 if(fabs(internTOBS[IndexFiltergram[j]]-internTOBS[i]) > 2.1*DataCadence) j=i;
02772 else
02773 {
02774 i=IndexFiltergram[j];
02775 j=IndexFiltergram[TargetWavelength];
02776 }
02777 }
02778
02779
02780 if(i == j)
02781 {
02782 QUALITY = QUALITY | QUAL_LOWKEYWORDNUM;
02783 if(TargetWavelength>=1 && TargetWavelength<nIndexFiltergram-1)
02784 {
02785 i=IndexFiltergram[TargetWavelength-1];
02786 j=IndexFiltergram[TargetWavelength+1];
02787 }
02788 if(TargetWavelength>=1 && TargetWavelength == nIndexFiltergram-1 )
02789 {
02790 i=IndexFiltergram[TargetWavelength-1];
02791 j=IndexFiltergram[TargetWavelength];
02792 }
02793 if(TargetWavelength<nIndexFiltergram-1 && TargetWavelength == 0)
02794 {
02795 i=IndexFiltergram[TargetWavelength];
02796 j=IndexFiltergram[TargetWavelength+1];
02797 }
02798 }
02799
02800
02801 if(KeywordMissing[j] == 1 || KeywordMissing[i] == 1)
02802 {
02803 printf("Error: some keywords are missing/corrupted to interpolate OBS_VR, OBS_VW, OBS_VN, CRLN_OBS, CROTA2, and CAR_ROT at target time %s\n",timeBegin2);
02804 if(SegmentRead[temp])
02805 {
02806 drms_free_array(Segments[temp]);
02807 Segments[temp]=NULL;
02808 if(Ierror[temp] != NULL) drms_free_array(Ierror[temp]);
02809 Ierror[temp]=NULL;
02810 SegmentRead[temp]=0;
02811 }
02812 QUALITY = QUALITY | QUAL_NOINTERPOLATEDKEYWORDS;
02813 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
02814 }
02815
02816
02817
02818
02819
02820 if(i != j)
02821 {
02822 DSUNOBSint=(DSUNOBS[j]-DSUNOBS[i])/(internTOBS[j]-internTOBS[i])*(TargetTime-internTOBS[i])+DSUNOBS[i];
02823 DSUNOBSint=DSUNOBSint/(double)AstroUnit;
02824 }
02825 else
02826 {
02827 DSUNOBSint=DSUNOBS[i];
02828 DSUNOBSint=DSUNOBSint/(double)AstroUnit;
02829 }
02830
02831
02832 tobs = TargetTime+(DSUNOBSint-1.0)/2.00398880422056639358e-03;
02833
02834
02835
02836 if(i != j)
02837 {
02838 printf("FSNs used for interpolation of OBS_VR, OBS_VW, OBS_VN, CRLN_OBS, CROTA2, and CAR_ROT: %d %d %d %d\n",FSN[j],FSN[i],HCAMID[j],HCAMID[i]);
02839
02840 DSUNOBSint=(DSUNOBS[j]-DSUNOBS[i])/(internTOBS[j]-internTOBS[i])*(tobs-internTOBS[i])+DSUNOBS[i];
02841 DSUNOBSint=DSUNOBSint/(double)AstroUnit;
02842 CRLTOBSint=(CRLTOBS[j]-CRLTOBS[i])/(internTOBS[j]-internTOBS[i])*(tobs-internTOBS[i])+CRLTOBS[i];
02843 CROTA2int =(CROTA2[j]-CROTA2[i])/(internTOBS[j]-internTOBS[i])*(tobs-internTOBS[i])+CROTA2[i];
02844 OBSVRint =(OBSVR[j]-OBSVR[i])/(internTOBS[j]-internTOBS[i])*(tobs-internTOBS[i])+OBSVR[i];
02845 OBSVWint =(OBSVW[j]-OBSVW[i])/(internTOBS[j]-internTOBS[i])*(tobs-internTOBS[i])+OBSVW[i];
02846 OBSVNint =(OBSVN[j]-OBSVN[i])/(internTOBS[j]-internTOBS[i])*(tobs-internTOBS[i])+OBSVN[i];
02847 ctime1 =-CRLNOBS[i];
02848 ctime2 =360.0*(float)(CARROT[j]-CARROT[i])-CRLNOBS[j];
02849 CRLNOBSint=(ctime2-ctime1)/(internTOBS[j]-internTOBS[i])*(tobs-internTOBS[i])+ctime1;
02850 if(CARROT[j] > CARROT[i])
02851 {
02852 if(CRLNOBSint > 0.0)
02853 {
02854 CRLNOBSint = 360.0 - CRLNOBSint;
02855 CARROTint = CARROT[j];
02856 }
02857 else
02858 {
02859 CRLNOBSint = -CRLNOBSint;
02860 CARROTint = CARROT[i];
02861 }
02862 }
02863 else
02864 {
02865 CRLNOBSint = -CRLNOBSint;
02866 CARROTint = CARROT[i];
02867 }
02868
02869 }
02870 else
02871 {
02872 printf("FSN used for interpolation of OBS_VR, OBS_VW, OBS_VN, CRLN_OBS, CROTA2, and CAR_ROT: %d %d\n",FSN[i],HCAMID[i]);
02873
02874 DSUNOBSint=DSUNOBS[i];
02875 DSUNOBSint=DSUNOBSint/(double)AstroUnit;
02876 CRLTOBSint=CRLTOBS[i];
02877 CROTA2int =CROTA2[i];
02878 OBSVRint =OBSVR[i];
02879 OBSVWint =OBSVW[i];
02880 OBSVNint =OBSVN[i];
02881 CRLNOBSint=CRLNOBS[i];
02882 CARROTint =CARROT[i];
02883
02884 QUALITY = QUALITY | QUAL_LOWKEYWORDNUM;
02885 }
02886
02887
02888
02889
02890
02891
02892 Mask = (unsigned char *)malloc(Nelem*sizeof(unsigned char));
02893 if(Mask == NULL)
02894 {
02895 printf("Error: cannot allocate memory for Mask\n");
02896 return 1;
02897 }
02898
02899
02900
02901
02902
02903
02904
02905
02906
02907
02908
02909
02910
02911
02912 if(SegmentRead[temp] == 0 && image != NULL)
02913 {
02914
02915 if(inRotationalFlat == 1)
02916 {
02917 printf("applying rotational flat field on record FSN=%d\n",FSN[temp]);
02918
02919 for(i=0;i<axisin[0]*axisin[1];++i) image[i]=image[i]*pztflat[i];
02920
02921 for(i=0;i<axisin[0]*axisin[1];++i) image[i]=image[i]/rotflat[i];
02922 }
02923
02924
02925 segin = drms_segment_lookup(recLev1->records[temp],"bad_pixel_list");
02926 printf("READ BAD PIXEL LIST OF TARGET FILTERGRAM FSN = %d\n",FSN[temp]);
02927 BadPixels = NULL;
02928 BadPixels = drms_segment_read(segin,segin->info->type,&status);
02929 if(status != DRMS_SUCCESS || BadPixels == NULL)
02930 {
02931 printf("Error: cannot read the list of bad pixels of level 1 filtergram FSN = %d at target time %s \n",FSN[temp],timeBegin2);
02932 drms_free_array(Segments[temp]);
02933 Segments[temp]=NULL;
02934 SegmentRead[temp]=-1;
02935
02936 }
02937 else
02938 {
02939
02940
02941 strcpy(HMISeriesTemp,CosmicRaySeries);
02942 strcat(HMISeriesTemp,"[][");
02943 sprintf(FSNtemps,"%d",FSN[temp]);
02944 strcat(HMISeriesTemp,FSNtemps);
02945 strcat(HMISeriesTemp,"]");
02946 rectemp=NULL;
02947
02948 rectemp=drms_open_records(drms_env,HMISeriesTemp,&statusA[0]);
02949 if(statusA[0] == DRMS_SUCCESS && rectemp != NULL && rectemp->n != 0)
02950 {
02951 segin = drms_segment_lookupnum(rectemp->records[0],0);
02952 CosmicRays = NULL;
02953
02954 CosmicRays = drms_segment_read(segin,segin->info->type,&status);
02955 if(status != DRMS_SUCCESS || CosmicRays == NULL)
02956 {
02957 printf("Error: the list of cosmic-ray hits could not be read for FSN %d\n",FSN[temp]);
02958
02959
02960
02961 QUALITY = QUALITY | QUAL_NOCOSMICRAY;
02962 QUALITYlev1[temp] = QUALITYlev1[temp] | QUAL_NOCOSMICRAY;
02963 CosmicRays = NULL;
02964
02965 }
02966
02967 COSMICCOUNT=drms_getkey_int(rectemp->records[0],COUNTS,&status);
02968 if(status != DRMS_SUCCESS || COSMICCOUNT == -1)
02969 {
02970 QUALITY = QUALITY | QUAL_NOCOSMICRAY;
02971 QUALITYlev1[temp] = QUALITYlev1[temp] | QUAL_NOCOSMICRAY;
02972 }
02973 }
02974 else
02975 {
02976 printf("Unable to open the series %s for FSN %d\n",HMISeriesTemp,FSN[temp]);
02977
02978
02979
02980 QUALITY = QUALITY | QUAL_NOCOSMICRAY;
02981 QUALITYlev1[temp] = QUALITYlev1[temp] | QUAL_NOCOSMICRAY;
02982 CosmicRays = NULL;
02983
02984 }
02985
02986 printf("CREATING MASK FOR GAP-FILLING OF TARGET FILTERGAM\n");
02987 status = MaskCreation(Mask,axisin[0],axisin[1],BadPixels,HIMGCFID[temp],image,CosmicRays,NBADPERM[temp]);
02988 if(status == 1) return 1;
02989
02990 if(BadPixels != NULL)
02991 {
02992 drms_free_array(BadPixels);
02993 BadPixels=NULL;
02994 }
02995 if(CosmicRays != NULL)
02996 {
02997 drms_free_array(CosmicRays);
02998 CosmicRays=NULL;
02999 }
03000 if(rectemp != NULL)
03001 {
03002 drms_close_records(rectemp,DRMS_FREE_RECORD);
03003 rectemp=NULL;
03004 }
03005
03006
03007 if(status != 0)
03008 {
03009 printf("Error: unable to create a mask for the gap filling function for level 1 filtergram FSN = %d at target time %s\n",FSN[temp],timeBegin2);
03010 drms_free_array(Segments[temp]);
03011 Segments[temp]=NULL;
03012 SegmentRead[temp]=-1;
03013
03014 }
03015 else
03016 {
03017 Ierror[temp] = drms_array_create(typeEr,2,axisout,NULL,&status);
03018 if(status != DRMS_SUCCESS || Ierror[temp] == NULL)
03019 {
03020 printf("Error: unable to create an array for Ierror at target time %s\n",timeBegin2);
03021 drms_free_array(Segments[temp]);
03022 Segments[temp]=NULL;
03023 Ierror[temp]=NULL;
03024 SegmentRead[temp]=-1;
03025
03026 }
03027 else
03028 {
03029 ierror = Ierror[temp]->data;
03030 printf("GAP FILLING THE TARGET FILTERGRAM\n");
03031 t0=dsecnd();
03032 status = do_gapfill(image,Mask,&const_param,ierror,axisin[0],axisin[1]);
03033 t1=dsecnd();
03034 printf("TIME ELAPSED TO GAPFILL: %f\n",t1-t0);
03035 if(status != 0)
03036 {
03037 printf("Error: gapfilling code did not work on the level 1 filtergram FSN = %d at target time %s\n",FSN[temp],timeBegin2);
03038 QUALITY = QUALITY | QUAL_NOGAPFILL;
03039 QUALITYlev1[temp] = QUALITYlev1[temp] | QUAL_NOGAPFILL;
03040 }
03041 SegmentRead[temp]=1;
03042 }
03043 }
03044 }
03045 }
03046
03047
03048
03049
03050
03051
03052 printf("GET INFORMATION ABOUT THE FRAMELIST\n");
03053 framelistSize = framelistInfo(TargetHFLID,TargetHPLTID,TargetHWLTID,WavelengthID,PHWPLPOS,WavelengthIndex,WavelengthLocation,&PolarizationType,CamId,&combine,&npol,MaxNumFiltergrams,&CadenceRead,CameraValues,FIDValues);
03054 if(framelistSize == 1) return 1;
03055
03056
03057
03058
03059
03060
03061
03062
03063
03064
03065
03066
03067
03068
03069
03070
03071
03072
03073
03074
03075
03076
03077
03078
03079
03080
03081
03082 printf("INFORMATION REGARDING THE FRAMELIST\n");
03083 printf("CamId= %d ; framelistSize= %d ; PolarizationType= %d ; combine cameras= %d ; npol= %d ; cadence=%f\n",CamId,framelistSize,PolarizationType,combine,npol,CadenceRead);
03084 printf("Framelist= ");
03085 for(i=0;i<framelistSize;++i) printf("I%d ",WavelengthIndex[i]);
03086 printf("\nLocation= ");
03087 for(i=0;i<framelistSize;++i) printf(" %d ",WavelengthLocation[i]);
03088 printf("\nPHWPLPOS \n");
03089 for(i=0;i<framelistSize;++i) printf("%d %d %d %d %d %d %d \n",PHWPLPOS[i*7],PHWPLPOS[i*7+1],PHWPLPOS[i*7+2],PHWPLPOS[i*7+3],PHWPLPOS[i*7+4],PHWPLPOS[i*7+5],PHWPLPOS[i*7+6]);
03090 printf("CAMERA ORDER\n");
03091 for(i=0;i<framelistSize;++i) printf(" %d ",CameraValues[i]);
03092 printf("\nFID ORDER\n");
03093 for(i=0;i<framelistSize;++i) printf(" %d ",FIDValues[i]);
03094 printf("\n");
03095
03096 if(framelistSize == 0)
03097 {
03098 printf("Error: cannot obtain information regarding the framelist for the level 1 filtergram FSN = %d at target time %s\n",FSN[temp],timeBegin2);
03099 QUALITY = QUALITY | QUAL_NOFRAMELISTINFO;
03100 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
03101 }
03102
03103 if(CadenceRead != DataCadence)
03104 {
03105 printf("Error: the cadence from the current framelist for the level 1 filtergram FSN = %d at target time %s is %f and does not match the cadence entered on the command line %f \n",FSN[temp],timeBegin2,CadenceRead,DataCadence);
03106 QUALITY = QUALITY | QUAL_WRONGCADENCE;
03107 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
03108 }
03109
03110
03111
03112
03113
03114 if (Lev1dWanted) recLev1d = drms_create_records(drms_env,framelistSize,HMISeriesLev1d,DRMS_PERMANENT,&status);
03115 else recLev1d = drms_create_records(drms_env,framelistSize,HMISeriesLev1d,DRMS_TRANSIENT,&status);
03116
03117 if(status != DRMS_SUCCESS || recLev1d == NULL || recLev1d->n < framelistSize)
03118 {
03119 printf("Error: cannot create records for the level 1d data %s at target time %s\n",HMISeriesLev1d,timeBegin2);
03120
03121
03122 return 1;
03123 }
03124 nRecs1d = framelistSize;
03125 Segments1d= 0;
03126 Segments1p= 0;
03127 TREC_EPOCH= drms_getkey_time(recLev1d->records[0],TRECEPOCHS,&status);
03128 if(status != DRMS_SUCCESS)
03129 {
03130 printf("Error: unable to read the %s keyword for level 1d data at target time %s\n",TRECEPOCHS,timeBegin2);
03131 return 1;
03132 }
03133 TREC_STEP= drms_getkey_time(recLev1d->records[0],TRECSTEPS,&status);
03134 if(status != DRMS_SUCCESS)
03135 {
03136 printf("Error: unable to read the keyword %s for level 1d data at target time %s\n",TRECSTEPS,timeBegin2);
03137 return 1;
03138 }
03139 if(TREC_STEP != DataCadence)
03140 {
03141 printf("Error: the cadence = %f is not equal to the T_REC_step = %f keyword of the level 1d data at target time %s\n",DataCadence,TREC_STEP,timeBegin2);
03142 return 1;
03143 }
03144 if(TREC_EPOCH != TREC_EPOCH0)
03145 {
03146 printf("Error: TREC_EPOCH= %f is not equal to the expected TREC_EPOCH = %f, at target time %s\n",TREC_EPOCH,TREC_EPOCH0,timeBegin2);
03147 return 1;
03148 }
03149
03150 trec = TargetTime;
03151
03152
03153
03154
03155 temp=-1;
03156 for(i=0;i<framelistSize;++i) if(WavelengthIndex[i] == WavelengthID && PHWPLPOS[i*7+4] == TargetHPLPOS[0] && PHWPLPOS[i*7+5] == TargetHPLPOS[1] && PHWPLPOS[i*7+6] == TargetHPLPOS[2]) temp=i;
03157 if(temp == -1)
03158 {
03159 printf("Error: the target filtergram does not match any frame of the corresponding framelist\n");
03160 QUALITY = QUALITY | QUAL_WRONGTARGET;
03161 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
03162 }
03163 temp=WavelengthLocation[temp];
03164
03165 OrganizeFramelist =(int *)malloc(framelistSize*sizeof(int));
03166 if(OrganizeFramelist == NULL)
03167 {
03168 printf("Error: memory could not be allocated to OrganizeFramelist\n");
03169 return 1;
03170 }
03171 OrganizeFramelist2=(int *)malloc(framelistSize*sizeof(int));
03172 if(OrganizeFramelist2 == NULL)
03173 {
03174 printf("Error: memory could not be allocated to OrganizeFramelist2\n");
03175 return 1;
03176 }
03177
03178 OrganizeFramelist[0]=WavelengthLocation[0]-temp;
03179 OrganizeFramelist2[0]=-1;
03180 for(i=1;i<framelistSize;++i)
03181 {
03182 OrganizeFramelist[i]=WavelengthLocation[i]-temp;
03183 OrganizeFramelist2[i]=signj(OrganizeFramelist[i]);
03184 for(k=0;k<i;++k) if(WavelengthIndex[i] == WavelengthIndex[k]) OrganizeFramelist2[i]= OrganizeFramelist2[k];
03185 }
03186
03187
03188 printf("ORGANIZE FRAMELIST: \n");
03189 for(i=0;i<framelistSize;++i) printf("%d ",OrganizeFramelist[i]);
03190 printf("\n ORGANIZE FRAMELIST2: \n");
03191 for(i=0;i<framelistSize;++i) printf("%d ",OrganizeFramelist2[i]);
03192 printf("\n");
03193
03194
03195
03196
03197
03198
03199
03200 if(TempIntNum > 2) MaxSearchDistanceL=(TIME)(TempIntNum/2-1)*DataCadence+TimeCaution;
03201 MaxSearchDistanceR=(TIME)(TempIntNum/2) *DataCadence+TimeCaution;
03202 FramelistArray = (int *)malloc(framelistSize*TempIntNum*sizeof(int));
03203 if(FramelistArray == NULL)
03204 {
03205 printf("Error: memory could not be allocated to FramelistArray\n");
03206 return 1;
03207 }
03208
03209 for(i=0;i<framelistSize;++i)
03210 {
03211
03212
03213
03214
03215
03216 camera=CameraValues[i];
03217 fidfilt=FIDValues[i];
03218 printf("FID value in framelist: %d\n",fidfilt);
03219
03220 FiltergramLocation=IndexFiltergram[TargetWavelength]+OrganizeFramelist[i];
03221 k=FiltergramLocation;
03222
03223 if(KeywordMissing[k] == 1)
03224 {
03225 FramelistArray[i]=-1;
03226 printf("Error: the filtergram index %d has missing keywords and the code cannot identify it",k);
03227 }
03228 else
03229 {
03230 if(HWL1POS[k] != PHWPLPOS[i*7] || HWL2POS[k] != PHWPLPOS[i*7+1] || HWL3POS[k] != PHWPLPOS[i*7+2] || HWL4POS[k] != PHWPLPOS[i*7+3] || HPL1POS[k] != PHWPLPOS[i*7+4] || HPL2POS[k] != PHWPLPOS[i*7+5] || HPL3POS[k] != PHWPLPOS[i*7+6] || HCAMID[k] != camera || CFINDEX[k] != TargetCFINDEX)
03231 {
03232 printf("Warning: a filtergram near the target location is not what it should be. Looking for the correct filtergram\n");
03233 FiltergramLocation=IndexFiltergram[TargetWavelength];
03234 for(ii=1;ii<framelistSize;++ii)
03235 {
03236 if(OrganizeFramelist[i] > 0)
03237 {
03238 k=FiltergramLocation+ii;
03239 if(k < nRecs1 && KeywordMissing[k] != 1)
03240 {
03241 if(HWL1POS[k] == PHWPLPOS[i*7] && HWL2POS[k] == PHWPLPOS[i*7+1] && HWL3POS[k] == PHWPLPOS[i*7+2] && HWL4POS[k] == PHWPLPOS[i*7+3] && HPL1POS[k] == PHWPLPOS[i*7+4] && HPL2POS[k] == PHWPLPOS[i*7+5] && HPL3POS[k] == PHWPLPOS[i*7+6] && HCAMID[k] == camera && CFINDEX[k] == TargetCFINDEX) break;
03242 }
03243 }
03244 else
03245 {
03246 k=FiltergramLocation-ii;
03247 if(k > 0 && KeywordMissing[k] != 1)
03248 {
03249 if(HWL1POS[k] == PHWPLPOS[i*7] && HWL2POS[k] == PHWPLPOS[i*7+1] && HWL3POS[k] == PHWPLPOS[i*7+2] && HWL4POS[k] == PHWPLPOS[i*7+3] && HPL1POS[k] == PHWPLPOS[i*7+4] && HPL2POS[k] == PHWPLPOS[i*7+5] && HPL3POS[k] == PHWPLPOS[i*7+6] && HCAMID[k] == camera && CFINDEX[k] == TargetCFINDEX) break;
03250 }
03251 }
03252 }
03253 if(ii == framelistSize)
03254 {
03255 FiltergramLocation=IndexFiltergram[TargetWavelength]+OrganizeFramelist[i];
03256 printf("Error: the filtergram FSN = %d is not what it should be\n",FSN[FiltergramLocation]);
03257 FramelistArray[i]=-1;
03258 }
03259 else
03260 {
03261 FramelistArray[i]=k;
03262 FiltergramLocation=k;
03263 }
03264 }
03265 else FramelistArray[i]=k;
03266 }
03267
03268
03269
03270
03271
03272 if(OrganizeFramelist2[i] <= 0)
03273 {
03274
03275 if(TempIntNum > 2)
03276 {
03277 k=FiltergramLocation-1;
03278 for(ii=1;ii<=TempIntNum/2-1;++ii)
03279 {
03280 while(k > 0)
03281 {
03282 if(KeywordMissing[k] != 1)
03283 {
03284 if(HWL1POS[k] != PHWPLPOS[i*7] || HWL2POS[k] != PHWPLPOS[i*7+1] || HWL3POS[k] != PHWPLPOS[i*7+2] || HWL4POS[k] != PHWPLPOS[i*7+3] || HPL1POS[k] != PHWPLPOS[i*7+4] || HPL2POS[k] != PHWPLPOS[i*7+5] || HPL3POS[k] != PHWPLPOS[i*7+6] || HCAMID[k] != camera || CFINDEX[k] != TargetCFINDEX)
03285 {
03286 if ((internTOBS[k]-TargetTime) >= -MaxSearchDistanceL) --k;
03287 else break;
03288 }
03289 else break;
03290 }
03291 else --k;
03292 }
03293 if(k < 0)
03294 {
03295 FramelistArray[i+framelistSize*ii]=-1;
03296 continue;
03297 }
03298 if(KeywordMissing[k] != 1)
03299 {
03300 if(HWL1POS[k] == PHWPLPOS[i*7] && HWL2POS[k] == PHWPLPOS[i*7+1] && HWL3POS[k] == PHWPLPOS[i*7+2] && HWL4POS[k] == PHWPLPOS[i*7+3] && HPL1POS[k] == PHWPLPOS[i*7+4] && HPL2POS[k] == PHWPLPOS[i*7+5] && HPL3POS[k] == PHWPLPOS[i*7+6] && HCAMID[k] == camera && CFINDEX[k] == TargetCFINDEX && (internTOBS[k]-TargetTime) >= -MaxSearchDistanceL)
03301 {
03302 FramelistArray[i+framelistSize*ii]=k;
03303 }
03304 else FramelistArray[i+framelistSize*ii]=-1;
03305 }
03306 else FramelistArray[i+framelistSize*ii]=-1;
03307 --k;
03308 }
03309 }
03310 k=FiltergramLocation+1;
03311 for(ii=TempIntNum/2;ii<TempIntNum;++ii)
03312 {
03313 while(k < nRecs1-1)
03314 {
03315 if(KeywordMissing[k] != 1)
03316 {
03317 if(HWL1POS[k] != PHWPLPOS[i*7] || HWL2POS[k] != PHWPLPOS[i*7+1] || HWL3POS[k] != PHWPLPOS[i*7+2] || HWL4POS[k] != PHWPLPOS[i*7+3] || HPL1POS[k] != PHWPLPOS[i*7+4] || HPL2POS[k] != PHWPLPOS[i*7+5] || HPL3POS[k] != PHWPLPOS[i*7+6] || HCAMID[k] != camera || CFINDEX[k] != TargetCFINDEX)
03318 {
03319 if ((internTOBS[k]-TargetTime) <= MaxSearchDistanceR) ++k;
03320 else break;
03321 }
03322 else break;
03323 }
03324 else ++k;
03325 }
03326 if(k > nRecs1-1)
03327 {
03328 FramelistArray[i+framelistSize*ii]=-1;
03329 continue;
03330 }
03331 if(KeywordMissing[k] != 1)
03332 {
03333 if(HWL1POS[k] == PHWPLPOS[i*7] && HWL2POS[k] == PHWPLPOS[i*7+1] && HWL3POS[k] == PHWPLPOS[i*7+2] && HWL4POS[k] == PHWPLPOS[i*7+3] && HPL1POS[k] == PHWPLPOS[i*7+4] && HPL2POS[k] == PHWPLPOS[i*7+5] && HPL3POS[k] == PHWPLPOS[i*7+6] && HCAMID[k] == camera && CFINDEX[k] == TargetCFINDEX && (internTOBS[k]-TargetTime) <= MaxSearchDistanceR)
03334 {
03335 FramelistArray[i+framelistSize*ii]=k;
03336 }
03337 else FramelistArray[i+framelistSize*ii]=-1;
03338 }
03339 else FramelistArray[i+framelistSize*ii]=-1;
03340 ++k;
03341 }
03342 }
03343 else
03344 {
03345
03346 k=FiltergramLocation-1;
03347 for(ii=1;ii<=TempIntNum/2;++ii)
03348 {
03349 while(k > 0)
03350 {
03351 if(KeywordMissing[k] != 1)
03352 {
03353 if(HWL1POS[k] != PHWPLPOS[i*7] || HWL2POS[k] != PHWPLPOS[i*7+1] || HWL3POS[k] != PHWPLPOS[i*7+2] || HWL4POS[k] != PHWPLPOS[i*7+3] || HPL1POS[k] != PHWPLPOS[i*7+4] || HPL2POS[k] != PHWPLPOS[i*7+5] || HPL3POS[k] != PHWPLPOS[i*7+6] || HCAMID[k] != camera || CFINDEX[k] != TargetCFINDEX)
03354 {
03355 if ((internTOBS[k]-TargetTime) >= -MaxSearchDistanceR) --k;
03356 else break;
03357 }
03358 else break;
03359 }
03360 else --k;
03361 }
03362 if(k < 0)
03363 {
03364 FramelistArray[i+framelistSize*ii]=-1;
03365 continue;
03366 }
03367 if(KeywordMissing[k] != 1)
03368 {
03369 if(HWL1POS[k] == PHWPLPOS[i*7] && HWL2POS[k] == PHWPLPOS[i*7+1] && HWL3POS[k] == PHWPLPOS[i*7+2] && HWL4POS[k] == PHWPLPOS[i*7+3] && HPL1POS[k] == PHWPLPOS[i*7+4] && HPL2POS[k] == PHWPLPOS[i*7+5] && HPL3POS[k] == PHWPLPOS[i*7+6] && HCAMID[k] == camera && CFINDEX[k] == TargetCFINDEX && (internTOBS[k]-TargetTime) >= -MaxSearchDistanceR)
03370 {
03371 FramelistArray[i+framelistSize*ii]=k;
03372 }
03373 else FramelistArray[i+framelistSize*ii]=-1;
03374 }
03375 else FramelistArray[i+framelistSize*ii]=-1;
03376 --k;
03377 }
03378 if( TempIntNum > 2)
03379 {
03380 k=FiltergramLocation+1;
03381 for(ii=TempIntNum/2+1;ii<TempIntNum;++ii)
03382 {
03383 while(k < nRecs1-1)
03384 {
03385 if(KeywordMissing[k] != 1)
03386 {
03387 if(HWL1POS[k] != PHWPLPOS[i*7] || HWL2POS[k] != PHWPLPOS[i*7+1] || HWL3POS[k] != PHWPLPOS[i*7+2] || HWL4POS[k] != PHWPLPOS[i*7+3] || HPL1POS[k] != PHWPLPOS[i*7+4] || HPL2POS[k] != PHWPLPOS[i*7+5] || HPL3POS[k] != PHWPLPOS[i*7+6] || HCAMID[k] != camera || CFINDEX[k] != TargetCFINDEX)
03388 {
03389 if ((internTOBS[k]-TargetTime) <= MaxSearchDistanceL) ++k;
03390 else break;
03391 }
03392 else break;
03393 }
03394 else ++k;
03395 }
03396 if(k > nRecs1-1)
03397 {
03398 FramelistArray[i+framelistSize*ii]=-1;
03399 continue;
03400 }
03401 if(KeywordMissing[k] != 1)
03402 {
03403 if(HWL1POS[k] == PHWPLPOS[i*7] && HWL2POS[k] == PHWPLPOS[i*7+1] && HWL3POS[k] == PHWPLPOS[i*7+2] && HWL4POS[k] == PHWPLPOS[i*7+3] && HPL1POS[k] == PHWPLPOS[i*7+4] && HPL2POS[k] == PHWPLPOS[i*7+5] && HPL3POS[k] == PHWPLPOS[i*7+6] && HCAMID[k] == camera && CFINDEX[k] == TargetCFINDEX && (internTOBS[k]-TargetTime) <= MaxSearchDistanceL)
03404 {
03405 FramelistArray[i+framelistSize*ii]=k;
03406 }
03407 else FramelistArray[i+framelistSize*ii]=-1;
03408 }
03409 else FramelistArray[i+framelistSize*ii]=-1;
03410 ++k;
03411 }
03412 }
03413 }
03414
03415 }
03416
03417
03418 free(OrganizeFramelist);
03419 free(OrganizeFramelist2);
03420 OrganizeFramelist=NULL;
03421 OrganizeFramelist2=NULL;
03422
03423
03424
03425
03426
03427 printf("LOOKING FOR FILTERGRAMS ALREADY IN MEMORY BUT THAT ARE NOT NEEDED ANYMORE\n");
03428 for(ii=0;ii<nRecs1;++ii)
03429 {
03430 if(SegmentRead[ii] == 1)
03431 {
03432 Needed=0;
03433
03434 for(k=0;k<framelistSize*TempIntNum;++k)
03435 {
03436 if (FramelistArray[k] == ii)
03437 {
03438 Needed=1;
03439 break;
03440 }
03441 }
03442
03443 if(Needed == 0)
03444 {
03445 drms_free_array(Segments[ii]);
03446 drms_free_array(Ierror[ii]);
03447 Segments[ii]= NULL;
03448 Ierror[ii] = NULL;
03449 SegmentRead[ii] = 0;
03450 }
03451
03452 }
03453
03454 }
03455
03456
03457
03458
03459
03460
03461
03462 arrLev1d = (DRMS_Array_t **)malloc(nRecs1d*sizeof(DRMS_Array_t *));
03463 if(arrLev1d == NULL)
03464 {
03465 printf("Error: memory could not be allocated to arrLev1d\n");
03466 return 1;
03467 }
03468
03469 for(k=0;k<nRecs1d;++k)
03470 {
03471 arrLev1d[k]= drms_array_create(type1d,2,axisout,NULL,&status);
03472 if(status != DRMS_SUCCESS || arrLev1d[k] == NULL)
03473 {
03474 printf("Error: cannot create a DRMS array for a level 1d filtergram with index %d at target time %s\n",k,timeBegin2);
03475
03476
03477
03478
03479 return 1;
03480 }
03481 }
03482 Segments1d=1;
03483
03484
03485
03486
03487
03488 printf("COMPUTING MEDIAN VALUE OF X0, Y0, AND RSUN\n");
03489 X0ARR =(float *)malloc(framelistSize*TempIntNum*sizeof(float));
03490 Y0ARR =(float *)malloc(framelistSize*TempIntNum*sizeof(float));
03491 RSUNARR =(float *)malloc(framelistSize*TempIntNum*sizeof(float));
03492
03493 for(k=0;k<framelistSize;++k)
03494 {
03495 for(i=0;i<TempIntNum;++i)
03496 {
03497 temp=FramelistArray[k+framelistSize*i];
03498 if(temp != -1)
03499 {
03500 X0ARR[k+framelistSize*i]=X0[temp];
03501 Y0ARR[k+framelistSize*i]=Y0[temp];
03502 RSUNARR[k+framelistSize*i]=RSUN[temp];
03503 }
03504 else
03505 {
03506 X0ARR[k+framelistSize*i] =MISSINGRESULT;
03507 Y0ARR[k+framelistSize*i] =MISSINGRESULT;
03508 RSUNARR[k+framelistSize*i]=MISSINGRESULT;
03509 }
03510
03511 }
03512 }
03513
03514
03515 status=fstats(framelistSize*TempIntNum,X0ARR,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
03516 printf("VALUES USED TO OBTAIN THE MEDIAN: %d\n",ngood);
03517 X0AVG=median;
03518 X0RMS=sigma;
03519 status=fstats(framelistSize*TempIntNum,Y0ARR,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
03520 Y0AVG=median;
03521 Y0RMS=sigma;
03522 status=fstats(framelistSize*TempIntNum,RSUNARR,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
03523 RSUNAVG=median;
03524 RSUNRMS=sigma;
03525
03526 free(X0ARR);
03527 free(Y0ARR);
03528 free(RSUNARR);
03529 X0ARR=NULL;
03530 Y0ARR=NULL;
03531 RSUNARR=NULL;
03532
03533 if(combine == 1)
03534 {
03535 printf("COMPUTING MEDIAN VALUE OF X0, Y0, AND RSUN, FOR FRONT CAMERA ONLY\n");
03536 X0ARR =(float *)malloc(framelistSize*TempIntNum*sizeof(float));
03537 Y0ARR =(float *)malloc(framelistSize*TempIntNum*sizeof(float));
03538 RSUNARR =(float *)malloc(framelistSize*TempIntNum*sizeof(float));
03539
03540 for(k=0;k<framelistSize;++k)
03541 {
03542 for(i=0;i<TempIntNum;++i)
03543 {
03544 temp=FramelistArray[k+framelistSize*i];
03545 if(temp != -1)
03546 {
03547 if(HCAMID[temp] == LIGHT_FRONT)
03548 {
03549 X0ARR[k+framelistSize*i]=X0[temp];
03550 Y0ARR[k+framelistSize*i]=Y0[temp];
03551 RSUNARR[k+framelistSize*i]=RSUN[temp];
03552 }
03553 else
03554 {
03555 X0ARR[k+framelistSize*i] =MISSINGRESULT;
03556 Y0ARR[k+framelistSize*i] =MISSINGRESULT;
03557 RSUNARR[k+framelistSize*i]=MISSINGRESULT;
03558 }
03559 }
03560 else
03561 {
03562 X0ARR[k+framelistSize*i] =MISSINGRESULT;
03563 Y0ARR[k+framelistSize*i] =MISSINGRESULT;
03564 RSUNARR[k+framelistSize*i]=MISSINGRESULT;
03565 }
03566
03567 }
03568 }
03569
03570
03571 status=fstats(framelistSize*TempIntNum,X0ARR,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
03572 printf("VALUES USED TO OBTAIN THE MEDIAN (FRONT CAMERA): %d\n",ngood);
03573 X0AVGF=median;
03574 status=fstats(framelistSize*TempIntNum,Y0ARR,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
03575 Y0AVGF=median;
03576 status=fstats(framelistSize*TempIntNum,RSUNARR,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
03577 RSUNAVGF=median;
03578
03579 free(X0ARR);
03580 free(Y0ARR);
03581 free(RSUNARR);
03582 X0ARR=NULL;
03583 Y0ARR=NULL;
03584 RSUNARR=NULL;
03585
03586 printf("COMPUTING MEDIAN VALUE OF X0, Y0, AND RSUN, FOR SIDE CAMERA ONLY\n");
03587 X0ARR =(float *)malloc(framelistSize*TempIntNum*sizeof(float));
03588 Y0ARR =(float *)malloc(framelistSize*TempIntNum*sizeof(float));
03589 RSUNARR =(float *)malloc(framelistSize*TempIntNum*sizeof(float));
03590
03591 for(k=0;k<framelistSize;++k)
03592 {
03593 for(i=0;i<TempIntNum;++i)
03594 {
03595 temp=FramelistArray[k+framelistSize*i];
03596 if(temp != -1)
03597 {
03598 if(HCAMID[temp] == LIGHT_SIDE)
03599 {
03600 X0ARR[k+framelistSize*i]=X0[temp];
03601 Y0ARR[k+framelistSize*i]=Y0[temp];
03602 RSUNARR[k+framelistSize*i]=RSUN[temp];
03603 }
03604 else
03605 {
03606 X0ARR[k+framelistSize*i] =MISSINGRESULT;
03607 Y0ARR[k+framelistSize*i] =MISSINGRESULT;
03608 RSUNARR[k+framelistSize*i]=MISSINGRESULT;
03609 }
03610 }
03611 else
03612 {
03613 X0ARR[k+framelistSize*i] =MISSINGRESULT;
03614 Y0ARR[k+framelistSize*i] =MISSINGRESULT;
03615 RSUNARR[k+framelistSize*i]=MISSINGRESULT;
03616 }
03617
03618 }
03619 }
03620
03621
03622 status=fstats(framelistSize*TempIntNum,X0ARR,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
03623 printf("VALUES USED TO OBTAIN THE MEDIAN (SIDE CAMERA): %d\n",ngood);
03624 X0AVGS=median;
03625 status=fstats(framelistSize*TempIntNum,Y0ARR,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
03626 Y0AVGS=median;
03627 status=fstats(framelistSize*TempIntNum,RSUNARR,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
03628 RSUNAVGS=median;
03629
03630 free(X0ARR);
03631 free(Y0ARR);
03632 free(RSUNARR);
03633 X0ARR=NULL;
03634 Y0ARR=NULL;
03635 RSUNARR=NULL;
03636
03637 }
03638
03639
03640
03641 for(k=0;k<framelistSize;++k)
03642 {
03643
03644 ActualTempIntNum=TempIntNum;
03645
03646
03647
03648
03649 for(i=0;i<TempIntNum;++i)
03650 {
03651 temp=FramelistArray[k+framelistSize*i];
03652 if(temp != -1)
03653 {
03654
03655
03656 if(SegmentRead[temp] != -1)
03657 {
03658
03659 if(SegmentRead[temp] == 0)
03660 {
03661 printf("segment needs to be read for FSN %d %d\n",FSN[temp],HCAMID[temp]);
03662 segin = drms_segment_lookupnum(recLev1->records[temp], 0);
03663 Segments[temp] = drms_segment_read(segin,type1d, &status);
03664 if (status != DRMS_SUCCESS || Segments[temp] == NULL)
03665 {
03666 printf("Error: could not read the segment of level 1 record FSN = %d at target time %s\n",FSN[temp],timeBegin2);
03667 ActualTempIntNum-=1;
03668 arrin[i] = NULL;
03669 arrerrors[i] = NULL;
03670 Segments[temp] = NULL;
03671 Ierror[temp] = NULL;
03672 SegmentRead[temp]=-1;
03673 }
03674 else
03675 {
03676 Ierror[temp] = drms_array_create(typeEr,2,axisout,NULL,&status);
03677 if(status != DRMS_SUCCESS || Ierror[temp] == NULL)
03678 {
03679 printf("Error: could not create an array for Ierror at target time %s\n",timeBegin2);
03680 drms_free_array(Segments[temp]);
03681 Segments[temp]=NULL;
03682 Ierror[temp]=NULL;
03683 SegmentRead[temp]=-1;
03684 ActualTempIntNum-=1;
03685 arrin[i] = NULL;
03686 arrerrors[i] = NULL;
03687 }
03688 else
03689 {
03690 arrin[i] = Segments[temp];
03691 arrerrors[i] = Ierror[temp];
03692 if(arrin[i]->axis[0] != axisin[0] || arrin[i]->axis[1] != axisin[1])
03693 {
03694 printf("Error: level 1 record FSN = %d at target time %s has a segment with dimensions %d x %d instead of %d x %d\n",FSN[temp],timeBegin2,arrin[i]->axis[0],arrin[i]->axis[1],axisin[0],axisin[1]);
03695 drms_free_array(Segments[temp]);
03696 drms_free_array(Ierror[temp]);
03697 ActualTempIntNum-=1;
03698 arrin[i] = NULL;
03699 arrerrors[i] = NULL;
03700 Segments[temp]=NULL;
03701 Ierror[temp]=NULL;
03702 SegmentRead[temp]=-1;
03703 }
03704 else
03705 {
03706 SegmentRead[temp]=1;
03707
03708
03709
03710
03711 segin = drms_segment_lookup(recLev1->records[temp],"bad_pixel_list");
03712 BadPixels = NULL;
03713 BadPixels = drms_segment_read(segin,segin->info->type,&status);
03714 if(status != DRMS_SUCCESS || BadPixels == NULL)
03715 {
03716 printf("Error: cannot read the list of bad pixels of level 1 filtergram FSN= %d\n",FSN[temp]);
03717 ActualTempIntNum-=1;
03718 drms_free_array(Segments[temp]);
03719 drms_free_array(Ierror[temp]);
03720 arrin[i] = NULL;
03721 arrerrors[i] = NULL;
03722 Segments[temp]=NULL;
03723 Ierror[temp]=NULL;
03724 SegmentRead[temp]=-1;
03725 }
03726 else
03727 {
03728
03729 strcpy(HMISeriesTemp,CosmicRaySeries);
03730 strcat(HMISeriesTemp,"[][");
03731 sprintf(FSNtemps,"%d",FSN[temp]);
03732 strcat(HMISeriesTemp,FSNtemps);
03733 strcat(HMISeriesTemp,"]");
03734 rectemp=NULL;
03735
03736 rectemp=drms_open_records(drms_env,HMISeriesTemp,&statusA[0]);
03737
03738 if(statusA[0] == DRMS_SUCCESS && rectemp != NULL && rectemp->n != 0)
03739 {
03740 segin = drms_segment_lookupnum(rectemp->records[0],0);
03741 CosmicRays = NULL;
03742
03743 CosmicRays = drms_segment_read(segin,segin->info->type,&status);
03744 if(status != DRMS_SUCCESS || CosmicRays == NULL)
03745 {
03746 printf("Error: the list of cosmic-ray hits could not be read for FSN %d\n",FSN[temp]);
03747
03748
03749
03750 QUALITY = QUALITY | QUAL_NOCOSMICRAY;
03751 QUALITYlev1[temp] = QUALITYlev1[temp] | QUAL_NOCOSMICRAY;
03752 CosmicRays = NULL;
03753
03754 }
03755
03756 COSMICCOUNT=drms_getkey_int(rectemp->records[0],COUNTS,&status);
03757 if(status != DRMS_SUCCESS || COSMICCOUNT == -1)
03758 {
03759 QUALITY = QUALITY | QUAL_NOCOSMICRAY;
03760 QUALITYlev1[temp] = QUALITYlev1[temp] | QUAL_NOCOSMICRAY;
03761 }
03762 }
03763 else
03764 {
03765 printf("Unable to open the series %s for FSN %d\n",HMISeriesTemp,FSN[temp]);
03766
03767
03768
03769 QUALITY = QUALITY | QUAL_NOCOSMICRAY;
03770 QUALITYlev1[temp] = QUALITYlev1[temp] | QUAL_NOCOSMICRAY;
03771 CosmicRays = NULL;
03772
03773 }
03774
03775 image = Segments[temp]->data;
03776
03777
03778
03779
03780
03781 if(inRotationalFlat == 1)
03782 {
03783
03784 HMIFlatField = (char *)malloc(MaxNString*sizeof(char *));
03785 if(HMIFlatField == NULL)
03786 {
03787 printf("Error: memory could not be allocated to HMIFlatField\n");
03788 return 1;
03789 }
03790 HMIFlatField = drms_getkey_string(recLev1->records[temp],FLATREC,&status);
03791 if (status != DRMS_SUCCESS)
03792 {
03793 printf("Error: could not read the FLAT_REC keyword for the target filtergram FSN= %d",FSN[temp]);
03794 return 1;
03795 }
03796
03797
03798 if(strcmp(HMIFlatField,HMIFlatField0) != 0)
03799 {
03800 printf("Warning: the hmi.flatfield record used to produce the level 1 records changed during the run of the observables code\n");
03801
03802 recflat = drms_open_records(drms_env,HMIFlatField,&status);
03803 if (status != DRMS_SUCCESS || recflat == NULL || recflat->n == 0)
03804 {
03805 printf("Error: record missing or corrupt for the flat field query %s\n",HMIFlatField);
03806 return 1;
03807 }
03808 drms_free_array(flatfield);
03809 strcpy(HMIFlatField0,HMIFlatField);
03810 segin = drms_segment_lookup(recflat->records[0],"flatfield");
03811 flatfield = drms_segment_read(segin,type1d,&status);
03812 if (status != DRMS_SUCCESS || flatfield == NULL)
03813 {
03814 printf("Error: could not read the data segment for the flat field query %s\n",HMIFlatField);
03815 return 1;
03816 }
03817 pztflat = flatfield->data;
03818 status=drms_close_records(recflat,DRMS_FREE_RECORD);
03819 }
03820
03821 printf("applying rotational flat field on record FSN=%d\n",FSN[temp]);
03822
03823 for(iii=0;iii<axisin[0]*axisin[1];++iii) image[iii]=image[iii]*pztflat[iii];
03824
03825 for(iii=0;iii<axisin[0]*axisin[1];++iii) image[iii]=image[iii]/rotflat[iii];
03826
03827 free(HMIFlatField);
03828
03829 }
03830
03831
03832
03833 status = MaskCreation(Mask,axisin[0],axisin[1],BadPixels,HIMGCFID[temp],image,CosmicRays,NBADPERM[temp]);
03834 if(status != 0)
03835 {
03836 printf("Error: unable to create a mask for the gap filling function\n");
03837 return 1;
03838 }
03839 if(BadPixels != NULL)
03840 {
03841 drms_free_array(BadPixels);
03842 BadPixels=NULL;
03843 }
03844 if(CosmicRays != NULL)
03845 {
03846 drms_free_array(CosmicRays);
03847 CosmicRays=NULL;
03848 }
03849 if(rectemp != NULL)
03850 {
03851 drms_close_records(rectemp,DRMS_FREE_RECORD);
03852 rectemp=NULL;
03853 }
03854
03855 image = arrin[i]->data;
03856 ierror = arrerrors[i]->data;
03857
03858 t0=dsecnd();
03859 printf("STARTING GAPFILL\n");
03860 status =do_gapfill(image,Mask,&const_param,ierror,axisin[0],axisin[1]);
03861 t1=dsecnd();
03862 printf("TIME ELAPSED TO GAPFILL: %f\n",t1-t0);
03863 if(status != 0)
03864 {
03865 printf("Error: gapfilling code did not work on level 1 filtergram FSN = %d at target time %s\n",FSN[temp],timeBegin2);
03866 QUALITY = QUALITY | QUAL_NOGAPFILL;
03867 QUALITYlev1[temp] = QUALITYlev1[temp] | QUAL_NOGAPFILL;
03868 }
03869 }
03870 }
03871 }
03872 }
03873 }
03874 else
03875 {
03876 printf("segment is already in memory for FSN %d %d\n",FSN[temp],HCAMID[temp]);
03877 arrin[i] = Segments[temp];
03878 arrerrors[i] = Ierror[temp];
03879 }
03880
03881 }
03882 else
03883 {
03884 printf("Error: the filtergram FSN = %d has corrupted/missing keyword(s), or a corrupted/missing data segment, at target time %s\n",FSN[temp],timeBegin2);
03885 ActualTempIntNum-=1;
03886 arrin[i] = NULL;
03887 arrerrors[i] = NULL;
03888 }
03889 }
03890 else
03891 {
03892 printf("a filtergram close to the target time is missing/corrupted\n");
03893 ActualTempIntNum-=1;
03894 arrin[i] = NULL;
03895 arrerrors[i] = NULL;
03896 }
03897 }
03898
03899
03900
03901
03902 if(QuickLook == 1) ThresholdPol=1;
03903 else ThresholdPol=2;
03904 if(ActualTempIntNum >= ThresholdPol)
03905 {
03906
03907 images = (float **)malloc(ActualTempIntNum*sizeof(float *));
03908 if(images == NULL)
03909 {
03910 printf("Error: memory could not be allocated to images\n");
03911 return 1;
03912 }
03913
03914 ierrors = (char **)malloc(ActualTempIntNum*sizeof(char *));
03915 if(ierrors == NULL)
03916 {
03917 printf("Error: memory could not be allocated to ierrors\n");
03918 return 1;
03919 }
03920
03921
03922 KeyInterp = (struct keyword *)malloc(ActualTempIntNum*sizeof(struct keyword));
03923 if(KeyInterp == NULL)
03924 {
03925 printf("Error: memory could not be allocated to KeyInterp\n");
03926 return 1;
03927 }
03928
03929
03930
03931 ii=0;
03932 for(i=0;i<TempIntNum;++i) if (arrin[i] != NULL && arrerrors[i] != NULL)
03933 {
03934 temp=FramelistArray[k+framelistSize*i];
03935 if(fabs(RSUN[temp]-RSUNAVG) > 1.82*RSUNerr || isnan(RSUN[temp]))
03936 {
03937 printf("Warning: image %d passed to do_interpolate has a RSUN value too different from the median value and will be rejected: %f %f\n",FSN[temp],RSUN[temp],RSUNAVG);
03938 ActualTempIntNum-=1;
03939 QUALITY = QUALITY | QUAL_LIMBFITISSUE;
03940 continue;
03941 }
03942 if(combine == 0)
03943 {
03944 if(fabs(X0[temp]-X0AVG) > RSUNerr || isnan(X0[temp]))
03945 {
03946 printf("Warning: image %d passed to do_interpolate has a X0 value too different from the median value and will be rejected: %f %f\n",FSN[temp],X0[temp],X0AVG);
03947 ActualTempIntNum-=1;
03948 QUALITY = QUALITY | QUAL_LIMBFITISSUE;
03949 continue;
03950 }
03951 if(fabs(Y0[temp]-Y0AVG) > RSUNerr || isnan(Y0[temp]))
03952 {
03953 printf("Warning: image %d passed to do_interpolate has a Y0 value too different from the median value and will be rejected: %f %f\n",FSN[temp],Y0[temp],Y0AVG);
03954 ActualTempIntNum-=1;
03955 QUALITY = QUALITY | QUAL_LIMBFITISSUE;
03956 continue;
03957 }
03958 }
03959 else
03960 {
03961 if(HCAMID[temp] == LIGHT_FRONT)
03962 {
03963 if(fabs(X0[temp]-X0AVGF) > RSUNerr || isnan(X0[temp]))
03964 {
03965 printf("Warning: image %d passed to do_interpolate has a X0 value too different from the median value and will be rejected: %f %f\n",FSN[temp],X0[temp],X0AVG);
03966 ActualTempIntNum-=1;
03967 QUALITY = QUALITY | QUAL_LIMBFITISSUE;
03968 continue;
03969 }
03970 if(fabs(Y0[temp]-Y0AVGF) > RSUNerr || isnan(Y0[temp]))
03971 {
03972 printf("Warning: image %d passed to do_interpolate has a Y0 value too different from the median value and will be rejected: %f %f\n",FSN[temp],Y0[temp],Y0AVG);
03973 ActualTempIntNum-=1;
03974 QUALITY = QUALITY | QUAL_LIMBFITISSUE;
03975 continue;
03976 }
03977 }
03978 else
03979 {
03980 if(fabs(X0[temp]-X0AVGS) > RSUNerr || isnan(X0[temp]))
03981 {
03982 printf("Warning: image %d passed to do_interpolate has a X0 value too different from the median value and will be rejected: %f %f\n",FSN[temp],X0[temp],X0AVG);
03983 ActualTempIntNum-=1;
03984 QUALITY = QUALITY | QUAL_LIMBFITISSUE;
03985 continue;
03986 }
03987 if(fabs(Y0[temp]-Y0AVGS) > RSUNerr || isnan(Y0[temp]))
03988 {
03989 printf("Warning: image %d passed to do_interpolate has a Y0 value too different from the median value and will be rejected: %f %f\n",FSN[temp],Y0[temp],Y0AVG);
03990 ActualTempIntNum-=1;
03991 QUALITY = QUALITY | QUAL_LIMBFITISSUE;
03992 continue;
03993 }
03994 }
03995 }
03996 images[ii]=arrin[i]->data;
03997 ierrors[ii]=arrerrors[i]->data;
03998 if(HCAMID[temp] == LIGHT_FRONT) KeyInterp[ii].camera=0;
03999 else KeyInterp[ii].camera=1;
04000 if(!strcmp(HWLTNSET[temp],"OPEN")) QUALITY = QUALITY | QUAL_ISSTARGET;
04001 if( (QUALITYin[temp] & Q_ACS_ECLP) == Q_ACS_ECLP) QUALITY = QUALITY | QUAL_ECLIPSE;
04002 KeyInterp[ii].rsun=RSUNAVG;
04003 printf("actual solar radius of image in = %f\n",RSUN[temp]);
04004 KeyInterp[ii].xx0=X0[temp];
04005 KeyInterp[ii].yy0=Y0[temp];
04006 KeyInterp[ii].dist=(float)(DSUNOBS[temp]/(double)AstroUnit);
04007 KeyInterp[ii].b0=CRLTOBS[temp]/180.*M_PI;
04008 KeyInterp[ii].p0=CROTA2[temp]/180.*M_PI;
04009 KeyInterp[ii].time=internTOBS[temp];
04010 KeyInterp[ii].focus=CFINDEX[temp];
04011
04012 rec=recLev1->records[temp];
04013 sprintf(recnums,"%ld",rec->recnum);
04014 if(k != 0 || ii != 0) strcat(source,",#");
04015 else strcat(source,"#");
04016 strcat(source,recnums);
04017
04018 QUALITYLEV1 = QUALITYLEV1 | QUALITYin[temp];
04019 QUALITY = QUALITY | QUALITYlev1[temp];
04020
04021 ii+=1;
04022 }
04023
04024
04025
04026 RSUNint=RSUNAVG;
04027 KeyInterpOut.rsun=RSUNint;
04028
04029 if(combine == 0)
04030 {
04031 KeyInterpOut.xx0 =X0AVG;
04032 KeyInterpOut.yy0 =Y0AVG;
04033 }
04034 else
04035 {
04036 if(CamId == LIGHT_FRONT)
04037 {
04038 KeyInterpOut.xx0 =X0AVGF;
04039 KeyInterpOut.yy0 =Y0AVGF;
04040 }
04041 else
04042 {
04043 KeyInterpOut.xx0 =X0AVGS;
04044 KeyInterpOut.yy0 =Y0AVGS;
04045 }
04046 }
04047
04048 KeyInterpOut.dist=(float)DSUNOBSint;
04049 KeyInterpOut.b0 =CRLTOBSint/180.*M_PI;
04050 KeyInterpOut.p0 =CROTA2int/180.*M_PI;
04051 tobs = TargetTime+(DSUNOBSint-1.0)/2.00398880422056639358e-03;
04052 KeyInterpOut.time=tobs;
04053 KeyInterpOut.focus=TargetCFINDEX;
04054
04055
04056
04057 printf("Calling temporal interpolation, de-rotation, and un-distortion code\n");
04058
04059
04060 printf("ACTUALTEMPINT %d\n",ActualTempIntNum);
04061 if(ActualTempIntNum != TempIntNum) QUALITY = QUALITY | QUAL_LOWINTERPNUM;
04062 else
04063 {
04064 minimum=KeyInterp[0].time;
04065 maximum=KeyInterp[0].time;
04066 for(ii=1;ii<TempIntNum;++ii)
04067 {
04068 if(KeyInterp[ii].time < minimum) minimum=KeyInterp[ii].time;
04069 if(KeyInterp[ii].time > maximum) maximum=KeyInterp[ii].time;
04070 }
04071 if((maximum-minimum) > (DataCadence*(double)(TempIntNum-1)+DataCadence/(double)framelistSize) ) QUALITY = QUALITY | QUAL_LOWINTERPNUM;
04072 }
04073
04074 if(ActualTempIntNum >= ThresholdPol)
04075 {
04076 for(ii=0;ii<ActualTempIntNum;++ii) printf("KEYWORDS IN: %f %f %f %f %f %f %f %d\n",KeyInterp[ii].rsun,KeyInterp[ii].xx0,KeyInterp[ii].yy0,KeyInterp[ii].dist,KeyInterp[ii].b0,KeyInterp[ii].p0,KeyInterp[ii].time,KeyInterp[ii].focus);
04077 printf("KEYWORDS OUT: %f %f %f %f %f %f %f %d\n",KeyInterpOut.rsun,KeyInterpOut.xx0,KeyInterpOut.yy0,KeyInterpOut.dist,KeyInterpOut.b0,KeyInterpOut.p0,KeyInterpOut.time,KeyInterpOut.focus);
04078
04079 t0=dsecnd();
04080 status=do_interpolate(images,ierrors,arrLev1d[k]->data,KeyInterp,&KeyInterpOut,&const_param,ActualTempIntNum,axisin[0],axisin[1],-1.0);
04081 t1=dsecnd();
04082 printf("TIME ELAPSED IN DO_INTERPOLATE: %f\n",t1-t0);
04083 }
04084 else
04085 {
04086 printf("Error: ActualTempIntNum <ThresholdPol\n");
04087 QUALITY = QUALITY | QUAL_NOTENOUGHINTERPOLANTS;
04088 status = 1;
04089 }
04090
04091 printf("End temporal interpolation, de-rotation, and un-distortion\n");
04092 if (status != 0)
04093 {
04094 printf("Error: temporal interpolation, de-rotation, and undistortion subroutine failed at target time %s\n",timeBegin2);
04095 drms_free_array(arrLev1d[k]);
04096 arrLev1d[k] = NULL;
04097 QUALITY = QUALITY | QUAL_INTERPOLATIONFAILED;
04098 }
04099 else
04100 {
04101 for(i=0;i<TempIntNum;++i)
04102 {
04103 temp=FramelistArray[k+framelistSize*i];
04104 if(temp != -1) break;
04105 }
04106
04107
04108
04109
04110 if(combine == 1)camera=3;
04111 if(combine == 0 && CamId == LIGHT_SIDE) camera=1;
04112 if(combine == 0 && CamId == LIGHT_FRONT) camera=2;
04113
04114
04115 cdelt1=1.0/RSUNint*asin(solar_radius/(DSUNOBSint*(double)AstroUnit))*180.*60.*60./M_PI;
04116
04117 statusA[0] = drms_setkey_time(recLev1d->records[k],TRECS,trec);
04118 statusA[1] = drms_setkey_time(recLev1d->records[k],TOBSS,tobs);
04119 statusA[2] = drms_setkey_int(recLev1d->records[k],CAMERAS,camera);
04120 statusA[3] = drms_setkey_string(recLev1d->records[k],HISTORYS,HISTORY);
04121 statusA[4] = drms_setkey_string(recLev1d->records[k],COMMENTS,COMMENT);
04122
04123 statusA[5]=0;
04124 statusA[6] = drms_setkey_int(recLev1d->records[k],HFLIDS,TargetHFLID);
04125 statusA[7] = drms_setkey_int(recLev1d->records[k],HCFTIDS,TargetCFINDEX);
04126
04127 statusA[8]=0;
04128 statusA[9] = drms_setkey_int(recLev1d->records[k],QUALITYS,QUALITY);
04129 statusA[10]= drms_setkey_double(recLev1d->records[k],DSUNOBSS,DSUNOBSint*(double)AstroUnit);
04130 statusA[11]= drms_setkey_float(recLev1d->records[k],CRLTOBSS,KeyInterpOut.b0*180./M_PI);
04131 statusA[12]= drms_setkey_float(recLev1d->records[k],CROTA2S,-KeyInterpOut.p0*180./M_PI);
04132 statusA[13]= drms_setkey_float(recLev1d->records[k],CDELT1S,cdelt1);
04133 statusA[14]= drms_setkey_float(recLev1d->records[k],CDELT2S,cdelt1);
04134 statusA[15]= drms_setkey_float(recLev1d->records[k],CRPIX1S,KeyInterpOut.xx0+1.);
04135 statusA[16]= drms_setkey_float(recLev1d->records[k],CRPIX2S,KeyInterpOut.yy0+1.);
04136 sprintf(jsocverss,"%d",jsoc_vers_num);
04137 statusA[17]= drms_setkey_string(recLev1d->records[k],BLDVERSS,jsocverss);
04138 statusA[18]= drms_setkey_double(recLev1d->records[k],OBSVRS,OBSVRint);
04139 statusA[19]= drms_setkey_double(recLev1d->records[k],OBSVWS,OBSVWint);
04140 statusA[20]= drms_setkey_double(recLev1d->records[k],OBSVNS,OBSVNint);
04141
04142 statusA[21]=0;
04143 statusA[22]= drms_setkey_float(recLev1d->records[k],CRLNOBSS,CRLNOBSint);
04144 statusA[23]= drms_setkey_int(recLev1d->records[k],CARROTS,CARROTint);
04145 statusA[24]= drms_setkey_int(recLev1d->records[k],TINTNUMS,ActualTempIntNum);
04146 statusA[25]= drms_setkey_int(recLev1d->records[k],SINTNUMS,const_param.order_int);
04147 statusA[26]= drms_setkey_string(recLev1d->records[k],CODEVER0S,CODEVERSION);
04148 statusA[27]= drms_setkey_string(recLev1d->records[k],DISTCOEFS,DISTCOEFPATH);
04149 statusA[28]= drms_setkey_string(recLev1d->records[k],ROTCOEFS,ROTCOEFPATH);
04150 statusA[29]= drms_setkey_int(recLev1d->records[k],ODICOEFFS,const_param.order_dist_coef);
04151 statusA[30]= drms_setkey_int(recLev1d->records[k],OROCOEFFS,2*const_param.order2_rot_coef);
04152 statusA[31]= drms_setkey_string(recLev1d->records[k],CODEVER1S,CODEVERSION1);
04153 statusA[32]= drms_setkey_string(recLev1d->records[k],CODEVER2S,CODEVERSION2);
04154 statusA[33]= drms_setkey_int(recLev1d->records[k],HPL1POSS,HPL1POS[temp]);
04155 statusA[34]= drms_setkey_int(recLev1d->records[k],HPL2POSS,HPL2POS[temp]);
04156 statusA[35]= drms_setkey_int(recLev1d->records[k],HPL3POSS,HPL3POS[temp]);
04157 statusA[36]= drms_setkey_int(recLev1d->records[k],HWL1POSS,HWL1POS[temp]);
04158 statusA[37]= drms_setkey_int(recLev1d->records[k],HWL2POSS,HWL2POS[temp]);
04159 statusA[38]= drms_setkey_int(recLev1d->records[k],HWL3POSS,HWL3POS[temp]);
04160 statusA[39]= drms_setkey_int(recLev1d->records[k],HWL4POSS,HWL4POS[temp]);
04161 statusA[40]= drms_setkey_int(recLev1d->records[k],FIDS,FID[temp]);
04162 printf("FID of written record: %d\n",FID[temp]);
04163 statusA[41]= drms_setkey_int(recLev1d->records[k],HCAMIDS,CamId);
04164 sprint_time(DATEOBS,tobs-DataCadence/2.0,"UTC",1);
04165 statusA[42]= drms_setkey_string(recLev1d->records[k],DATEOBSS,DATEOBS);
04166 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
04167 statusA[43]= drms_setkey_string(recLev1d->records[k],DATES,DATEOBS);
04168 statusA[44]= drms_setkey_int(recLev1d->records[k],HWLTIDS,TargetHWLTID);
04169 statusA[45]= drms_setkey_int(recLev1d->records[k],HPLTIDS,TargetHPLTID);
04170 statusA[46]= drms_setkey_int(recLev1d->records[k],WavelengthIDS,WavelengthID);
04171 if(camera == 2) strcpy(DATEOBS,"HMI_FRONT2");
04172 if(camera == 1) strcpy(DATEOBS,"HMI_SIDE1");
04173 if(camera == 3) strcpy(DATEOBS,"HMI_COMBINED");
04174 statusA[47]= drms_setkey_string(recLev1d->records[k],INSTRUMES,DATEOBS);
04175 statusA[48]= drms_setkey_string(recLev1d->records[k],ROTFLAT,QueryFlatField);
04176
04177 TotalStatus=0;
04178 for(i=0;i<49;++i) TotalStatus+=statusA[i];
04179 if(TotalStatus != 0)
04180 {
04181 for(ii=0;ii<49;++ii) if(statusA[ii] != 0) printf(" %d ",ii);
04182 printf("\n");
04183 printf("WARNING: could not set some of the keywords modified by the temporal interpolation subroutine for the level 1d data at target time %s\n",timeBegin2);
04184 }
04185
04186
04187 if (Lev1dWanted)
04188 {
04189 segout = drms_segment_lookupnum(recLev1d->records[k],0);
04190 arrLev1d[k]->bzero=segout->bzero;
04191 arrLev1d[k]->bscale=segout->bscale;
04192 arrLev1d[k]->israw=0;
04193 drms_segment_write(segout,arrLev1d[k],0);
04194 }
04195
04196 }
04197
04198 free(images);
04199 free(ierrors);
04200 images=NULL;
04201 ierrors=NULL;
04202 free(KeyInterp);
04203 KeyInterp=NULL;
04204 }
04205 else
04206 {
04207 printf("Error: not enough valid level 1 filtergrams to produce a level 1d filtergram at target time %s\n",timeBegin2);
04208 drms_free_array(arrLev1d[k]);
04209 arrLev1d[k]=NULL;
04210 QUALITY = QUALITY | QUAL_NOTENOUGHINTERPOLANTS;
04211
04212 }
04213
04214
04215 }
04216
04217
04218 free(FramelistArray);
04219 FramelistArray=NULL;
04220 }
04221
04222
04223
04224
04225
04226
04227
04228
04229
04230
04231
04232
04233
04234 if (TestLevIn[2]==1)
04235 {
04236
04237
04238
04239
04240
04241
04242
04243
04244
04245
04246 PolarizationType=1;
04247 if(CamId == LIGHT_SIDE) camera=1;
04248 if(CamId == LIGHT_FRONT) camera=2;
04249 if(camera == 2)
04250 {
04251 printf("Error: the use of time-averaged IQUV data as input requires CamId = side camera\n");
04252 return 1;
04253 }
04254
04255
04256
04257
04258 if(PolarizationType ==3 || PolarizationType ==2) strcpy(HMISeries,HMISeriesLev1pb);
04259 else strcpy(HMISeries,HMISeriesLev1pa);
04260 strcat(HMISeries,"[");
04261 strcat(HMISeries,timeBegin2);
04262 strcat(HMISeries,"][");
04263
04264 sprintf(CamIds,"%d",camera);
04265 strcat(HMISeries,CamIds);
04266 strcat(HMISeries,"]");
04267
04268
04269 recLev1p = drms_open_records(drms_env,HMISeries,&status);
04270
04271 if (status == DRMS_SUCCESS && recLev1p != NULL && recLev1p->n > 0)
04272 {
04273 if(recLev1p->n > 1)
04274 {
04275 printf("Too many records: %d\n",recLev1p->n);
04276 return 1;
04277 }
04278 printf("Number of level 1p records opened= %d\n",recLev1p->n);
04279
04280 nSegs1p = drms_record_numsegments(recLev1p->records[0]);
04281 printf("NUMBER OF DATA SEGMENTS: %d\n",nSegs1p);
04282
04283
04284 if (PolarizationType ==1)
04285 {
04286 nRecs1p = nSegs1p/4;
04287 }
04288 if (PolarizationType ==2 || PolarizationType ==3)
04289 {
04290 nRecs1p = nSegs1p/2;
04291 }
04292 printf("NUMBER OF WAVELENGTHS: %d\n",nRecs1p);
04293
04294 trec = drms_getkey_time(recLev1p->records[0],TRECS,&status);
04295 if(status != DRMS_SUCCESS)
04296 {
04297 printf("Error: cannot read the keyword %s\n",TRECS);
04298 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04299 }
04300 if(trec != TargetTime)
04301 {
04302 printf("Error: %s of a level 1p record is not equal to the target time %s\n",TRECS,timeBegin2);
04303 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04304 }
04305 tobs = drms_getkey_time(recLev1p->records[0],TOBSS,&status);
04306 if(status != DRMS_SUCCESS)
04307 {
04308 printf("Error: cannot read the keyword %s\n",TOBSS);
04309 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04310 }
04311
04312
04313
04314
04315
04316 TREC_STEP= drms_getkey_time(recLev1p->records[0],TRECSTEPS,&status);
04317 if(status != DRMS_SUCCESS)
04318 {
04319 printf("Error: cannot read the keyword %s\n",TRECSTEPS);
04320 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04321 }
04322
04323 if(TREC_STEP != DataCadence)
04324 {
04325 printf("Error: the cadence is not equal to the T_REC_step keyword of the level 1p data\n");
04326 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04327 }
04328
04329 DRMS_SegmentDimInfo_t di;
04330 segin = drms_segment_lookupnum(recLev1p->records[0],0);
04331 status = drms_segment_getdims(segin,&di);
04332 if(status != DRMS_SUCCESS)
04333 {
04334 printf("Error: cannot read the dimensions of the data segment of level 1p data\n");
04335 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04336 }
04337 axisin[0]= di.axis[0];
04338 axisin[1]= di.axis[1];
04339 axisout[0]=axisin[0];
04340 axisout[1]=axisin[1];
04341
04342 QUALITYLEV1 = drms_getkey_int(recLev1p->records[0],QUALLEV1S,&status);
04343 if(status != DRMS_SUCCESS)
04344 {
04345 printf("Error: cannot read the keyword %s\n",QUALLEV1S);
04346 }
04347 strcpy(source,drms_getkey_string(recLev1p->records[0],SOURCES,&status));
04348 if(status != DRMS_SUCCESS)
04349 {
04350 printf("Error: cannot read the keyword %s\n",SOURCES);
04351 }
04352 QUALITY = drms_getkey_int(recLev1p->records[0],QUALITYS,&status);
04353 if(status != DRMS_SUCCESS)
04354 {
04355 printf("Error: cannot read the keyword %s\n",QUALITYS);
04356 }
04357
04358
04359
04360 Segments1p=0;
04361
04362 }
04363 else
04364 {
04365 printf("Unable to open the series %s, or no record found, for time %s\n",HMISeries,timeBegin2);
04366 if(status == DRMS_SUCCESS) QUALITY = QUALITY | QUAL_MISSINGLEV1P;
04367 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04368 }
04369 }
04370
04371
04372
04373
04374
04375
04376
04377
04378
04379
04380
04381
04382
04383
04384
04385
04386 if (Lev1pWanted == 1 || (Lev15Wanted == 1 && TestLevIn[2]!=1))
04387 {
04388
04389
04390 printf("PRODUCING LEVEL 1P RECORDS\n");
04391
04392 int ActualnSegs1d=0;
04393
04394
04395
04396 if(Segments1d == 0)
04397 {
04398 printf("READING THE LEVEL 1d DATA SEGMENTS\n");
04399
04400
04401 statusA[0]=1;
04402 i=0;
04403 while(statusA[0] != DRMS_SUCCESS)
04404 {
04405 TargetHFLID = drms_getkey_int(recLev1d->records[i],HFLIDS,&statusA[0]);
04406 if(i < nRecs1d-1) ++i;
04407 else break;
04408 }
04409 statusA[1]=1;
04410 i=0;
04411 while(statusA[1] != DRMS_SUCCESS)
04412 {
04413 TargetHPLTID = drms_getkey_int(recLev1d->records[i],HPLTIDS,&statusA[1]);
04414 if(i < nRecs1d-1) ++i;
04415 else break;
04416 }
04417 statusA[2]=1;
04418 i=0;
04419 while(statusA[2] != DRMS_SUCCESS)
04420 {
04421 TargetHWLTID = drms_getkey_int(recLev1d->records[i],HWLTIDS,&statusA[2]);
04422 if(i < nRecs1d-1) ++i;
04423 else break;
04424 }
04425 statusA[3]=1;
04426 i=0;
04427 while(statusA[3] != DRMS_SUCCESS)
04428 {
04429 WavelengthID2 = drms_getkey_int(recLev1d->records[i],WavelengthIDS,&statusA[3]);
04430 if(i < nRecs1d-1) ++i;
04431 else break;
04432 }
04433 if(WavelengthID2 != WavelengthID)
04434 {
04435 printf("Error: WavelengthID2 != WavelengthID\n");
04436 return 1;
04437 }
04438
04439
04440 if( (statusA[0]+statusA[1]+statusA[2]+statusA[3]) != 0)
04441 {
04442 printf("Error: some needed keywords cannot be read on level 1d data at target time %s\n",timeBegin2);
04443 Segments1d = 0;
04444 QUALITY = QUALITY | QUAL_MISSINGKEYWORDLEV1D;
04445 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04446 }
04447 else
04448 {
04449 framelistSize = framelistInfo(TargetHFLID,TargetHPLTID,TargetHWLTID,WavelengthID,PHWPLPOS,WavelengthIndex,WavelengthLocation,&PolarizationType,CamId,&combine,&npol,MaxNumFiltergrams,&CadenceRead,CameraValues,FIDValues);
04450 if(framelistSize == 1) return 1;
04451 }
04452
04453 if(framelistSize != nRecs1d)
04454 {
04455 printf("Error: there are %d level 1d records open at target time %s, but the expected number of records is %d\n",nRecs1d,timeBegin2,framelistSize);
04456 QUALITY = QUALITY | QUAL_MISSINGKEYWORDLEV1D;
04457 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04458 }
04459
04460 ActualnSegs1d=nRecs1d;
04461
04462 arrLev1d = (DRMS_Array_t **)malloc(nRecs1d*sizeof(DRMS_Array_t *));
04463 if(arrLev1d == NULL)
04464 {
04465 printf("Error: memory could not be allocated to arrLev1d\n");
04466 return 1;
04467 }
04468
04469 for(i=0;i<nRecs1d;++i)
04470 {
04471
04472
04473
04474
04475
04476
04477
04478
04479 if(recLev1d->records[i] != NULL)
04480 {
04481 segin = drms_segment_lookupnum(recLev1d->records[i], 0);
04482 arrLev1d[i] = drms_segment_read(segin,type1d,&status);
04483 if(status != DRMS_SUCCESS || arrLev1d[i] == NULL)
04484 {
04485 printf("Error: could not read the segment for level 1d data index %d at target time %s \n",i,timeBegin2);
04486 arrLev1d[i] = NULL;
04487 ActualnSegs1d-=1;
04488 }
04489 }
04490 else
04491 {
04492 arrLev1d[i] = NULL;
04493 ActualnSegs1d-=1;
04494 }
04495 }
04496
04497 Segments1d = 1;
04498 }
04499 else
04500 {
04501 ActualnSegs1d=nRecs1d;
04502 for(i=0;i<nRecs1d;++i)
04503 {
04504 if(arrLev1d[i] == NULL) ActualnSegs1d-=1;
04505 }
04506 }
04507
04508
04509
04510
04511
04512
04513 if(ActualnSegs1d != framelistSize)
04514 {
04515 printf("Error: some records for the level 1d filtergrams are missing to produce level 1p data at target time %s\n",timeBegin2);
04516 Segments1d=0;
04517 QUALITY = QUALITY | QUAL_MISSINGLEV1D;
04518 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04519 }
04520
04521
04522
04523
04524
04525 if(PolarizationType ==3)
04526 {
04527 nRecs1p=ActualnSegs1d/2;
04528 npolout=2;
04529 nSegs1p=npolout*nRecs1p;
04530 Lev1pOffset=40;
04531 if (Lev1pWanted) recLev1p = drms_create_records(drms_env,1,HMISeriesLev1pb,DRMS_PERMANENT,&status);
04532 else recLev1p = drms_create_records(drms_env,1,HMISeriesLev1pb,DRMS_TRANSIENT,&status);
04533 }
04534 if(PolarizationType ==2)
04535 {
04536 nRecs1p=ActualnSegs1d/4;
04537 npolout=2;
04538 nSegs1p=npolout*nRecs1p;
04539 Lev1pOffset=40;
04540 if (Lev1pWanted) recLev1p = drms_create_records(drms_env,1,HMISeriesLev1pb,DRMS_PERMANENT,&status);
04541 else recLev1p = drms_create_records(drms_env,1,HMISeriesLev1pb,DRMS_TRANSIENT,&status);
04542 }
04543 if(PolarizationType ==1)
04544 {
04545 nRecs1p=ActualnSegs1d/npol;
04546 npolout=4;
04547 nSegs1p=npolout*nRecs1p;
04548 Lev1pOffset=0;
04549 if (Lev15Wanted)
04550 {
04551 printf("Warning: you asked for a lev1.5 output but the code can only produce lev1p data from this framelist and the camera selected\n");
04552 Lev15Wanted = 0;
04553 }
04554 if (Lev1pWanted) recLev1p = drms_create_records(drms_env,1,HMISeriesLev1pa,DRMS_PERMANENT,&status);
04555 else recLev1p = drms_create_records(drms_env,1,HMISeriesLev1pa,DRMS_TRANSIENT,&status);
04556 }
04557
04558 if (status != DRMS_SUCCESS || recLev1p == NULL)
04559 {
04560 printf("Could not create a record for the level 1p series at target time %s\n",timeBegin2);
04561
04562
04563
04564 return 1;
04565 }
04566
04567
04568
04569
04570 arrLev1p = (DRMS_Array_t **)malloc(nSegs1p*sizeof(DRMS_Array_t *));
04571 if(arrLev1p == NULL)
04572 {
04573 printf("Error: memory could not be allocated to arrLev1p\n");
04574 return 1;
04575 }
04576
04577 for(i=0;i<nSegs1p;++i)
04578 {
04579 arrLev1p[i] = drms_array_create(type1p,2,axisout,NULL,&status);
04580 if(status != DRMS_SUCCESS || arrLev1p[i] == NULL)
04581 {
04582 printf("Error: cannot create an array for a level 1p data at target time %s\n",timeBegin2);
04583
04584
04585
04586
04587
04588
04589
04590
04591
04592
04593
04594
04595 return 1;
04596 }
04597 }
04598
04599 images = (float **)malloc(npol*sizeof(float *));
04600 if(images == NULL)
04601 {
04602 printf("Error: memory could not be allocated to images\n");
04603 return 1;
04604 }
04605
04606 imagesout = (float **)malloc(npolout*sizeof(float *));
04607 if(imagesout == NULL)
04608 {
04609 printf("Error: memory could not be allocated to imagesout\n");
04610 return 1;
04611 }
04612
04613 ps1 = (int *)malloc(npol*sizeof(int *));
04614 if(ps1 == NULL)
04615 {
04616 printf("Error: memory could not be allocated to ps1\n");
04617 return 1;
04618 }
04619 ps2 = (int *)malloc(npol*sizeof(int *));
04620 if(ps2 == NULL)
04621 {
04622 printf("Error: memory could not be allocated to ps2\n");
04623 return 1;
04624 }
04625 ps3 = (int *)malloc(npol*sizeof(int *));
04626 if(ps3 == NULL)
04627 {
04628 printf("Error: memory could not be allocated to ps3\n");
04629 return 1;
04630 }
04631
04632 fid = (int *)malloc(nRecs1d*sizeof(int *));
04633 if(fid == NULL)
04634 {
04635 printf("Error: memory could not be allocated to fid\n");
04636 return 1;
04637 }
04638
04639 Wavelengths = (int *)malloc(nRecs1p*sizeof(int *));
04640 if(Wavelengths == NULL)
04641 {
04642 printf("Error: memory could not be allocated to Wavelengths\n");
04643 return 1;
04644 }
04645
04646
04647
04648 printf("Looking for available level 1d arrays. Number expected: %d\n",nRecs1d);
04649 ii=0;
04650 for(i=0;i<nRecs1d;++i)
04651 {
04652 if (arrLev1d[i] != NULL)
04653 {
04654 fid[i]=drms_getkey_int(recLev1d->records[i],FIDS,&status);
04655 if( status != 0)
04656 {
04657 printf("Error: unable to read the keyword fid in a level 1d record\n");
04658 QUALITY = QUALITY | QUAL_MISSINGKEYWORDLEV1D;
04659 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04660 }
04661 printf("FID of lev1d data: %d\n",fid[i]);
04662 if(i == 0) Wavelengths[0]=WhichWavelength(fid[0]);
04663 temp=0;
04664 for(k=0;k<=ii;++k) if(WhichWavelength(fid[i]) == Wavelengths[k]) temp=1;
04665 if(temp == 0)
04666 {
04667 ii+=1;
04668 Wavelengths[ii]=WhichWavelength(fid[i]);
04669 }
04670 }
04671 else
04672 {
04673 printf("Error: wavelength index %d is missing\n",i);
04674 QUALITY = QUALITY | QUAL_MISSINGLEV1D;
04675 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04676 }
04677 }
04678
04679 if(ii+1 != nRecs1p)
04680 {
04681 printf("Error: the number of wavelengths in level 1d data: %d; is not what it should be: %d\n",ii+1,nRecs1p);
04682 QUALITY = QUALITY | QUAL_WRONGWAVELENGTHNUM;
04683 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04684 }
04685
04686
04687 if(QuickLook == 1)
04688 {
04689 TSEL=20.;
04690 TFRONT=20.;
04691 }
04692 else
04693 {
04694
04695 strcpy(HMISeriesTemp,HMISeriesTemperature);
04696 strcat(HMISeriesTemp,"[");
04697 strcat(HMISeriesTemp,timeBegin2);
04698 strcat(HMISeriesTemp,"]");
04699 rectemp=NULL;
04700 rectemp = drms_open_records(drms_env,HMISeriesTemp,&statusA[0]);
04701 printf("TEMPERATURE QUERY = %s\n",HMISeriesTemp);
04702 if(statusA[0] == DRMS_SUCCESS && rectemp != NULL && rectemp->n != 0) TSEL=drms_getkey_float(rectemp->records[0],TS08,&status);
04703 else status = 1;
04704 if(status != DRMS_SUCCESS || isnan(TSEL))
04705 {
04706 printf("Error: the temperature keyword %s could not be read\n",TS08);
04707 QUALITY = QUALITY | QUAL_NOTEMP;
04708 TSEL=20.;
04709 }
04710 statusA[1]=1;
04711 if(statusA[0] == DRMS_SUCCESS && rectemp != NULL && rectemp->n != 0) TFRONT=(drms_getkey_float(rectemp->records[0],TS01,&statusA[0])+drms_getkey_float(rectemp->records[0],TS02,&statusA[1]))/2.0;
04712 if(statusA[0] != DRMS_SUCCESS || statusA[1] != DRMS_SUCCESS || isnan(TFRONT))
04713 {
04714 printf("Error: temperature keyword %s and/or %s could not be read\n",TS01,TS02);
04715 QUALITY = QUALITY | QUAL_NOTEMP;
04716 TFRONT=20.;
04717 }
04718 printf("TEMPERATURES = %f %f\n",TSEL,TFRONT);
04719 if(rectemp != NULL)
04720 {
04721 drms_close_records(rectemp,DRMS_FREE_RECORD);
04722 rectemp=NULL;
04723 }
04724 }
04725
04726 printf("WAVELENGTH ORDER = ");
04727 for(i=0;i<nRecs1p;++i) printf(" %d ",Wavelengths[i]);
04728 printf("\n");
04729
04730 printf("CALIBRATION OF POLARIZATION AND SORTING OF THE IMAGES IN THE ORDER I0, I1, I2, I3, I4... \n");
04731
04732 for(k=0;k<nRecs1p;++k)
04733 {
04734
04735 i=0;
04736 for(ii=0;ii<nRecs1d;++ii) if (WhichWavelength(fid[ii]) == k)
04737 {
04738 printf("wavelength=%d, polarization %d\n",k,i);
04739 images[i]=arrLev1d[ii]->data;
04740 ps1[i]=drms_getkey_int(recLev1d->records[ii],HPL1POSS,&statusA[0]);
04741 ps2[i]=drms_getkey_int(recLev1d->records[ii],HPL2POSS,&statusA[1]);
04742 ps3[i]=drms_getkey_int(recLev1d->records[ii],HPL3POSS,&statusA[2]);
04743 if( (statusA[0]+statusA[1]+statusA[2]) != 0)
04744 {
04745 printf("Error: unable to read one or several keyword(s) in level 1d data\n");
04746 QUALITY = QUALITY | QUAL_MISSINGKEYWORDLEV1D;
04747 if(Lev15Wanted) CreateEmptyRecord=1; goto NextTargetTime;
04748 }
04749 i+=1;
04750 }
04751
04752 for(i=0;i<npolout;++i) imagesout[i]=arrLev1p[k*npolout+i]->data;
04753
04754
04755
04756
04757
04758 printf("Producing level 1p data\n");
04759 t0=dsecnd();
04760 polcal(&pars,npol,PolarizationType,images,imagesout,ps1,ps2,ps3,TSEL,TFRONT,axisout[0],axisout[1],axisout[1]);
04761 t1=dsecnd();
04762 printf("TIME ELAPSED IN POLCAL: %f\n",t1-t0);
04763
04764
04765
04766
04767 if(Lev1pWanted)
04768 {
04769 t0=dsecnd();
04770 for(i=0;i<npolout;++i)
04771 {
04772 segout = drms_segment_lookup(recLev1p->records[0],Lev1pSegName[i+Lev1pOffset+k*npolout]);
04773 arrLev1p[k*npolout+i]->bzero=segout->bzero;
04774 arrLev1p[k*npolout+i]->bscale=segout->bscale;
04775 arrLev1p[k*npolout+i]->israw=0;
04776 drms_segment_write(segout,arrLev1p[k*npolout+i], 0);
04777 }
04778 t1=dsecnd();
04779 printf("TIME ELAPSED TO WRITE THE LEVEL 1p SEGMENTS: %f\n",t1-t0);
04780 }
04781
04782 }
04783
04784
04785
04786
04787
04788 drms_copykeys(recLev1p->records[0],recLev1d->records[0],0, kDRMS_KeyClass_Explicit);
04789
04790
04791 statusA[0] = drms_setkey_float(recLev1p->records[0] ,TFRONTS,TFRONT);
04792 statusA[1] = drms_setkey_float(recLev1p->records[0] ,TSELS,TSEL);
04793 statusA[2] = drms_setkey_int(recLev1p->records[0] ,POLCALMS,method);
04794 statusA[3] = drms_setkey_string(recLev1p->records[0],CODEVER3S,CODEVERSION3);
04795 DSUNOBSint = drms_getkey_double(recLev1d->records[0],DSUNOBSS,&status);
04796 if(status != DRMS_SUCCESS || isnan(DSUNOBSint))
04797 {
04798 printf("Error: %s keyword cannot be read on level 1d data at target time %s\n",DSUNOBSS,timeBegin2);
04799 statusA[4]=1;
04800 }
04801 else
04802 {
04803 ctime1=asin((double)solar_radius/DSUNOBSint)*180.*60.*60./M_PI;
04804 printf("RSUN_OBS = %f\n",ctime1);
04805 statusA[4] = drms_setkey_float(recLev1p->records[0],RSUNOBSS,ctime1);
04806 }
04807 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
04808 statusA[5]= drms_setkey_string(recLev1p->records[0],DATES,DATEOBS);
04809 strcat(source,"]");
04810 statusA[6]= drms_setkey_string(recLev1p->records[0],SOURCES,source);
04811 statusA[7]= drms_setkey_int(recLev1p->records[0],QUALLEV1S,QUALITYLEV1);
04812
04813 TotalStatus=0;
04814 for(i=0;i<8;++i) TotalStatus+=statusA[i];
04815 if(TotalStatus != 0)
04816 {
04817 for(i=0;i<8;++i) printf(" %d ",statusA[i]);
04818 printf("\n");
04819 printf("WARNING: could not set some of the keywords for the level 1p data at target time %s\n",timeBegin2);
04820 }
04821
04822
04823
04824
04825 if(Lev1pWanted)
04826 {
04827
04828
04829 X0AVG = (float)drms_getkey_double(recLev1d->records[0],CRPIX1S,&status);
04830 if(status != DRMS_SUCCESS || isnan(X0AVG))
04831 {
04832 printf("Error: %s keyword cannot be read on level 1d data at target time %s\n",CRPIX1S,timeBegin2);
04833 }
04834 X0AVG=X0AVG-1;
04835 printf("X0AVG= %f\n",X0AVG);
04836
04837 Y0AVG = (float)drms_getkey_double(recLev1d->records[0],CRPIX2S,&status);
04838 if(status != DRMS_SUCCESS || isnan(Y0AVG))
04839 {
04840 printf("Error: %s keyword cannot be read on level 1d data at target time %s\n",CRPIX2S,timeBegin2);
04841 }
04842 Y0AVG=Y0AVG-1;
04843 printf("Y0AVG= %f\n",Y0AVG);
04844
04845
04846 for(k=0;k<nRecs1p;++k)
04847 {
04848 for(i=0;i<npolout;++i)
04849 {
04850 strcpy(TOTVALSSS[k*npolout+i],"TOTVALS[");
04851 sprintf(query,"%d",k*npolout+i);
04852 strcat(TOTVALSSS[k*npolout+i],query);
04853 strcat(TOTVALSSS[k*npolout+i],"]");
04854
04855 strcpy(MISSVALSSS[k*npolout+i],"MISSVALS[");
04856 sprintf(query,"%d",k*npolout+i);
04857 strcat(MISSVALSSS[k*npolout+i],query);
04858 strcat(MISSVALSSS[k*npolout+i],"]");
04859
04860 strcpy(DATAVALSSS[k*npolout+i],"DATAVALS[");
04861 sprintf(query,"%d",k*npolout+i);
04862 strcat(DATAVALSSS[k*npolout+i],query);
04863 strcat(DATAVALSSS[k*npolout+i],"]");
04864
04865 strcpy(DATAMEANSS[k*npolout+i],"DATAMEA2[");
04866 sprintf(query,"%d",k*npolout+i);
04867 strcat(DATAMEANSS[k*npolout+i],query);
04868 strcat(DATAMEANSS[k*npolout+i],"]");
04869
04870 strcpy(DATAMINSS[k*npolout+i],"DATAMIN2[");
04871 sprintf(query,"%d",k*npolout+i);
04872 strcat(DATAMINSS[k*npolout+i],query);
04873 strcat(DATAMINSS[k*npolout+i],"]");
04874
04875 strcpy(DATAMAXSS[k*npolout+i],"DATAMAX2[");
04876 sprintf(query,"%d",k*npolout+i);
04877 strcat(DATAMAXSS[k*npolout+i],query);
04878 strcat(DATAMAXSS[k*npolout+i],"]");
04879
04880 strcpy(DATAMEDNSS[k*npolout+i],"DATAMED2[");
04881 sprintf(query,"%d",k*npolout+i);
04882 strcat(DATAMEDNSS[k*npolout+i],query);
04883 strcat(DATAMEDNSS[k*npolout+i],"]");
04884
04885 strcpy(DATARMSSS[k*npolout+i],"DATARMS2[");
04886 sprintf(query,"%d",k*npolout+i);
04887 strcat(DATARMSSS[k*npolout+i],query);
04888 strcat(DATARMSSS[k*npolout+i],"]");
04889
04890 strcpy(DATASKEWSS[k*npolout+i],"DATASKE2[");
04891 sprintf(query,"%d",k*npolout+i);
04892 strcat(DATASKEWSS[k*npolout+i],query);
04893 strcat(DATASKEWSS[k*npolout+i],"]");
04894
04895 strcpy(DATAKURTSS[k*npolout+i],"DATAKUR2[");
04896 sprintf(query,"%d",k*npolout+i);
04897 strcat(DATAKURTSS[k*npolout+i],query);
04898 strcat(DATAKURTSS[k*npolout+i],"]");
04899
04900 strcpy(DATAMINS2S[k*npolout+i],"DATAMIN[");
04901 sprintf(query,"%d",k*npolout+i);
04902 strcat(DATAMINS2S[k*npolout+i],query);
04903 strcat(DATAMINS2S[k*npolout+i],"]");
04904
04905 strcpy(DATAMAXS2S[k*npolout+i],"DATAMAX[");
04906 sprintf(query,"%d",k*npolout+i);
04907 strcat(DATAMAXS2S[k*npolout+i],query);
04908 strcat(DATAMAXS2S[k*npolout+i],"]");
04909
04910 strcpy(DATAMEANS2S[k*npolout+i],"DATAMEAN[");
04911 sprintf(query,"%d",k*npolout+i);
04912 strcat(DATAMEANS2S[k*npolout+i],query);
04913 strcat(DATAMEANS2S[k*npolout+i],"]");
04914
04915 strcpy(DATAMEDNS2S[k*npolout+i],"DATAMEDN[");
04916 sprintf(query,"%d",k*npolout+i);
04917 strcat(DATAMEDNS2S[k*npolout+i],query);
04918 strcat(DATAMEDNS2S[k*npolout+i],"]");
04919
04920 strcpy(DATARMSS2S[k*npolout+i],"DATARMS[");
04921 sprintf(query,"%d",k*npolout+i);
04922 strcat(DATARMSS2S[k*npolout+i],query);
04923 strcat(DATARMSS2S[k*npolout+i],"]");
04924
04925 strcpy(DATASKEWS2S[k*npolout+i],"DATASKEW[");
04926 sprintf(query,"%d",k*npolout+i);
04927 strcat(DATASKEWS2S[k*npolout+i],query);
04928 strcat(DATASKEWS2S[k*npolout+i],"]");
04929
04930 strcpy(DATAKURTS2S[k*npolout+i],"DATAKURT[");
04931 sprintf(query,"%d",k*npolout+i);
04932 strcat(DATAKURTS2S[k*npolout+i],query);
04933 strcat(DATAKURTS2S[k*npolout+i],"]");
04934
04935 image=(float *)arrLev1p[k*npolout+i]->data;
04936
04937
04938 status=fstats(axisout[0]*axisout[1],image,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
04939 if(status != 0)
04940 {
04941 printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
04942 }
04943 statusA[6]= drms_setkey_float(recLev1p->records[0],DATAMINSS[k*npolout+i],(float)minimum);
04944 statusA[7]= drms_setkey_float(recLev1p->records[0],DATAMAXSS[k*npolout+i],(float)maximum);
04945 statusA[8]= drms_setkey_float(recLev1p->records[0],DATAMEDNSS[k*npolout+i],(float)median);
04946 statusA[9]= drms_setkey_float(recLev1p->records[0],DATAMEANSS[k*npolout+i],(float)mean);
04947 statusA[10]= drms_setkey_float(recLev1p->records[0],DATARMSSS[k*npolout+i],(float)sigma);
04948 statusA[11]= drms_setkey_float(recLev1p->records[0],DATASKEWSS[k*npolout+i],(float)skewness);
04949 statusA[12]= drms_setkey_float(recLev1p->records[0],DATAKURTSS[k*npolout+i],(float)kurtosis);
04950 statusA[13]= drms_setkey_int(recLev1p->records[0],TOTVALSSS[k*npolout+i],axisout[0]*axisout[1]);
04951 statusA[14]= drms_setkey_int(recLev1p->records[0],DATAVALSSS[k*npolout+i],ngood);
04952 statusA[15]= drms_setkey_int(recLev1p->records[0],MISSVALSSS[k*npolout+i],axisout[0]*axisout[1]-ngood);
04953
04954 for(ii=0;ii<axisout[0]*axisout[1];++ii)
04955 {
04956 row =ii / axisout[0];
04957 column=ii % axisout[0];
04958 distance = sqrt(((float)row-Y0AVG)*((float)row-Y0AVG)+((float)column-X0AVG)*((float)column-X0AVG));
04959 if(distance > 0.99*RSUNint) image[ii]=NAN;
04960 }
04961
04962 status=fstats(axisout[0]*axisout[1],image,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
04963 if(status != 0)
04964 {
04965 printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
04966 }
04967 statusA[16]= drms_setkey_float(recLev1p->records[0],DATAMINS2S[k*npolout+i],(float)minimum);
04968 statusA[17]= drms_setkey_float(recLev1p->records[0],DATAMAXS2S[k*npolout+i],(float)maximum);
04969 statusA[18]= drms_setkey_float(recLev1p->records[0],DATAMEDNS2S[k*npolout+i],(float)median);
04970 statusA[19]= drms_setkey_float(recLev1p->records[0],DATAMEANS2S[k*npolout+i],(float)mean);
04971 statusA[20]= drms_setkey_float(recLev1p->records[0],DATARMSS2S[k*npolout+i],(float)sigma);
04972 statusA[21]= drms_setkey_float(recLev1p->records[0],DATASKEWS2S[k*npolout+i],(float)skewness);
04973 statusA[22]= drms_setkey_float(recLev1p->records[0],DATAKURTS2S[k*npolout+i],(float)kurtosis);
04974
04975 TotalStatus=0;
04976 for(ii=6;ii<23;++ii) TotalStatus+=statusA[ii];
04977 if(TotalStatus != 0)
04978 {
04979 for(ii=6;ii<23;++ii) printf(" %d ",statusA[ii]);
04980 printf("\n");
04981 printf("WARNING: could not set some of the keywords for the level 1p data at target time %s\n",timeBegin2);
04982 }
04983
04984 }
04985 }
04986 }
04987
04988 Segments1p=1;
04989
04990 }
04991
04992
04993
04994
04995
04996
04997
04998
04999
05000
05001
05002
05003 if (Lev15Wanted)
05004 {
05005
05006
05007
05008
05009
05010 int ActualnSegs1p=nSegs1p;
05011
05012 if(Segments1p == 0)
05013 {
05014 if(PolarizationType ==2 || PolarizationType ==3) nSegs1p=nRecs1p*2;
05015 if(PolarizationType ==1) nSegs1p=nRecs1p*4;
05016
05017 arrLev1p = (DRMS_Array_t **)malloc(nSegs1p*sizeof(DRMS_Array_t *));
05018 if(arrLev1p == NULL)
05019 {
05020 printf("Error: memory could not be allocated to arrLev1p\n");
05021 return 1;
05022 }
05023
05024
05025
05026
05027
05028
05029
05030
05031
05032
05033 printf("READING DATA SEGMENTS OF LEVEL 1p RECORDS\n");
05034 for(i=0;i<nSegs1p;++i)
05035 {
05036 segin = drms_segment_lookupnum(recLev1p->records[0],i);
05037 arrLev1p[i] = drms_segment_read(segin,type1p,&status);
05038 if(status != DRMS_SUCCESS || arrLev1p[i] == NULL)
05039 {
05040 printf("Error: could not read the segment for level 1p data index %d at target time %s \n",i,timeBegin2);
05041 arrLev1p[i] = NULL;
05042 ActualnSegs1p-=1;
05043 }
05044 }
05045
05046
05047
05048 if(PolarizationType == 1)
05049 {
05050 if(ActualnSegs1p != nSegs1p)
05051 {
05052 printf("Error: some level 1p data are missing to produce level 1.5 data: %d\n",ActualnSegs1p);
05053 QUALITY = QUALITY | QUAL_MISSINGLEV1P;
05054 CreateEmptyRecord=1; goto NextTargetTime;
05055 }
05056 else
05057 {
05058 printf("PRODUCING LCP+RCP FROM I,Q,U,V DATA\n");
05059 LCP=(float *)malloc(axisin[0]*axisin[1]*sizeof(float));
05060 RCP=(float *)malloc(axisin[0]*axisin[1]*sizeof(float));
05061 if(LCP == NULL || RCP == NULL)
05062 {
05063 printf("Error: could not create LCP and/or RCP array \n");
05064 return 1;
05065 }
05066 j=0;
05067 for (i=0;i<nRecs1p;++i)
05068 {
05069 temparr1=(float *)arrLev1p[i*4]->data;
05070 temparr2=(float *)arrLev1p[i*4+3]->data;
05071 temparr3=(float *)arrLev1p[j*2]->data;
05072 temparr4=(float *)arrLev1p[j*2+1]->data;
05073 for(ii=0;ii<axisin[0]*axisin[1];++ii)
05074 {
05075 LCP[ii]=temparr1[ii]+temparr2[ii];
05076 RCP[ii]=temparr1[ii]-temparr2[ii];
05077 temparr3[ii]=LCP[ii];
05078 temparr4[ii]=RCP[ii];
05079 }
05080 j+=1;
05081 }
05082 free(LCP);
05083 free(RCP);
05084 LCP=NULL;
05085 RCP=NULL;
05086 PolarizationType=3;
05087 nSegs1p=nRecs1p*2;
05088 ActualnSegs1p = nSegs1p;
05089 for(i=nSegs1p;i<nRecs1p*4;++i) if(arrLev1p[i] != NULL)
05090 {
05091 drms_free_array(arrLev1p[i]);
05092 arrLev1p[i]=NULL;
05093 }
05094 }
05095 }
05096
05097 printf("READING DONE\n");
05098
05099
05100 TargetHFLID = drms_getkey_int(recLev1p->records[0],HFLIDS,&status);
05101 if(status != DRMS_SUCCESS)
05102 {
05103 printf("Error: HFLID keyword cannot be read on level 1p data at target time %s\n",timeBegin2);
05104 CreateEmptyRecord=1; goto NextTargetTime;
05105 }
05106 TargetHPLTID = drms_getkey_int(recLev1p->records[0],HPLTIDS,&status);
05107 if(status != DRMS_SUCCESS)
05108 {
05109 printf("Error: HPLTID keyword cannot be read on level 1p data at target time %s\n",timeBegin2);
05110 CreateEmptyRecord=1; goto NextTargetTime;
05111 }
05112 TargetHWLTID = drms_getkey_int(recLev1p->records[0],HWLTIDS,&status);
05113 if(status != DRMS_SUCCESS)
05114 {
05115 printf("Error: HWLTID keyword cannot be read on level 1p data at target time %s\n",timeBegin2);
05116 CreateEmptyRecord=1; goto NextTargetTime;
05117 }
05118 WavelengthID2 = drms_getkey_int(recLev1p->records[0],WavelengthIDS,&status);
05119 if(status != DRMS_SUCCESS)
05120 {
05121 printf("Error: WavelengthID2 keyword cannot be read on level 1p data at target time %s\n",timeBegin2);
05122 CreateEmptyRecord=1; goto NextTargetTime;
05123 }
05124 if(WavelengthID2 != WavelengthID)
05125 {
05126 printf("Error: WavelengthID2 != WavelengthID\n");
05127 return 1;
05128 }
05129
05130 Segments1p=1;
05131 framelistSize = framelistInfo(TargetHFLID,TargetHPLTID,TargetHWLTID,WavelengthID,PHWPLPOS,WavelengthIndex,WavelengthLocation,&PolarizationType,CamId,&combine,&npol,MaxNumFiltergrams,&CadenceRead,CameraValues,FIDValues);
05132 if(framelistSize == 1) return 1;
05133 }
05134
05135
05136
05137
05138 if(ActualnSegs1p != nSegs1p)
05139 {
05140 printf("Error: some level 1p data are missing to produce level 1.5 data: %d\n",ActualnSegs1p);
05141 QUALITY = QUALITY | QUAL_MISSINGLEV1P;
05142 CreateEmptyRecord=1; goto NextTargetTime;
05143 }
05144
05145
05146
05147
05148 i=0;
05149 while(WavelengthIndex[i] != 2) i++;
05150 int HCMNBT,HCMWBT,HCMPOLT,HCME1T;
05151
05152 if( (framelistSize/npol == 6) || (framelistSize/npol == 8) || (framelistSize/npol == 10) )
05153 {
05154 HCMNBT = PHWPLPOS[7*i+3]+12;
05155 HCMWBT = PHWPLPOS[7*i+1]+6;
05156 HCMPOLT= PHWPLPOS[7*i+2];
05157 HCME1T = PHWPLPOS[7*i+0]-3;
05158 }
05159 if(framelistSize/npol == 5)
05160 {
05161 HCMNBT = PHWPLPOS[7*i+3];
05162 HCMWBT = PHWPLPOS[7*i+1];
05163 HCMPOLT= PHWPLPOS[7*i+2];
05164 HCME1T = PHWPLPOS[7*i+0];
05165 }
05166 printf("KEYWORD VALUES: %d %d %d %d\n",HCMNBT,HCMWBT,HCME1T,HCMPOLT);
05167
05168 char keylist[]="HWL4POS,HWL3POS,HWL2POS,HWL1POS,NWL,HCAMID,FSN_REC";
05169 int unique=0;
05170 int n0,n1;
05171
05172 int NBC,WBC,E1C,POLC,NC,CAMERAUSED,FSNLOOKUP;
05173
05174 arrayL0 = drms_record_getvector(drms_env,HMISeriesLookup, keylist, typeLO, unique, &status);
05175 if(status != DRMS_SUCCESS)
05176 {
05177 printf("Error: cannot read a list of keywords in the look-up table series\n");
05178 QUALITY = QUALITY | QUAL_NOLOOKUPKEYWORD;
05179 CreateEmptyRecord=1; goto NextTargetTime;
05180 }
05181 printf("DIMENSIONS= %d %d\n",arrayL0->axis[0],arrayL0->axis[1]);
05182 keyL=arrayL0->data;
05183 arrayL1 = drms_record_getvector(drms_env,HMISeriesLookup,TRECS,DRMS_TYPE_DOUBLE , unique, &status);
05184 if(status != DRMS_SUCCESS)
05185 {
05186 printf("Error: cannot read a list of keywords in the look-up table series\n");
05187 QUALITY = QUALITY | QUAL_NOLOOKUPKEYWORD;
05188 CreateEmptyRecord=1; goto NextTargetTime;
05189 }
05190
05191 timeL=arrayL1->data;
05192
05193 n1=arrayL0->axis[1];
05194 n0=arrayL0->axis[0];
05195 if(n1 != arrayL1->axis[1])
05196 {
05197 printf("Error: The number of look-up table records identified by T_REC is not the same as the number of records identified by HWL4POS,HWL3POS,HWL2POS,HWL1POS, and NWL\n");
05198 QUALITY = QUALITY | QUAL_NOLOOKUPRECORD;
05199 CreateEmptyRecord=1; goto NextTargetTime;
05200 }
05201 printf("%d RECORDS FOUND FOR THE LOOK-UP TABLES\n",n1);
05202
05203 count = (double *)malloc(n1*sizeof(double));
05204 if(count == NULL)
05205 {
05206 printf("Error: memory could not be allocated to count\n");
05207 return 1;
05208 }
05209
05210 temptime = 1000000000.0;
05211 temp = 0;
05212
05213 for(i=0;i<n1;++i)
05214 {
05215 count[i] = fabs(timeL[i]-TargetTime);
05216 NBC = abs(keyL[ i]-HCMNBT);
05217 POLC= abs(keyL[ n1+i]-HCMPOLT);
05218 WBC = abs(keyL[2*n1+i]-HCMWBT);
05219 E1C = abs(keyL[3*n1+i]-HCME1T);
05220 NC = abs(keyL[4*n1+i]-framelistSize/npol);
05221 CAMERAUSED = abs(keyL[5*n1+i]-CamId);
05222 if(count[i] < temptime && NBC+WBC+POLC+E1C+NC+CAMERAUSED == 0)
05223 {
05224 temptime=count[i];
05225 temp=i;
05226 }
05227 }
05228 if(temptime == 1000000000.0)
05229 {
05230 printf("Error: could not find a look-up table with the correct keywords to produce level 1.5 data\n");
05231 QUALITY = QUALITY | QUAL_NOLOOKUPRECORD;
05232 CreateEmptyRecord=1; goto NextTargetTime;
05233 }
05234
05235 printf("Index of the retrieved look-up table %d %d\n",temp,n1);
05236 printf("Keyword values of the retrieved look-up table: %d %d %d %d %d %d %d\n",keyL[temp]-HCMNBT,keyL[n1+temp]-HCMPOLT,keyL[2*n1+temp]-HCMWBT,keyL[3*n1+temp]-HCME1T,keyL[4*n1+temp]-framelistSize/npol,keyL[5*n1+temp]-CamId,keyL[6*n1+temp]);
05237 FSNLOOKUP = keyL[6*n1+temp];
05238
05239
05240 sprintf(query,"%d",FSNLOOKUP);
05241 strcpy(HMILookup,HMISeriesLookup);
05242 strcat(HMILookup,"[");
05243 strcat(HMILookup,query);
05244 strcat(HMILookup,"][][");
05245 if(CamId == 0) strcat(HMILookup,"0");
05246 if(CamId == 1) strcat(HMILookup,"1");
05247 if(CamId == 2) strcat(HMILookup,"2");
05248 if(CamId == 3) strcat(HMILookup,"3");
05249 strcat(HMILookup,"][");
05250 sprintf(query,"%d",HCME1T);
05251 strcat(HMILookup,query);
05252 strcat(HMILookup,"][");
05253 sprintf(query,"%d",HCMWBT);
05254 strcat(HMILookup,query);
05255 strcat(HMILookup,"][");
05256 sprintf(query,"%d",HCMPOLT);
05257 strcat(HMILookup,query);
05258 strcat(HMILookup,"][");
05259 sprintf(query,"%d",HCMNBT);
05260 strcat(HMILookup,query);
05261 strcat(HMILookup,"][");
05262 NC=framelistSize/npol;
05263 sprintf(query,"%d",NC);
05264 strcat(HMILookup,query);
05265 strcat(HMILookup,"]");
05266
05267 printf("QUERY= %s\n",HMILookup);
05268
05269 drms_free_array(arrayL0);
05270 arrayL0=NULL;
05271 drms_free_array(arrayL1);
05272 arrayL1=NULL;
05273 free(count);
05274 count=NULL;
05275
05276 lookup = drms_open_records(drms_env,HMILookup,&status);
05277 if (status == DRMS_SUCCESS && lookup != NULL)
05278 {
05279 if (lookup->n > 1)
05280 {
05281 printf("Error: more than 1 lookup table record was downloaded.\n");
05282 QUALITY = QUALITY | QUAL_NOLOOKUPRECORD;
05283 CreateEmptyRecord=1; goto NextTargetTime;
05284 }
05285 if (lookup->n <= 0)
05286 {
05287 printf("Error:no record for the look-up tables were downloaded.\n");
05288 QUALITY = QUALITY | QUAL_NOLOOKUPRECORD;
05289 CreateEmptyRecord=1; goto NextTargetTime;
05290 }
05291 }
05292 else
05293 {
05294 printf("Error: can't open the look-up table series.\n");
05295 QUALITY = QUALITY | QUAL_NOLOOKUPRECORD;
05296 CreateEmptyRecord=1; goto NextTargetTime;
05297 }
05298
05299 segin = drms_segment_lookupnum(lookup->records[0], 0);
05300 arrintable= drms_segment_read(segin, segin->info->type, &status);
05301 if (status != DRMS_SUCCESS)
05302 {
05303 printf("Error: unable to read the data segment of the look-up table record\n");
05304 QUALITY = QUALITY | QUAL_NOLOOKUPRECORD;
05305 CreateEmptyRecord=1; goto NextTargetTime;
05306 }
05307 else printf("look-up table record read\n");
05308
05309
05310
05311
05312
05313 int FSNDIFF;
05314 char keylistCoeff[]="COEFF0,COEFF1,COEFF2,COEFF3";
05315
05316 arrayL0 = drms_record_getvector(drms_env,HMISeriesCoeffs, keylistCoeff, DRMS_TYPE_DOUBLE, unique, &status);
05317 if(status != DRMS_SUCCESS)
05318 {
05319 printf("Error: cannot read a list of keywords in the look-up table series\n");
05320 QUALITY = QUALITY | QUAL_NOCOEFFKEYWORD;
05321 CreateEmptyRecord=1; goto NextTargetTime;
05322 }
05323 TIME *keyt=NULL;
05324 keyt=arrayL0->data;
05325
05326 arrayL1 = drms_record_getvector(drms_env,HMISeriesCoeffs,TRECS,DRMS_TYPE_DOUBLE, unique, &status);
05327 if(status != DRMS_SUCCESS)
05328 {
05329 printf("Error: cannot read a list of keywords in the look-up table series\n");
05330 QUALITY = QUALITY | QUAL_NOCOEFFKEYWORD;
05331 CreateEmptyRecord=1; goto NextTargetTime;
05332 }
05333 timeL=arrayL1->data;
05334
05335 arrayL2 = drms_record_getvector(drms_env,HMISeriesCoeffs,CALFSNS,DRMS_TYPE_INT, unique, &status);
05336 if(status != DRMS_SUCCESS)
05337 {
05338 printf("Error: cannot read a list of keywords in the look-up table series\n");
05339 QUALITY = QUALITY | QUAL_NOCOEFFKEYWORD;
05340 CreateEmptyRecord=1; goto NextTargetTime;
05341 }
05342 FSNL=arrayL2->data;
05343
05344 n1=arrayL0->axis[1];
05345 n0=arrayL0->axis[0];
05346 if(n1 != arrayL1->axis[1] || n1 != arrayL2->axis[1])
05347 {
05348 printf("Error: The number of polynomial coefficient records identified by T_REC is not the same as the number of records identified by COEFF0, COEFF1, COEFF2, and COEFF3, or by FSN_REC\n");
05349 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05350 CreateEmptyRecord=1; goto NextTargetTime;
05351 }
05352 printf("%d RECORDS FOUND FOR THE POLYNOMIAL COEFFICIENTS\n",n1);
05353
05354 count = (double *)malloc(n1*sizeof(double));
05355 if(count == NULL)
05356 {
05357 printf("Error: memory could not be allocated to count\n");
05358 return 1;
05359 }
05360
05361
05362 if(QuickLook != 1)
05363 {
05364 temptime = 1000000000.0;
05365 temp = 0;
05366 temp2= 0;
05367
05368 for(i=0;i<n1;++i)
05369 {
05370 count[i] = fabs(timeL[i]-TargetTime);
05371 FSNDIFF = FSNL[i]-FSNLOOKUP;
05372 if(count[i] < temptime && FSNDIFF == 0 && timeL[i]<TargetTime)
05373 {
05374 temptime=count[i];
05375 temp=i;
05376 }
05377 }
05378 if(temptime > 86400.0)
05379 {
05380 printf("Error: could not find polynomial coefficients with the correct keywords and within 24 hours of the target time to produce level 1.5 data\n");
05381 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05382 CreateEmptyRecord=1; goto NextTargetTime;
05383 }
05384
05385 temptime = 1000000000.0;
05386
05387 for(i=0;i<n1;++i)
05388 {
05389 count[i] = fabs(timeL[i]-TargetTime);
05390 FSNDIFF = FSNL[i]-FSNLOOKUP;
05391 if(count[i] < temptime && FSNDIFF == 0 && timeL[i]>=TargetTime)
05392 {
05393 temptime=count[i];
05394 temp2=i;
05395 }
05396 }
05397 if(temptime > 86400.0)
05398 {
05399 printf("Error: could not find polynomial coefficients with the correct keywords and within 12 hours of the target time to produce level 1.5 data\n");
05400 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05401 CreateEmptyRecord=1; goto NextTargetTime;
05402 }
05403
05404 }
05405 else
05406 {
05407
05408 temptime = 1000000000.0;
05409 temp = 0;
05410
05411 for(i=0;i<n1;++i)
05412 {
05413 count[i] = fabs(timeL[i]-TargetTime);
05414 if(count[i] < temptime)
05415 {
05416 temptime=count[i];
05417 temp=i;
05418 }
05419 }
05420 temp2=temp;
05421 }
05422
05423 printf("Indeces of the retrieved polynomial record %d %d %d\n",temp,temp2,n1);
05424 printf("Keyword values of the retrieved polynomial record: %d\n",FSNL[temp]-FSNLOOKUP);
05425
05426
05427 sprint_time(query,timeL[temp],"TAI",1);
05428 strcpy(HMICoeffs,HMISeriesCoeffs);
05429 strcat(HMICoeffs,"[");
05430 strcat(HMICoeffs,query);
05431 strcat(HMICoeffs,"]");
05432
05433 printf("QUERY= %s\n",HMICoeffs);
05434
05435 recpoly = drms_open_records(drms_env,HMICoeffs,&status);
05436 if (status == DRMS_SUCCESS && recpoly != NULL && recpoly->n != 0)
05437 {
05438 coeff[0]=drms_getkey_double(recpoly->records[0],COEFF0S,&status);
05439 if(status != DRMS_SUCCESS)
05440 {
05441 printf("Error: could not read a polynomial coefficient\n");
05442 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05443 CreateEmptyRecord=1; goto NextTargetTime;
05444 }
05445 coeff[1]=drms_getkey_double(recpoly->records[0],COEFF1S,&status);
05446 if(status != DRMS_SUCCESS)
05447 {
05448 printf("Error: could not read a polynomial coefficient\n");
05449 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05450 CreateEmptyRecord=1; goto NextTargetTime;
05451 }
05452 coeff[2]=drms_getkey_double(recpoly->records[0],COEFF2S,&status);
05453 if(status != DRMS_SUCCESS)
05454 {
05455 printf("Error: could not read a polynomial coefficient\n");
05456 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05457 CreateEmptyRecord=1; goto NextTargetTime;
05458 }
05459 coeff[3]=drms_getkey_double(recpoly->records[0],COEFF3S,&status);
05460 if(status != DRMS_SUCCESS)
05461 {
05462 printf("Error: could not read a polynomial coefficient\n");
05463 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05464 CreateEmptyRecord=1; goto NextTargetTime;
05465 }
05466 }
05467 else
05468 {
05469 printf("Error: can't open the polynomial coefficients series.\n");
05470 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05471 CreateEmptyRecord=1; goto NextTargetTime;
05472 }
05473
05474
05475 if(QuickLook != 1)
05476 {
05477
05478 sprint_time(query,timeL[temp2],"TAI",1);
05479 strcpy(HMICoeffs,HMISeriesCoeffs);
05480 strcat(HMICoeffs,"[");
05481 strcat(HMICoeffs,query);
05482 strcat(HMICoeffs,"]");
05483
05484 printf("QUERY= %s\n",HMICoeffs);
05485
05486 recpoly2 = drms_open_records(drms_env,HMICoeffs,&status);
05487 if (status == DRMS_SUCCESS && recpoly2 != NULL && recpoly2->n != 0)
05488 {
05489 coeff2[0]=drms_getkey_double(recpoly2->records[0],COEFF0S,&status);
05490 if(status != DRMS_SUCCESS)
05491 {
05492 printf("Error: could not read a polynomial coefficient\n");
05493 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05494 CreateEmptyRecord=1; goto NextTargetTime;
05495 }
05496 coeff2[1]=drms_getkey_double(recpoly2->records[0],COEFF1S,&status);
05497 if(status != DRMS_SUCCESS)
05498 {
05499 printf("Error: could not read a polynomial coefficient\n");
05500 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05501 CreateEmptyRecord=1; goto NextTargetTime;
05502 }
05503 coeff2[2]=drms_getkey_double(recpoly2->records[0],COEFF2S,&status);
05504 if(status != DRMS_SUCCESS)
05505 {
05506 printf("Error: could not read a polynomial coefficient\n");
05507 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05508 CreateEmptyRecord=1; goto NextTargetTime;
05509 }
05510 coeff2[3]=drms_getkey_double(recpoly2->records[0],COEFF3S,&status);
05511 if(status != DRMS_SUCCESS)
05512 {
05513 printf("Error: could not read a polynomial coefficient\n");
05514 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05515 CreateEmptyRecord=1; goto NextTargetTime;
05516 }
05517 }
05518 else
05519 {
05520 printf("Error: can't open the polynomial coefficients series.\n");
05521 QUALITY = QUALITY | QUAL_NOCOEFFPRECORD;
05522 CreateEmptyRecord=1; goto NextTargetTime;
05523 }
05524 }
05525 else recpoly2 = NULL;
05526
05527
05528 printf("Polynomial coefficient values: %e %e %e %e\n",coeff[0],coeff[1],coeff[2],coeff[3]);
05529
05530 if(QuickLook != 1)
05531 {
05532 printf("Polynomial coefficient values: %e %e %e %e\n",coeff2[0],coeff2[1],coeff2[2],coeff2[3]);
05533 coeff[0]=(TargetTime-timeL[temp])*(coeff2[0]-coeff[0])/(timeL[temp2]-timeL[temp])+coeff[0];
05534 coeff[1]=(TargetTime-timeL[temp])*(coeff2[1]-coeff[1])/(timeL[temp2]-timeL[temp])+coeff[1];
05535 coeff[2]=(TargetTime-timeL[temp])*(coeff2[2]-coeff[2])/(timeL[temp2]-timeL[temp])+coeff[2];
05536 coeff[3]=(TargetTime-timeL[temp])*(coeff2[3]-coeff[3])/(timeL[temp2]-timeL[temp])+coeff[3];
05537 }
05538
05539
05540 drms_free_array(arrayL0);
05541 arrayL0=NULL;
05542 drms_free_array(arrayL1);
05543 arrayL1=NULL;
05544 drms_free_array(arrayL2);
05545 arrayL2=NULL;
05546 free(count);
05547 count=NULL;
05548
05549 strcpy(HISTORY,"Polynomial Coefficients used for Doppler velocity correction: ");
05550 sprintf(query,"%e",coeff[0]);
05551 strcat(HISTORY,query);
05552 strcat(HISTORY," ");
05553 sprintf(query,"%e",coeff[1]);
05554 strcat(HISTORY,query);
05555 strcat(HISTORY," ");
05556 sprintf(query,"%e",coeff[2]);
05557 strcat(HISTORY,query);
05558 strcat(HISTORY," ");
05559 sprintf(query,"%e",coeff[3]);
05560 strcat(HISTORY,query);
05561 strcat(HISTORY," ");
05562
05563
05564
05565
05566
05567 nRecs15 = 6;
05568
05569 recLev15a = drms_create_records(drms_env,1,HMISeriesLev15a,DRMS_PERMANENT,&statusA[0]);
05570 recLev15b = drms_create_records(drms_env,1,HMISeriesLev15b,DRMS_PERMANENT,&statusA[1]);
05571 recLev15c = drms_create_records(drms_env,1,HMISeriesLev15c,DRMS_PERMANENT,&statusA[2]);
05572 recLev15d = drms_create_records(drms_env,1,HMISeriesLev15d,DRMS_PERMANENT,&statusA[3]);
05573 recLev15e = drms_create_records(drms_env,1,HMISeriesLev15e,DRMS_PERMANENT,&statusA[4]);
05574 printf("Observables will be saved in the following series:\n %s \n %s \n %s \n %s \n %s \n",HMISeriesLev15a,HMISeriesLev15b,HMISeriesLev15c,HMISeriesLev15d,HMISeriesLev15e);
05575
05576
05577 if ( (statusA[0]+statusA[1]+statusA[2]+statusA[3]+statusA[4]) != DRMS_SUCCESS || recLev15a == NULL || recLev15b == NULL || recLev15c == NULL|| recLev15d == NULL || recLev15e == NULL)
05578 {
05579 printf("Could not create a record for one or several level 1.5 data series, at target time %s\n",timeBegin2);
05580
05581
05582
05583
05584
05585
05586
05587 return 1;
05588 }
05589 if(recLev15a->n == 0 || recLev15b->n == 0 || recLev15c->n == 0 || recLev15d->n == 0 || recLev15e->n == 0)
05590 {
05591 printf("Could not create a record for one or several level 1.5 data series, at target time %s\n",timeBegin2);
05592
05593
05594
05595
05596
05597
05598
05599 return 1;
05600 }
05601
05602
05603 arrLev15 = (DRMS_Array_t **)malloc(nRecs15*sizeof(DRMS_Array_t *));
05604 if(arrLev15 == NULL)
05605 {
05606 printf("Error: memory could not be allocated to arrLev15\n");
05607 return 1;
05608 }
05609 for (i=0;i<nRecs15;++i)
05610 {
05611 arrLev15[i] = drms_array_create(type15,2,axisout,NULL,&status);
05612 if(status != DRMS_SUCCESS || arrLev15[i] == NULL)
05613 {
05614 printf("Error: cannot create an array for a level 1.5 data at target time %s\n",timeBegin2);
05615
05616
05617 return 1;
05618 }
05619
05620 }
05621
05622
05623
05624 X0AVG = (float)drms_getkey_double(recLev1p->records[0],CRPIX1S,&status);
05625 if(status != DRMS_SUCCESS || isnan(X0AVG))
05626 {
05627 printf("Error: %s keyword cannot be read on level 1p data at target time %s\n",CRPIX1S,timeBegin2);
05628 QUALITY = QUALITY | QUAL_MISSINGKEYWORDLEV1P;
05629 CreateEmptyRecord=1; goto NextTargetTime;
05630 }
05631 X0AVG=X0AVG-1;
05632
05633 Y0AVG = (float)drms_getkey_double(recLev1p->records[0],CRPIX2S,&status);
05634 if(status != DRMS_SUCCESS || isnan(Y0AVG))
05635 {
05636 printf("Error: %s keyword cannot be read on level 1p data at target time %s\n",CRPIX2S,timeBegin2);
05637 QUALITY = QUALITY | QUAL_MISSINGKEYWORDLEV1P;
05638 CreateEmptyRecord=1; goto NextTargetTime;
05639 }
05640 Y0AVG=Y0AVG-1;
05641
05642 DSUNOBSint = drms_getkey_double(recLev1p->records[0],DSUNOBSS,&status);
05643 if(status != DRMS_SUCCESS || isnan(DSUNOBSint))
05644 {
05645 printf("Error: %s keyword cannot be read on level 1p data at target time %s\n",DSUNOBSS,timeBegin2);
05646 QUALITY = QUALITY | QUAL_MISSINGKEYWORDLEV1P;
05647 CreateEmptyRecord=1; goto NextTargetTime;
05648 }
05649
05650 cdelt1 = (float)drms_getkey_double(recLev1p->records[0],CDELT1S,&status);
05651 if(status != DRMS_SUCCESS || isnan(cdelt1))
05652 {
05653 printf("Error: %s keyword cannot be read on level 1p data at target time %s\n",CDELT1S,timeBegin2);
05654 QUALITY = QUALITY | QUAL_MISSINGKEYWORDLEV1P;
05655 CreateEmptyRecord=1; goto NextTargetTime;
05656 }
05657
05658
05659
05660
05661
05662
05663
05664
05665
05666
05667 RSUNint=1.0/cdelt1*asin((double)solar_radius/DSUNOBSint)*180.*60.*60./M_PI;
05668
05669
05670 DopplerParameters.FSRNB=FSR[0];
05671 DopplerParameters.FSRWB=FSR[1];
05672 DopplerParameters.FSRE1=FSR[2];
05673 DopplerParameters.FSRE2=FSR[3];
05674 DopplerParameters.FSRE3=FSR[4];
05675 DopplerParameters.FSRE4=FSR[5];
05676 DopplerParameters.FSRE5=FSR[6];
05677 DopplerParameters.dlamdv=dlamdv;
05678 DopplerParameters.maxVtest=ntest*2;
05679 DopplerParameters.maxNx=maxNx;
05680 DopplerParameters.ntest=ntest;
05681 DopplerParameters.dvtest=dvtest;
05682 DopplerParameters.MISSINGDATA=MISSINGDATA;
05683 DopplerParameters.MISSINGRESULT=MISSINGRESULT;
05684 DopplerParameters.coeff0=coeff[0];
05685 DopplerParameters.coeff1=coeff[1];
05686 DopplerParameters.coeff2=coeff[2];
05687 DopplerParameters.coeff3=coeff[3];
05688 DopplerParameters.QuickLook=QuickLook;
05689
05690 t0=dsecnd();
05691
05692 Dopplergram(arrLev1p,arrLev15,nSegs1p,arrintable,RSUNint,X0AVG,Y0AVG,DopplerParameters,MISSVALS,&SATVALS,cdelt1,TargetTime);
05693
05694 t1=dsecnd();
05695 printf("TIME ELAPSED IN DOPPLERGRAM(): %f\n",t1-t0);
05696
05697 printf("KEYWORDS OF Dopplergram() %f %f %f\n",RSUNint,X0AVG,Y0AVG);
05698 printf("%f %f %f %f %f %f %f %f %d %d %d %f %f %f \n", DopplerParameters.FSRNB,DopplerParameters.FSRWB,DopplerParameters.FSRE1,DopplerParameters.FSRE2,DopplerParameters.FSRE3,DopplerParameters.FSRE4,DopplerParameters.FSRE5,DopplerParameters.dlamdv,DopplerParameters.maxVtest,DopplerParameters.maxNx,DopplerParameters.ntest,DopplerParameters.dvtest,DopplerParameters.MISSINGDATA,DopplerParameters.MISSINGRESULT);
05699
05700
05701 t0=dsecnd();
05702 segout = drms_segment_lookupnum(recLev15a->records[0], 0);
05703 arrLev15[0]->bzero=segout->bzero;
05704 arrLev15[0]->bscale=segout->bscale;
05705 arrLev15[0]->israw=0;
05706 drms_segment_write(segout,arrLev15[0], 0);
05707
05708 segout = drms_segment_lookupnum(recLev15b->records[0], 0);
05709 arrLev15[1]->bzero=segout->bzero;
05710 arrLev15[1]->bscale=segout->bscale;
05711 arrLev15[1]->israw=0;
05712 drms_segment_write(segout,arrLev15[1], 0);
05713
05714 segout = drms_segment_lookupnum(recLev15c->records[0], 0);
05715 arrLev15[2]->bzero=segout->bzero;
05716 arrLev15[2]->bscale=segout->bscale;
05717 arrLev15[2]->israw=0;
05718 drms_segment_write(segout,arrLev15[2], 0);
05719
05720 segout = drms_segment_lookupnum(recLev15d->records[0], 0);
05721 arrLev15[3]->bzero=segout->bzero;
05722 arrLev15[3]->bscale=segout->bscale;
05723 arrLev15[3]->israw=0;
05724 drms_segment_write(segout,arrLev15[3], 0);
05725
05726 segout = drms_segment_lookupnum(recLev15e->records[0], 0);
05727 arrLev15[4]->bzero=segout->bzero;
05728 arrLev15[4]->bscale=segout->bscale;
05729 arrLev15[4]->israw=0;
05730 drms_segment_write(segout,arrLev15[4], 0);
05731 t1=dsecnd();
05732 printf("TIME ELAPSED TO WRITE THE LEVEL 1.5 SEGMENTS: %f\n",t1-t0);
05733
05734
05735
05736 image0=(float *)arrLev15[5]->data;
05737
05738 for(i=0;i<axisout[0]*axisout[1];++i)
05739 {
05740 row =i / axisout[0];
05741 column=i % axisout[0];
05742 distance = sqrt(((float)row-Y0AVG)*((float)row-Y0AVG)+((float)column-X0AVG)*((float)column-X0AVG));
05743 if(distance > 0.99*RSUNint)
05744 {
05745 image0[i]=NAN;
05746 }
05747 }
05748 status=fstats(axisout[0]*axisout[1],arrLev15[5]->data,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
05749 if(status != 0) printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
05750 statusA[2]= drms_setkey_float(recLev15a->records[0],RAWMEDNS,(float)median);
05751 if(statusA[2] != 0)
05752 {
05753 printf("WARNING: could not set some of the keywords modified by the temporal interpolation subroutine for the Dopplergram at target time %s %f\n",timeBegin2,(float)median);
05754 }
05755
05756
05757
05758 t0=dsecnd();
05759 drms_copykeys(recLev15a->records[0],recLev1p->records[0],1,kDRMS_KeyClass_Explicit);
05760
05761 status=fstats(axisout[0]*axisout[1],arrLev15[0]->data,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
05762 if(status != 0)
05763 {
05764 printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
05765 }
05766
05767 image=arrLev15[0]->data;
05768 for(i=0;i<axisout[0]*axisout[1];i++) if (image[i] > (32767.*arrLev15[0]->bscale+arrLev15[0]->bzero) || image[i] < (-32768.*arrLev15[0]->bscale+arrLev15[0]->bzero))
05769 {
05770 MISSVALS[0]+=1;
05771 ngood -= 1;
05772 }
05773
05774
05775
05776 statusA[0]= drms_setkey_int(recLev15a->records[0],TOTVALSS,ngood+MISSVALS[0]);
05777 statusA[1]= drms_setkey_int(recLev15a->records[0],DATAVALSS,ngood);
05778
05779 statusA[2]= drms_setkey_int(recLev15a->records[0],MISSVALSS,MISSVALS[0]);
05780 statusA[3]= drms_setkey_float(recLev15a->records[0],DATAMINS,(float)minimum);
05781 statusA[4]= drms_setkey_float(recLev15a->records[0],DATAMAXS,(float)maximum);
05782 statusA[5]= drms_setkey_float(recLev15a->records[0],DATAMEDNS,(float)median);
05783 statusA[6]= drms_setkey_float(recLev15a->records[0],DATAMEANS,(float)mean);
05784 statusA[7]= drms_setkey_float(recLev15a->records[0],DATARMSS,(float)sigma);
05785 statusA[8]= drms_setkey_float(recLev15a->records[0],DATASKEWS,(float)skewness);
05786 statusA[9]= drms_setkey_float(recLev15a->records[0],DATAKURTS,(float)kurtosis);
05787 statusA[10]=drms_setkey_int(recLev15a->records[0],CALFSNS,FSNLOOKUP);
05788 statusA[11]=drms_setkey_string(recLev15a->records[0],LUTQUERYS,HMILookup);
05789 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
05790 statusA[12]= drms_setkey_string(recLev15a->records[0],DATES,DATEOBS);
05791 statusA[13]= drms_setkey_int(recLev15a->records[0],QUALITYS,QUALITY);
05792 statusA[14]= drms_setkey_int(recLev15a->records[0],SATVALSS,SATVALS);
05793 statusA[15]= drms_setkey_string(recLev15a->records[0],SOURCES,source);
05794 statusA[16]= drms_setkey_int(recLev15a->records[0],QUALLEV1S,QUALITYLEV1);
05795
05796 TotalStatus=0;
05797 for(i=0;i<17;++i) TotalStatus+=statusA[i];
05798 if(TotalStatus != 0)
05799 {
05800 printf("WARNING: could not set some of the keywords modified by the temporal interpolation subroutine for the Dopplergram at target time %s\n",timeBegin2);
05801 }
05802
05803
05804
05805 drms_copykeys(recLev15b->records[0],recLev1p->records[0],1,kDRMS_KeyClass_Explicit);
05806
05807 status=fstats(axisout[0]*axisout[1],arrLev15[1]->data,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
05808 if(status != 0)
05809 {
05810 printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
05811 }
05812
05813 image=arrLev15[1]->data;
05814 for(i=0;i<axisout[0]*axisout[1];i++) if (image[i] > (2147483647.*arrLev15[1]->bscale+arrLev15[1]->bzero) || image[i] < (-2147483648.*arrLev15[1]->bscale+arrLev15[1]->bzero))
05815 {
05816 MISSVALS[1]+=1;
05817 ngood -= 1;
05818 }
05819
05820
05821 statusA[0]= drms_setkey_int(recLev15b->records[0],TOTVALSS,ngood+MISSVALS[1]);
05822 statusA[1]= drms_setkey_int(recLev15b->records[0],DATAVALSS,ngood);
05823 statusA[2]= drms_setkey_int(recLev15b->records[0],MISSVALSS,MISSVALS[1]);
05824 statusA[3]= drms_setkey_float(recLev15b->records[0],DATAMINS,(float)minimum);
05825 statusA[4]= drms_setkey_float(recLev15b->records[0],DATAMAXS,(float)maximum);
05826 statusA[5]= drms_setkey_float(recLev15b->records[0],DATAMEDNS,(float)median);
05827 statusA[6]= drms_setkey_float(recLev15b->records[0],DATAMEANS,(float)mean);
05828 statusA[7]= drms_setkey_float(recLev15b->records[0],DATARMSS,(float)sigma);
05829 statusA[8]= drms_setkey_float(recLev15b->records[0],DATASKEWS,(float)skewness);
05830 statusA[9]= drms_setkey_float(recLev15b->records[0],DATAKURTS,(float)kurtosis);
05831 statusA[10]=drms_setkey_int(recLev15b->records[0],CALFSNS,FSNLOOKUP);
05832 statusA[11]=drms_setkey_string(recLev15b->records[0],LUTQUERYS,HMILookup);
05833 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
05834 statusA[12]= drms_setkey_string(recLev15b->records[0],DATES,DATEOBS);
05835 statusA[13]= drms_setkey_int(recLev15b->records[0],QUALITYS,QUALITY);
05836 statusA[14]= drms_setkey_int(recLev15b->records[0],SATVALSS,SATVALS);
05837 statusA[15]= drms_setkey_string(recLev15b->records[0],SOURCES,source);
05838 statusA[16]= drms_setkey_int(recLev15b->records[0],QUALLEV1S,QUALITYLEV1);
05839
05840 TotalStatus=0;
05841 for(i=0;i<17;++i) TotalStatus+=statusA[i];
05842 if(TotalStatus != 0)
05843 {
05844 printf("WARNING: could not set some of the keywords modified by the temporal interpolation subroutine for the magnetogram at target time %s\n",timeBegin2);
05845 }
05846
05847
05848
05849 drms_copykeys(recLev15c->records[0],recLev1p->records[0],1,kDRMS_KeyClass_Explicit);
05850
05851 status=fstats(axisout[0]*axisout[1],arrLev15[2]->data,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
05852 if(status != 0)
05853 {
05854 printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
05855 }
05856
05857 image=arrLev15[2]->data;
05858 for(i=0;i<axisout[0]*axisout[1];i++) if (image[i] > (32767.*arrLev15[2]->bscale+arrLev15[2]->bzero) || image[i] < (-32768.*arrLev15[2]->bscale+arrLev15[2]->bzero))
05859 {
05860 MISSVALS[2]+=1;
05861 ngood -= 1;
05862 }
05863
05864
05865 statusA[0]= drms_setkey_int(recLev15c->records[0],TOTVALSS,ngood+MISSVALS[2]);
05866 statusA[1]= drms_setkey_int(recLev15c->records[0],DATAVALSS,ngood);
05867 statusA[2]= drms_setkey_int(recLev15c->records[0],MISSVALSS,MISSVALS[2]);
05868 statusA[3]= drms_setkey_float(recLev15c->records[0],DATAMINS,(float)minimum);
05869 statusA[4]= drms_setkey_float(recLev15c->records[0],DATAMAXS,(float)maximum);
05870 statusA[5]= drms_setkey_float(recLev15c->records[0],DATAMEDNS,(float)median);
05871 statusA[6]= drms_setkey_float(recLev15c->records[0],DATAMEANS,(float)mean);
05872 statusA[7]= drms_setkey_float(recLev15c->records[0],DATARMSS,(float)sigma);
05873 statusA[8]= drms_setkey_float(recLev15c->records[0],DATASKEWS,(float)skewness);
05874 statusA[9]= drms_setkey_float(recLev15c->records[0],DATAKURTS,(float)kurtosis);
05875 statusA[10]=drms_setkey_int(recLev15c->records[0],CALFSNS,FSNLOOKUP);
05876 statusA[11]=drms_setkey_string(recLev15c->records[0],LUTQUERYS,HMILookup);
05877 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
05878 statusA[12]= drms_setkey_string(recLev15c->records[0],DATES,DATEOBS);
05879 statusA[13]= drms_setkey_int(recLev15c->records[0],QUALITYS,QUALITY);
05880 statusA[14]= drms_setkey_int(recLev15c->records[0],SATVALSS,SATVALS);
05881 statusA[15]= drms_setkey_string(recLev15c->records[0],SOURCES,source);
05882 statusA[16]= drms_setkey_int(recLev15c->records[0],QUALLEV1S,QUALITYLEV1);
05883
05884 TotalStatus=0;
05885 for(i=0;i<17;++i) TotalStatus+=statusA[i];
05886 if(TotalStatus != 0)
05887 {
05888 printf("WARNING: could not set some of the keywords modified by the temporal interpolation subroutine for the linedepth at target time %s\n",timeBegin2);
05889 }
05890
05891
05892
05893
05894 drms_copykeys(recLev15d->records[0],recLev1p->records[0],1,kDRMS_KeyClass_Explicit);
05895
05896 status=fstats(axisout[0]*axisout[1],arrLev15[3]->data,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
05897 if(status != 0)
05898 {
05899 printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
05900 }
05901
05902
05903 image=arrLev15[3]->data;
05904 for(i=0;i<axisout[0]*axisout[1];i++) if (image[i] > (32767.*arrLev15[3]->bscale+arrLev15[3]->bzero) || image[i] < (-32768.*arrLev15[3]->bscale+arrLev15[3]->bzero))
05905 {
05906 MISSVALS[3]+=1;
05907 ngood -= 1;
05908 }
05909
05910
05911 statusA[0]= drms_setkey_int(recLev15d->records[0],TOTVALSS,ngood+MISSVALS[3]);
05912 statusA[1]= drms_setkey_int(recLev15d->records[0],DATAVALSS,ngood);
05913 statusA[2]= drms_setkey_int(recLev15d->records[0],MISSVALSS,MISSVALS[3]);
05914 statusA[3]= drms_setkey_float(recLev15d->records[0],DATAMINS,(float)minimum);
05915 statusA[4]= drms_setkey_float(recLev15d->records[0],DATAMAXS,(float)maximum);
05916 statusA[5]= drms_setkey_float(recLev15d->records[0],DATAMEDNS,(float)median);
05917 statusA[6]= drms_setkey_float(recLev15d->records[0],DATAMEANS,(float)mean);
05918 statusA[7]= drms_setkey_float(recLev15d->records[0],DATARMSS,(float)sigma);
05919 statusA[8]= drms_setkey_float(recLev15d->records[0],DATASKEWS,(float)skewness);
05920 statusA[9]= drms_setkey_float(recLev15d->records[0],DATAKURTS,(float)kurtosis);
05921 statusA[10]=drms_setkey_int(recLev15d->records[0],CALFSNS,FSNLOOKUP);
05922 statusA[11]=drms_setkey_string(recLev15d->records[0],LUTQUERYS,HMILookup);
05923 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
05924 statusA[12]= drms_setkey_string(recLev15d->records[0],DATES,DATEOBS);
05925 statusA[13]= drms_setkey_int(recLev15d->records[0],QUALITYS,QUALITY);
05926 statusA[14]= drms_setkey_int(recLev15d->records[0],SATVALSS,SATVALS);
05927 statusA[15]= drms_setkey_string(recLev15d->records[0],SOURCES,source);
05928 statusA[16]= drms_setkey_int(recLev15d->records[0],QUALLEV1S,QUALITYLEV1);
05929
05930 TotalStatus=0;
05931 for(i=0;i<17;++i) TotalStatus+=statusA[i];
05932 if(TotalStatus != 0)
05933 {
05934 printf("WARNING: could not set some of the keywords modified by the temporal interpolation subroutine for the linewidth at target time %s\n",timeBegin2);
05935 }
05936
05937
05938
05939 drms_copykeys(recLev15e->records[0],recLev1p->records[0],1,kDRMS_KeyClass_Explicit);
05940
05941 status=fstats(axisout[0]*axisout[1],arrLev15[4]->data,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
05942 if(status != 0)
05943 {
05944 printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
05945 }
05946
05947 image=arrLev15[4]->data;
05948 for(i=0;i<axisout[0]*axisout[1];i++) if (image[i] > (32767.*arrLev15[4]->bscale+arrLev15[4]->bzero) || image[i] < (-32768.*arrLev15[4]->bscale+arrLev15[4]->bzero))
05949 {
05950 MISSVALS[4]+=1;
05951 ngood -= 1;
05952 }
05953
05954
05955 statusA[0]= drms_setkey_int(recLev15e->records[0],TOTVALSS,ngood+MISSVALS[4]);
05956 statusA[1]= drms_setkey_int(recLev15e->records[0],DATAVALSS,ngood);
05957 statusA[2]= drms_setkey_int(recLev15e->records[0],MISSVALSS,MISSVALS[4]);
05958 statusA[3]= drms_setkey_float(recLev15e->records[0],DATAMINS,(float)minimum);
05959 statusA[4]= drms_setkey_float(recLev15e->records[0],DATAMAXS,(float)maximum);
05960 statusA[5]= drms_setkey_float(recLev15e->records[0],DATAMEDNS,(float)median);
05961 statusA[6]= drms_setkey_float(recLev15e->records[0],DATAMEANS,(float)mean);
05962 statusA[7]= drms_setkey_float(recLev15e->records[0],DATARMSS,(float)sigma);
05963 statusA[8]= drms_setkey_float(recLev15e->records[0],DATASKEWS,(float)skewness);
05964 statusA[9]= drms_setkey_float(recLev15e->records[0],DATAKURTS,(float)kurtosis);
05965 statusA[10]=drms_setkey_int(recLev15e->records[0],CALFSNS,FSNLOOKUP);
05966 statusA[11]=drms_setkey_string(recLev15e->records[0],LUTQUERYS,HMILookup);
05967 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
05968 statusA[12]= drms_setkey_string(recLev15e->records[0],DATES,DATEOBS);
05969 statusA[13]= drms_setkey_int(recLev15e->records[0],QUALITYS,QUALITY);
05970 statusA[14]= drms_setkey_int(recLev15e->records[0],SATVALSS,SATVALS);
05971 statusA[15]= drms_setkey_string(recLev15e->records[0],SOURCES,source);
05972 statusA[16]= drms_setkey_int(recLev15e->records[0],QUALLEV1S,QUALITYLEV1);
05973
05974 TotalStatus=0;
05975 for(i=0;i<17;++i) TotalStatus+=statusA[i];
05976 if(TotalStatus != 0)
05977 {
05978 printf("WARNING: could not set some of the keywords modified by the temporal interpolation subroutine for the continuum intensity at target time %s\n",timeBegin2);
05979 }
05980
05981
05982 image0=(float *)arrLev15[0]->data;
05983 image1=(float *)arrLev15[1]->data;
05984 image2=(float *)arrLev15[2]->data;
05985 image3=(float *)arrLev15[3]->data;
05986 image4=(float *)arrLev15[4]->data;
05987
05988 for(i=0;i<axisout[0]*axisout[1];++i)
05989 {
05990 row =i / axisout[0];
05991 column=i % axisout[0];
05992 distance = sqrt(((float)row-Y0AVG)*((float)row-Y0AVG)+((float)column-X0AVG)*((float)column-X0AVG));
05993 if(distance > 0.99*RSUNint)
05994 {
05995 image0[i]=NAN;
05996 image1[i]=NAN;
05997 image2[i]=NAN;
05998 image3[i]=NAN;
05999 image4[i]=NAN;
06000 }
06001 }
06002
06003
06004
06005
06006 status=fstats(axisout[0]*axisout[1],arrLev15[0]->data,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
06007 if(status != 0) printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
06008 statusA[0]= drms_setkey_float(recLev15a->records[0],DATAMINS2,(float)minimum);
06009 statusA[1]= drms_setkey_float(recLev15a->records[0],DATAMAXS2,(float)maximum);
06010 statusA[2]= drms_setkey_float(recLev15a->records[0],DATAMEDNS2,(float)median);
06011 statusA[3]= drms_setkey_float(recLev15a->records[0],DATAMEANS2,(float)mean);
06012 statusA[4]= drms_setkey_float(recLev15a->records[0],DATARMSS2,(float)sigma);
06013 statusA[5]= drms_setkey_float(recLev15a->records[0],DATASKEWS2,(float)skewness);
06014 statusA[6]= drms_setkey_float(recLev15a->records[0],DATAKURTS2,(float)kurtosis);
06015 status=fstats(axisout[0]*axisout[1],arrLev15[1]->data,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
06016 if(status != 0) printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
06017 statusA[7]= drms_setkey_float(recLev15b->records[0],DATAMINS2,(float)minimum);
06018 statusA[8]= drms_setkey_float(recLev15b->records[0],DATAMAXS2,(float)maximum);
06019 statusA[9]= drms_setkey_float(recLev15b->records[0],DATAMEDNS2,(float)median);
06020 statusA[10]= drms_setkey_float(recLev15b->records[0],DATAMEANS2,(float)mean);
06021 statusA[11]= drms_setkey_float(recLev15b->records[0],DATARMSS2,(float)sigma);
06022 statusA[12]= drms_setkey_float(recLev15b->records[0],DATASKEWS2,(float)skewness);
06023 statusA[13]= drms_setkey_float(recLev15b->records[0],DATAKURTS2,(float)kurtosis);
06024 status=fstats(axisout[0]*axisout[1],arrLev15[2]->data,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
06025 if(status != 0) printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
06026 statusA[14]= drms_setkey_float(recLev15c->records[0],DATAMINS2,(float)minimum);
06027 statusA[15]= drms_setkey_float(recLev15c->records[0],DATAMAXS2,(float)maximum);
06028 statusA[16]= drms_setkey_float(recLev15c->records[0],DATAMEDNS2,(float)median);
06029 statusA[17]= drms_setkey_float(recLev15c->records[0],DATAMEANS2,(float)mean);
06030 statusA[18]= drms_setkey_float(recLev15c->records[0],DATARMSS2,(float)sigma);
06031 statusA[19]= drms_setkey_float(recLev15c->records[0],DATASKEWS2,(float)skewness);
06032 statusA[20]= drms_setkey_float(recLev15c->records[0],DATAKURTS2,(float)kurtosis);
06033 status=fstats(axisout[0]*axisout[1],arrLev15[3]->data,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
06034 if(status != 0) printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
06035 statusA[21]= drms_setkey_float(recLev15d->records[0],DATAMINS2,(float)minimum);
06036 statusA[22]= drms_setkey_float(recLev15d->records[0],DATAMAXS2,(float)maximum);
06037 statusA[23]= drms_setkey_float(recLev15d->records[0],DATAMEDNS2,(float)median);
06038 statusA[24]= drms_setkey_float(recLev15d->records[0],DATAMEANS2,(float)mean);
06039 statusA[25]= drms_setkey_float(recLev15d->records[0],DATARMSS2,(float)sigma);
06040 statusA[26]= drms_setkey_float(recLev15d->records[0],DATASKEWS2,(float)skewness);
06041 statusA[27]= drms_setkey_float(recLev15d->records[0],DATAKURTS2,(float)kurtosis);
06042 status=fstats(axisout[0]*axisout[1],arrLev15[4]->data,&minimum,&maximum,&median,&mean,&sigma,&skewness,&kurtosis,&ngood);
06043 if(status != 0) printf("Error: the statistics function did not run properly at target time %s\n",timeBegin2);
06044 statusA[28]= drms_setkey_float(recLev15e->records[0],DATAMINS2,(float)minimum);
06045 statusA[29]= drms_setkey_float(recLev15e->records[0],DATAMAXS2,(float)maximum);
06046 statusA[30]= drms_setkey_float(recLev15e->records[0],DATAMEDNS2,(float)median);
06047 statusA[31]= drms_setkey_float(recLev15e->records[0],DATAMEANS2,(float)mean);
06048 statusA[32]= drms_setkey_float(recLev15e->records[0],DATARMSS2,(float)sigma);
06049 statusA[33]= drms_setkey_float(recLev15e->records[0],DATASKEWS2,(float)skewness);
06050 statusA[34]= drms_setkey_float(recLev15e->records[0],DATAKURTS2,(float)kurtosis);
06051
06052 statusA[35]= drms_setkey_string(recLev15a->records[0],HISTORYS,HISTORY);
06053 statusA[36]= drms_setkey_string(recLev15b->records[0],HISTORYS,HISTORY);
06054 statusA[37]= drms_setkey_string(recLev15c->records[0],HISTORYS,HISTORY);
06055 statusA[38]= drms_setkey_string(recLev15d->records[0],HISTORYS,HISTORY);
06056 statusA[39]= drms_setkey_string(recLev15e->records[0],HISTORYS,HISTORY);
06057
06058 TotalStatus=0;
06059 for(i=0;i<35;++i) TotalStatus+=statusA[i];
06060 if(TotalStatus != 0)
06061 {
06062 printf("WARNING: could not set some of the keywords modified by the temporal interpolation subroutine for the Dopplergram at target time %s\n",timeBegin2);
06063 }
06064 t1=dsecnd();
06065 printf("TIME ELAPSED TO SET LEV 1.5 KEYWORDS AND CALCULATE STATISTICS KEYWORDS: %f\n",t1-t0);
06066
06067 }
06068
06069
06070
06071 NextTargetTime:
06072
06073
06074
06075
06076
06077
06078
06079
06080
06081
06082 printf("FREEING RECORD\n");
06083
06084 if(Mask != NULL)
06085 {
06086 free(Mask);
06087 Mask = NULL;
06088 }
06089
06090 if(recLev1d != NULL)
06091 {
06092 printf("recLev1d != NULL\n");
06093 if(recLev1d->n > 0)
06094 {
06095
06096 if (Lev1dWanted) status=drms_close_records(recLev1d,DRMS_INSERT_RECORD);
06097 else status=drms_close_records(recLev1d,DRMS_FREE_RECORD);
06098 recLev1d=NULL;
06099 if(arrLev1d != NULL)
06100 {
06101 for(i=0;i<nRecs1d;++i) if(arrLev1d[i] != NULL)
06102 {
06103 drms_free_array(arrLev1d[i]);
06104 arrLev1d[i]=NULL;
06105 }
06106 free(arrLev1d);
06107 }
06108 arrLev1d=NULL;
06109 Segments1d=0;
06110 }
06111 }
06112
06113
06114 if(recLev1p != NULL)
06115 {
06116 printf("recLev1p != NULL\n");
06117 if(recLev1p->n > 0)
06118 {
06119 if (Lev1pWanted) status=drms_close_records(recLev1p,DRMS_INSERT_RECORD);
06120 else status=drms_close_records(recLev1p,DRMS_FREE_RECORD);
06121 recLev1p=NULL;
06122 for(i=0;i<nSegs1p;++i) if(arrLev1p[i] != NULL)
06123 {
06124 drms_free_array(arrLev1p[i]);
06125 arrLev1p[i]=NULL;
06126 }
06127 if(arrLev1p != NULL) free(arrLev1p);
06128 arrLev1p=NULL;
06129 if(ps1 != NULL) free(ps1);
06130 if(ps2 != NULL) free(ps2);
06131 if(ps3 != NULL) free(ps3);
06132 if(fid != NULL) free(fid);
06133 if(Wavelengths != NULL) free(Wavelengths);
06134 if(images != NULL) free(images);
06135 if(imagesout != NULL) free(imagesout);
06136 Wavelengths=NULL;
06137 images=NULL;
06138 imagesout=NULL;
06139 ps1=NULL;
06140 ps2=NULL;
06141 ps3=NULL;
06142 fid=NULL;
06143 Segments1p=0;
06144 }
06145 }
06146
06147
06148 if(Lev15Wanted)
06149 {
06150
06151 if(arrintable != NULL)
06152 {
06153 drms_free_array(arrintable);
06154 }
06155 arrintable=NULL;
06156 if(lookup != NULL) status=drms_close_records(lookup,DRMS_FREE_RECORD);
06157 lookup = NULL;
06158 if(recpoly != NULL) status=drms_close_records(recpoly,DRMS_FREE_RECORD);
06159 recpoly= NULL;
06160 if(recpoly2 != NULL) status=drms_close_records(recpoly2,DRMS_FREE_RECORD);
06161 recpoly2= NULL;
06162 if(arrayL0 != NULL)
06163 {
06164 drms_free_array(arrayL0);
06165 arrayL0=NULL;
06166 }
06167 if(arrayL1 != NULL)
06168 {
06169 drms_free_array(arrayL1);
06170 arrayL1=NULL;
06171 }
06172 if(arrayL2 != NULL)
06173 {
06174 drms_free_array(arrayL2);
06175 arrayL2=NULL;
06176 }
06177 if(count != NULL)
06178 {
06179 free(count);
06180 count=NULL;
06181 }
06182
06183 if(recLev15a != NULL)
06184 {
06185 printf("recLev15a != NULL\n");
06186 if(CreateEmptyRecord != 1)
06187 {
06188 printf("Inserting record for the observables\n");
06189 status=drms_close_records(recLev15a,DRMS_INSERT_RECORD);
06190 status=drms_close_records(recLev15b,DRMS_INSERT_RECORD);
06191 status=drms_close_records(recLev15c,DRMS_INSERT_RECORD);
06192 status=drms_close_records(recLev15d,DRMS_INSERT_RECORD);
06193 status=drms_close_records(recLev15e,DRMS_INSERT_RECORD);
06194 recLev15a=NULL;
06195 recLev15b=NULL;
06196 recLev15c=NULL;
06197 recLev15d=NULL;
06198 recLev15e=NULL;
06199 }
06200 else
06201 {
06202 printf("Warning: creating empty lev1.5 record\n");
06203
06204 if(CamId == LIGHT_SIDE) camera=1;
06205 if(CamId == LIGHT_FRONT) camera=2;
06206
06207 QUALITY= QUALITY | QUAL_NODATA;
06208 statusA[0] = drms_setkey_time(recLev15a->records[0],TRECS,TargetTime);
06209
06210 statusA[2] = drms_setkey_int(recLev15a->records[0],CAMERAS,camera);
06211 statusA[3] = drms_setkey_int(recLev15a->records[0],QUALITYS,QUALITY);
06212 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
06213 statusA[4]= drms_setkey_string(recLev15a->records[0],DATES,DATEOBS);
06214
06215 statusA[0] = drms_setkey_time(recLev15b->records[0],TRECS,TargetTime);
06216
06217 statusA[2] = drms_setkey_int(recLev15b->records[0],CAMERAS,camera);
06218 statusA[3] = drms_setkey_int(recLev15b->records[0],QUALITYS,QUALITY);
06219 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
06220 statusA[4]= drms_setkey_string(recLev15b->records[0],DATES,DATEOBS);
06221
06222 statusA[0] = drms_setkey_time(recLev15c->records[0],TRECS,TargetTime);
06223
06224 statusA[2] = drms_setkey_int(recLev15c->records[0],CAMERAS,camera);
06225 statusA[3] = drms_setkey_int(recLev15c->records[0],QUALITYS,QUALITY);
06226 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
06227 statusA[4]= drms_setkey_string(recLev15c->records[0],DATES,DATEOBS);
06228
06229 statusA[0] = drms_setkey_time(recLev15d->records[0],TRECS,TargetTime);
06230
06231 statusA[2] = drms_setkey_int(recLev15d->records[0],CAMERAS,camera);
06232 statusA[3] = drms_setkey_int(recLev15d->records[0],QUALITYS,QUALITY);
06233 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
06234 statusA[4]= drms_setkey_string(recLev15d->records[0],DATES,DATEOBS);
06235
06236 statusA[0] = drms_setkey_time(recLev15e->records[0],TRECS,TargetTime);
06237
06238 statusA[2] = drms_setkey_int(recLev15e->records[0],CAMERAS,camera);
06239 statusA[3] = drms_setkey_int(recLev15e->records[0],QUALITYS,QUALITY);
06240 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
06241 statusA[4]= drms_setkey_string(recLev15e->records[0],DATES,DATEOBS);
06242
06243
06244 status=drms_close_records(recLev15a,DRMS_INSERT_RECORD);
06245 status=drms_close_records(recLev15b,DRMS_INSERT_RECORD);
06246 status=drms_close_records(recLev15c,DRMS_INSERT_RECORD);
06247 status=drms_close_records(recLev15d,DRMS_INSERT_RECORD);
06248 status=drms_close_records(recLev15e,DRMS_INSERT_RECORD);
06249 recLev15a=NULL;
06250 recLev15b=NULL;
06251 recLev15c=NULL;
06252 recLev15d=NULL;
06253 recLev15e=NULL;
06254 }
06255 for (i=0;i<nRecs15;++i) if(arrLev15[i] != NULL)
06256 {
06257 drms_free_array(arrLev15[i]);
06258 arrLev15[i]=NULL;
06259 }
06260 if(arrLev15 != NULL) free(arrLev15);
06261 arrLev15=NULL;
06262 }
06263 else
06264 {
06265
06266 recLev15a = drms_create_records(drms_env,1,HMISeriesLev15a,DRMS_PERMANENT,&statusA[0]);
06267 recLev15b = drms_create_records(drms_env,1,HMISeriesLev15b,DRMS_PERMANENT,&statusA[1]);
06268 recLev15c = drms_create_records(drms_env,1,HMISeriesLev15c,DRMS_PERMANENT,&statusA[2]);
06269 recLev15d = drms_create_records(drms_env,1,HMISeriesLev15d,DRMS_PERMANENT,&statusA[3]);
06270 recLev15e = drms_create_records(drms_env,1,HMISeriesLev15e,DRMS_PERMANENT,&statusA[4]);
06271
06272 if(CamId == LIGHT_SIDE) camera=1;
06273 if(CamId == LIGHT_FRONT) camera=2;
06274
06275 printf("Warning: creating empty lev1.5 record\n");
06276 QUALITY= QUALITY | QUAL_NODATA;
06277 statusA[0] = drms_setkey_time(recLev15a->records[0],TRECS,TargetTime);
06278
06279 statusA[2] = drms_setkey_int(recLev15a->records[0],CAMERAS,camera);
06280 statusA[3] = drms_setkey_int(recLev15a->records[0],QUALITYS,QUALITY);
06281 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
06282 statusA[4]= drms_setkey_string(recLev15a->records[0],DATES,DATEOBS);
06283
06284 statusA[0] = drms_setkey_time(recLev15b->records[0],TRECS,TargetTime);
06285
06286 statusA[2] = drms_setkey_int(recLev15b->records[0],CAMERAS,camera);
06287 statusA[3] = drms_setkey_int(recLev15b->records[0],QUALITYS,QUALITY);
06288 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
06289 statusA[4]= drms_setkey_string(recLev15b->records[0],DATES,DATEOBS);
06290
06291 statusA[0] = drms_setkey_time(recLev15c->records[0],TRECS,TargetTime);
06292
06293 statusA[2] = drms_setkey_int(recLev15c->records[0],CAMERAS,camera);
06294 statusA[3] = drms_setkey_int(recLev15c->records[0],QUALITYS,QUALITY);
06295 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
06296 statusA[4]= drms_setkey_string(recLev15c->records[0],DATES,DATEOBS);
06297
06298 statusA[0] = drms_setkey_time(recLev15d->records[0],TRECS,TargetTime);
06299
06300 statusA[2] = drms_setkey_int(recLev15d->records[0],CAMERAS,camera);
06301 statusA[3] = drms_setkey_int(recLev15d->records[0],QUALITYS,QUALITY);
06302 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
06303 statusA[4]= drms_setkey_string(recLev15d->records[0],DATES,DATEOBS);
06304
06305 statusA[0] = drms_setkey_time(recLev15e->records[0],TRECS,TargetTime);
06306
06307 statusA[2] = drms_setkey_int(recLev15e->records[0],CAMERAS,camera);
06308 statusA[3] = drms_setkey_int(recLev15e->records[0],QUALITYS,QUALITY);
06309 sprint_time(DATEOBS,CURRENT_SYSTEM_TIME,"UTC",1);
06310 statusA[4]= drms_setkey_string(recLev15e->records[0],DATES,DATEOBS);
06311
06312 status=drms_close_records(recLev15a,DRMS_INSERT_RECORD);
06313 status=drms_close_records(recLev15b,DRMS_INSERT_RECORD);
06314 status=drms_close_records(recLev15c,DRMS_INSERT_RECORD);
06315 status=drms_close_records(recLev15d,DRMS_INSERT_RECORD);
06316 status=drms_close_records(recLev15e,DRMS_INSERT_RECORD);
06317 recLev15a=NULL;
06318 recLev15b=NULL;
06319 recLev15c=NULL;
06320 recLev15d=NULL;
06321 recLev15e=NULL;
06322 }
06323
06324 QUALITY=0;
06325 CreateEmptyRecord=0;
06326 }
06327
06328 PreviousTargetTime=TargetTime;
06329 TargetTime+=DataCadence;
06330 initialrun=0;
06331 }
06332
06333
06334 if (TestLevIn[0]==1)
06335 {
06336 printf("FREEING GENERAL ARRAYS\n");
06337 if(recLev1->n > 0)
06338 {
06339 status=drms_close_records(recLev1,DRMS_FREE_RECORD);
06340 recLev1=NULL;
06341 free(internTOBS);
06342 free(HWL1POS);
06343 free(HWL2POS);
06344 free(HWL3POS);
06345 free(HWL4POS);
06346 free(HPL1POS);
06347 free(HPL2POS);
06348 free(HPL3POS);
06349 free(HWLTID);
06350 free(HPLTID);
06351 free(FID);
06352 free(HFLID);
06353 free(HCAMID);
06354 free(RSUN);
06355 free(CROTA2);
06356 free(CRLTOBS);
06357 free(DSUNOBS);
06358 free(X0);
06359 free(Y0);
06360 free(SegmentRead);
06361 free(KeywordMissing);
06362 free(Segments);
06363
06364 if(inRotationalFlat == 1)
06365 {
06366 drms_free_array(flatfield);
06367 drms_free_array(flatfieldrot);
06368 status=drms_close_records(recflatrot,DRMS_FREE_RECORD);
06369 }
06370 free(Ierror);
06371 free(IndexFiltergram);
06372 free(FSN);
06373 free(CFINDEX);
06374 free(HIMGCFID);
06375 for(i=0;i<nRecs1;++i) free(IMGTYPE[i]);
06376 free(IMGTYPE);
06377 free(CDELT1);
06378 for(i=0;i<nRecs1;++i) free(HWLTNSET[i]);
06379 free(HWLTNSET);
06380 free(NBADPERM);
06381 free(QUALITYin);
06382 free(QUALITYlev1);
06383 }
06384 }
06385
06386 if(TestLevIn[0]==1)
06387 {
06388 free_interpol(&const_param);
06389 }
06390
06391 if(Lev1pWanted || (Lev15Wanted && TestLevIn[2]==0))
06392 {
06393 status = free_polcal(&pars);
06394 }
06395
06396 free(CODEVERSION);
06397 free(CODEVERSION1);
06398 free(CODEVERSION2);
06399 free(CODEVERSION3);
06400 free(DISTCOEFFILEF);
06401 free(DISTCOEFFILES);
06402 free(ROTCOEFFILE);
06403 free(DISTCOEFPATH);
06404 free(ROTCOEFPATH);
06405
06406
06407 status=0;
06408 t1=dsecnd();
06409 printf("TOTAL TIME ELAPSED IN OBSERVABLES CODE: %f\n",t1-tstart);
06410 return status;
06411
06412 }
06413
06414