; ;************************************************************************************************** ; ;Program: HalphaFlareReports ; ;Build H-alpha Flare Reports ; ;Written by Bill Denig ; ;Version History ; 17 Mar 13 - Initial development ; 22 Jun 14 - Updated I/O ; ;Notes: ; 1. Source data is the monthly FlareData_yyyy-mm-dd.txt files ; ;************************************************************************************************** ; ;************************************************************************************************** ; Pro CMPDate,iYear,iMonth,iDay,xTime,iLat,iLon,iYrCMP,iMnCMP,iDyCMP,xTmCMP ; ;************************************************************************************************** ; ;Program: CMPDate ; ;Routine to calculate the central meridian passage for the Solar Region Report ; ;Written by Bill Denig ; ;Version History ; 05 Mar 12 - Initial development (adopted from Fortran routine by Dan Wilkinson ; ;Notes: TBD ; ;************************************************************************************************** ; xLatRadians = Float(iLat)*3.14159/180. xLonRadians = Float(iLon)*3.14159/180. ; ;Solar rotation speed (degrees/day) calculated using "Astrophysical Quantities" by C.W. Allen ; daysPerDegree = 1./(13.39-2.7*sin(xLatRadians)^2) daysToCMP = Float(iLon) * daysPerDegree ; If (iYear MOD 4) NE 0 Then Begin dayOffset = [ 0, 31, 59, 90,120,151,181,212,243,273,304,334,365] ;Not a leap year EndIf Else Begin dayOffset = [ 0, 31, 60, 91,121,152,182,213,244,274,305,335,366] ;Leap year EndElse ; iYrCMP = iYear iDayOfYear = dayOffset(iMonth-1) + iDay iMinutesOfYear = Long(iDayOfYear-1)*Long(1440) + Long(Fix(1440.*xtime/24. + 0.5)) iMinutesToCMP = Long(daysToCMP*1440) iMinutesOfYearToCMP = iMinutesOfYear + iMinutesToCMP iMaximumMinutesInYear = Long(dayOffset(12))*Long(1440) ;print,'0 ',iYrCMP,iMinutesOfYear,iMinutesToCMP,iMinutesOfYearToCMP ; ;Check if CMP is next year (positive time offset late in year) ; iMinuteOfYearToCMP is the carryover into next year ; If iMinutesOfYearToCMP GE iMaximumMinutesInYear Then Begin iYrCMP = iYrCMP + 1 iMinutesOfYearToCMP = iMinutesOfYearToCMP - iMaximumMinutesInYear EndIf ;print,'1 ',iYrCMP,iMinutesOfYear,iMinutesToCMP,iMinutesOfYearToCMP ; ;Check if CMP is last year (negative time offset early in year) ; iMinuteOfYearToCMP ; If iMinutesOfYearToCMP LT 0 Then Begin iYrCMP = iYrCMP - 1 If iYrCMP MOD 4 NE 0 Then Begin iMinutesOfYearToCMP = Long(365)*Long(1440) + iMinutesOfYear + iMinutesToCMP EndIf Else Begin iMinutesOfYearToCMP = Long(366)*Long(1440) + iMinutesOfYear + iMinutesToCMP EndElse EndIf ;print,'2 ',iYrCMP,iMinutesOfYear,iMinutesToCMP,iMinutesOfYearToCMP ; ;Reset dayOffset for current iYrCMP ; If iYrCMP MOD 4 NE 0 Then Begin dayOffset = [ 0, 31, 59, 90,120,151,181,212,243,273,304,334,365] EndIf Else Begin dayOffset = [ 0, 31, 60, 91,121,152,182,213,244,274,305,335,366] EndElse ; For i=1,12 Do Begin i1 = Long(dayOffset(i-1))*Long(1440) i2 = Long(dayOffset(i ))*Long(1440) If (iMinutesOfYearToCMP GE i1) AND (iMinutesOfYearToCMP LE i2) Then Begin iMinutesOfMonthToCMP = iMinutesOfYearToCMP - i1 xDaysOfMonthToCMP = Float(iMinutesOfMonthToCMP)/1440. iMnCMP = i iDyCMP = Fix(xDaysOfMonthToCMP) + 1 xTmCMP = 24.*(Float(iMinutesOfMonthToCMP)/1440. MOD 1.) EndIf EndFor ;print,'3 ',iYrCMP,iMinutesOfYear,iMinutesToCMP,iMinutesOfYearToCMP ;Print,'3 ',iYear,iMonth,iDay,xTime,iLat,iLon,iYrCMP,iMnCMP,iDyCMP,xTmCMP,Float(iMinutesToCMP/1440.) ; Return End ; ;************************************************************************************************** ; Pro HalphaFlareReports ; ;Input Record Specification: ; ; -------------------------------------------------------------------- ; COLUMNS FMT DESCRIPTION ; -------------------------------------------------------------------- ; 1- 2 I2 Year (yy) ; 3- 4 I2 Month (mm) ; 5- 6 I2 Day (dd) ; 7- 9 A3 Three letter station identifier (sss) ; "HOL" = Station #403; Holloman AFB, NM (KHMN) decode to "HOLL" ; "LEA" = Station #649; Learmonth, Australia (APLM) - decode to "LEAR" ; "SVI" = Station #526; San Vito, Italy (LISS) - decode to "SVTO" ; " " = Station #547; (KANZ) - decode to "KANZ" ; 10 I1 Qualifier - time - TBD ; 11-14 I4 Start time (hhmm) ; 15 I1 Qualifier - time - TBD ; 16-19 I4 Maximim time (hhmm) ; 20 I1 Qualifier - time - TBD ; 21-24 I4 End time (hhmm) ; 25 I1 Quadrant (q = [1/NE, 2/SE , 3/SW , 4/NW) ; 26-27 I2 Central meridian distance (heliographic longitude) ; 28-29 I2 Heliographic latitude ; 30-32 I3 Unknown - TBD ; 33-36 I4 NOAA Region Number ; 37-xx Unknown ; ; -------------------------------------------------------------------- ; ;Output File Specification ; Column Format Description ; ; 1- 2 I2 Data code: always 31 for flares ; 3- 5 I3 Station Code ; "649" = Holloman, NM ; "403" = Learmonth, Australia ; "526" = San Vito, Italy ; 6- 7 I2 Year (yy) ; 8- 9 I2 Month (mm) ; 10-11 I2 Day (dd) ; 12-13 A2 Asterisks mark record with unconfirmed change. ; 14-17 A4 Start time (hhmm) ; 18 A1 Qualifier: D=after,E=before,U=uncertain. ; 19-22 A4 End time (hhmm) ; 23 A1 Qualifier: D=after,E=before,U=uncertain. ; 24-27 A4 Max time (hhmm) ; 28 A1 Qualifier: D=after,E=before,U=uncertain. ; 29 A1 N/North or S/South ; 30-31 I2 Heliographic Latitude ; 32 A1 E/East or W/West ; 33-34 I2 Central meridian distance ; 35 A1 Importance based on flare area = S,1,2 or 3. ; 36 A1 Brightness: F=faint, N=normal, B=bright. ; 37 A1 Completeness = C,P,V or S; indicates the kind of ; observation and the completeness of it. ; 38-41 F4.0 Time flare area was measured. ; 42-46 F5.0 Area; apparent area in millionths of solar disk. ; 47-50 F4.1 Area; corrected area in square degrees; implied ; decimal point between columns 49 and 50. (x10) ; 51 1x -blank- ; 52-56 F5.2 Line width; width in Angstroms of H-alpha line. ; 57-59 A3 Intensity; brightness of H-alpha emission line ; expressed as percentage above continuum. ; 60 A1 X-ray class: C,M,X code the maximum power of 10 ; the 1-8 Angstrom flux attains. ; 61-63 F3.1 X-ray intensity: a number from 1.0 to 9.9 that ; multiplies the X-ray class. ; 64-67 I4 Calcium plage region in which flare occurred. ; 68-71 A4 4-letter station abbreviation ; 72 I1 Seeing: atmospheric stability during observations. ; 73-75 A3 USAF flare reports contain letter-coded remarks ; 76 A1 plus status codes ; 77-80 A4 plus NOAA/USAF region numbers. ; 81-85 I5 NOAA/USAF sunspot region number. ; 86 A1 -blank-; may be used to add a letter to a region. ; 87-88 I2 Central meridian passage - Year (yy) - TBD ; 89-90 I2 Central meridian passage - Month (mm) - TBD ; 91-94 F4.1 Central meridian passage - Day (dd.f) - TBD ; 95 1x -blank- ; 96-100 I5 Grouping number assigned by WDC-A Boulder. ; ; --------------------------------------------------------------------- ; validRecordCounter = 0 rejectRecordCounter = 0 ; ;Get filename and open file ; openFile$: inputFileName = '' On_IOerror,openFileError$ If !version.OS_Name EQ 'Microsoft Windows' Then Begin dataDirectory = 'C:\bdlatitude\Denig - Science\dmsp_sw\ngdc_sw\HalphaFlareReports\' inputFileName = dialog_PickFile(title='Request Input FileName', $ path=dataDirectory,filter = '*.txt',/read, $ get_Path = dataDirectory) EndIf Else Begin Read,Prompt=' Enter input filename ( to end): ',inputFileName EndElse If strLen(inputFileName) EQ 0 Then goTo, endJob$ Print,' Input File: ',inputFileName Get_LUN,u_dpd Openr,u_dpd,inputFileName goTo,readFile$ ; ;Error handler (open) ; openFileError$: On_IOerror,null Free_LUN,u_dpd Print,' ***** File not found - Try again *****' goTo,openFile$ ; readFile$: lineIn = '' On_IOerror,fileError$ Readf,u_dpd,lineIn Print,' Header: ',lineIn ; ;Open outFile ; printFile$: printFileName = '' On_IOerror,printFileError$ If !version.OS_Name EQ 'Microsoft Windows' Then Begin printFileName=dialog_PickFile(title=' Request Print FileName', $ path=dataDirectory,/OverWrite_Prompt,/Write) EndIf Else Begin Read,Prompt=' Enter print filename ( to end): ',printFileName EndElse If strLen(printFileName) EQ 0 Then goTo, endJob$ Print,' Print File: ',printFileName If !version.OS_Name NE 'Microsoft Windows' Then Begin If File_Test(printFileName) Then Begin charAnswer='' Read,Prompt=' *****Overwrite existing output file? Y/N ( to end): ',charAnswer If StrLen(charAnswer) EQ 0 Then goTo,endJob$ charAnswerStrip=StrMid(charAnswer,0,1) If charAnswerStrip EQ 'n' OR charAnswerStrip EQ 'N' Then goTo,printFile$ If charAnswerStrip NE 'y' AND charAnswerStrip NE 'Y' Then Begin Print,' ***** Answer not understood - return to output file request' goTo,printFile$ EndIf EndIf EndIf ; Get_LUN,u_out Openw,u_out,printFileName goTo,processFile$ ; printFileError$: On_IOerror,null Free_LUN,u_out Print,' ***** Error - terminating *****' ; processFile$: While ~ EOF(u_dpd) Do Begin ; Print,'' ;Blank line Readf,u_dpd,lineIn Print,'' If strLen(lineIn) EQ 51 Then Begin Print,'Accepted: ',lineIn EndIf Else Begin Print,'Rejected: ',lineIn rejectRecordCounter = rejectRecordCounter + 1 goTo,NextRecord$ EndElse ; ;Process Record ; ;Catch,theError ; ; 1 2 3 4 5 6 7 8 9 10 lineOut='1234567890123456789012345678901234567890123456789012345678901234567890123456789012345678901234567890' lineOut='zzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzz' ; x01 = strMid(lineIn, 0, 2) ;year (yy) x02 = strMid(lineIn, 2, 2) ;month (mm) x03 = strMid(lineIn, 4, 2) ;day (dd x04 = strMid(lineIn, 6, 3) ;station name x05q= strMid(lineIn, 9, 1) ;qualifier - start time x05 = strMid(lineIn,10, 2) ;start hour (hh) x06 = strMid(lineIn,12, 2) ;start minute (mm) x07q= strMid(lineIn,14, 1) ;qualifier - max time x07 = strMid(lineIn,15, 2) ;maximum hour (hh) x08 = strMid(lineIn,17, 2) ;maximum minute (mm) x09q= strMid(lineIn,19, 1) ;qualifier - end time x09 = strMid(lineIn,20, 2) ;end hour (hh) x10 = strMid(lineIn,22, 2) ;end minute (mm) x11 = strMid(lineIn,24, 1) ;quadrant x12 = strMid(lineIn,25, 2) ;central meridian distance - degrees x13 = strMid(lineIn,27, 2) ;heliographic latitude x14 = strMid(lineIn,29, 2) ;Unknown - TBD x15 = strMid(lineIn,31, 4) ;NOAA region number lineLength = strLen(lineIn) ;Length of the input record x16 = strMid(lineIn,35,lineLength-35) ;Remaining characters ; iYear = Fix(x01) iYear = iYear + 2000 iMonth = Fix(x02) iDay = Fix(x03) stationName = x04 strthh = Fix(x05) strtmm = Fix(x06) peakhh = Fix(x07) peakmm = Fix(x08) stophh = Fix(x09) stopmm = Fix(x10) ; xTime = Float(peakhh) + Float(peakmm)/60. ; ;Determine which solar quadrant the h-alpha flare occurred ; iQuad = fix(x11) ;Convert string to integer xLat = 'X' ;Set default xLon = 'X' ;Set default if iQuad EQ 1 then begin ;Quadrant #1 - N/E xLat = 'N' xLon = 'E' endIf if iQuad EQ 2 then begin ;Quadrant #2 - S/E xLat = 'S' xLon = 'E' endIf if iQuad EQ 3 then begin ;Quadrant #3 - S/W xLat = 'S' xLon = 'W' endIf if iQuad EQ 4 then begin ;Quadrant #4 - N/W xLat = 'N' xLon = 'W' endIf ; iLon = Fix(x12) iLat = Fix(x13) If iQuad EQ 2 OR iQuad EQ 3 Then iLat = -1*iLat If iQuad EQ 3 OR iQuad EQ 4 Then iLon = -1*iLon ; ;Determine station 4-letter designation for observation ; stationID = 'XXXX' ;Default character string starionNm = 'XXX' ;Default character string if x04 EQ 'HOL' then Begin ;Holloman AFB, NM (KHLM) stationID = 'HOLL' stationNm = '649' endIf if x04 EQ 'LEA' then Begin ;Learmonth, Australia (LEAR) stationID = 'LEAR' stationNm = '403' endIf if x04 EQ 'SVI' then Begin ;San Vito, Italy (LISS) stationID = 'SVTO' stationNm = '526' endIf ; if x04 EQ 'KAN' then Begin ;Kan ; stationID = 'KANZ' ; stationNm = '547' ; endIf ; ;NOAA/USAF Region Number (5 digits) ; xRegionNm = 'xxxxx' iRegionNm = fix(x15) if x15 EQ '0000' Then xRegionNm = ' ' ;Blank if no region number if iregionNm GT 0 AND iRegionNm LT 9999 Then Begin strPut,xRegionNm,'1',0 strPut,xRegionNm,x15,1 endIf ; ;Calculate Central Meridian Passage ; CMPDate,iYear,iMonth,iDay,xTime,iLat,iLon,iYrCMP,iMnCMP,iDyCMP,xTmCMP ;Individual CMP ;print,iYear,iMonth,iDay,xTime,iLat,iLon,iYrCMP,iMnCMP,iDyCMP,xTmCMP ;Print ; xYrCMP = 'xx' ;CMP Year If iYrCMP GE 2000 Then iTemp = iYrCMP - 2000 If iYrCMP LT 2000 Then iTemp = iYrCMP - 1900 xTemp=strTrim(string(iTemp),1) if strLen(xTemp) EQ 1 Then Begin strPut,xYrCMP,xTemp,1 strPut,xYrCMP, '0',0 endIf if strLen(xTemp) EQ 2 Then strPut,xYrCMP,xTemp,0 print,'Year: ',iYrCMP,' ',xYrCMP ; xMnCMP = 'xx' ;CMP Month xTemp=strTrim(string(iMnCMP),1) if strLen(xTemp) EQ 1 Then Begin strPut,xMnCMP,xTemp,1 strPut,xMnCMP, '0',0 endIf if strLen(xTemp) EQ 2 Then strPut,xMnCMP,xTemp,0 ;print,'Month: ',iMnCMP,' ',xMnCMP ; xDyCMP = 'xx.x' ;CMP Day xTemp=strTrim(string(iDyCMP),1) if strLen(xTemp) EQ 1 Then Begin strPut,xDyCMP,xTemp,1 strPut,xDyCMP, '0',0 endIf if strLen(xTemp) EQ 2 Then strPut,xDyCMP,xTemp,0 ; iFrac = Fix(10.*xTmCMP/24.) ;xTmCMP is day fraction - uses truncation xTemp=strTrim(string(iFrac),1) if strLen(xTemp) EQ 1 Then strPut,xDyCMP,xTemp,3 ;print,'Day: ',iDyCMP,xTmCMP,' ',xDyCMP ; strPut,lineOut, '31' , 0 ;Data code; always '31' for flare events strPut,lineOut,stationNm,2 ;Station code strPut,lineOut, x01 , 5 ;Year (A2) strPut,lineOut, x02 , 7 ;Month (A2) strPut,lineOut, x03 , 9 ;Day (A2) strPut,lineOut, ' ' , 11 ;Unconfirmed change (TBD) strPut,lineOut, x05 , 13 ;Start hour (A2) strPut,lineOut, x06 , 15 ;Start minute (A2) strPut,lineOut, '-' , 17 ;-qualifier--start time - TBD (X05q) strPut,lineOut, x09 , 18 ;End hour (A2) strPut,lineOut, x10 , 20 ;End minute (A2) strPut,lineOut, '-' , 22 ;-qualifier----end time - TBD (x09q) strPut,lineOut, x07 , 23 ;Peak hour (A2) strPut,lineOut, x08 , 25 ;Peak minute (A2) strPut,lineOut, '-' , 27 ;-qualifier---peak time - TBD (x07q) strPut,lineOut, xLat , 28 ;Latitude - N/S; Default = 'X' strPut,lineOut, x13 , 29 ;Latitude - value strPut,lineOut, xLon , 31 ;Longitude - E/W; Default = 'X' strPut,lineOut, x12 , 32 ;Longitude - value strPut,lineOut, '-' , 34 ;Flare importance - TBD strPut,lineOut, '-' , 35 ;Flare brightness - TBD strPut,lineOut, '-' , 36 ;Flare completeness - TBD strPut,lineOut,'----', 37 ;Time when flare area measured - TBD strPut,lineOut,'-----',41 ;Flare area - measured - millionths of solar disk strPut,lineOut,'----', 46 ;Flare area - corrected - square degrees strPut,lineOut, ' ', 50 ;-blank- strPut,lineOut,'-----',51 ;Line width of H-alpha emission - angstroms strPut,lineOut, '---', 56 ;Intensity of H-alpha emission above background - percent strPut,lineOut, '-', 59 ;X-ray class - C/M/X strPut,lineOut, '---', 60 ;X-ray intensity - 0.0 -to- 9.9 strPut,lineOut,'----', 63 ;Calcium plage region strPut,lineOut,stationID,67 ;Station 4-letter strPut,lineOut, '-', 71 ;Seeing: Atmospheric Stability strPut,lineOut, '---', 72 ;Remarts: USAF flare report strPut,lineOut, '-', 75 ; Status strPut,lineOut,'----', 76 ; NOAA/USAF region number strPut,lineOut,xRegionNm,80 ;NOAA/USAF region number - 1st digit = 1 (after xx-xx-20xx) strPut,lineOut, ' ' , 85 ;-blank-; may be used to add a letter to a region strPut,lineOut,xYrCMP, 86 ;Central meridian passage - year (yy) - I2 strPut,lineOut,xMnCMP, 88 ;Central meridian passage - month (mm) - I2 strPut,lineOut,xDyCMP, 90 ;Central meridian passage - day (dd) - F4.1 strPut,lineOut, ' ' , 94 ;-blank- strPut,lineOut,'-----',95 ;Grouping Number (assigned by WDC-A) Print,lineOut ; validRecordCounter = validRecordCounter + 1 ; ;Print record ; print,'OutF: ',lineOut printf,u_Out,lineOut nextRecord$: EndWhile ; Close,u_dpd Close,u_out goTo,organizeJob$ ; fileError$: On_IOerror,null Print,' ***** Unexpected File error - terminating' Free_LUN,u_dpd Free_LUN,u_out goTo,endJob$ ; organizeJob$: Print,' Number of Valid Records: ',validRecordCounter Print,' Number of Rejected Records: ',rejectRecordCounter Print,'' ; continueJob$: charAnswer='' ; endJob$: Print,' End of Job' End